GLpApp_GLpYUEpetraDataPool.hpp

Go to the documentation of this file.
00001 #ifndef GLPAPP_GLPYUEPETRADATAPOOL_H
00002 #define GLPAPP_GLPYUEPETRADATAPOOL_H
00003 
00004 //#include "Epetra_config.h"
00005 
00006 #include <iostream>
00007 
00008 #include "Epetra_Map.h"
00009 #include "Epetra_MultiVector.h"
00010 #include "Epetra_Vector.h"
00011 #include "Epetra_Import.h"
00012 #include "Epetra_Export.h"
00013 #include "Epetra_CrsMatrix.h"
00014 #include "Epetra_FECrsMatrix.h"
00015 #include "Epetra_LinearProblem.h"
00016 #include "Epetra_LAPACK.h"
00017 #include "Epetra_FEVector.h"
00018 #include "Epetra_IntSerialDenseVector.h"
00019 #include "Epetra_SerialDenseMatrix.h"
00020 #include "Epetra_SerialDenseVector.h"
00021 #include "GenSQP_DataPool.hpp"
00022 #include "GenSQP_YUEpetraVector.hpp"
00023 #include "Epetra_SerialDenseMatrix.h"
00024 #include "Epetra_SerialDenseVector.h"
00025 
00026 #ifdef HAVE_MPI
00027 #  include "Epetra_MpiComm.h"
00028 #else
00029 #  include "Epetra_SerialComm.h"
00030 #endif
00031 
00032 namespace GLpApp {
00033 
00034 class GLpYUEpetraDataPool : public GenSQP::DataPool
00035 {
00036 public:
00037 
00038   GLpYUEpetraDataPool(
00039     Teuchos::RefCountPtr<const Epetra_Comm>    const& commptr
00040     ,const double                              beta
00041     ,const double                              len_x     // Ignored if myfile is *not* empty
00042     ,const double                              len_y     // Ignored if myfile is *not* empty
00043     ,const int                                 local_nx  // Ignored if myfile is *not* empty
00044     ,const int                                 local_ny  // Ignored if myfile is *not* empty
00045     ,const char                                myfile[]
00046     ,const bool                                trace
00047     );
00048 
00052   void computeAll( const GenSQP::Vector &x );
00053 
00055   int  solveAugsys( const Teuchos::RefCountPtr<const Epetra_MultiVector> & rhsy,
00056                     const Teuchos::RefCountPtr<const Epetra_MultiVector> & rhsu,
00057                     const Teuchos::RefCountPtr<const Epetra_MultiVector> & rhsp,
00058                     const Teuchos::RefCountPtr<Epetra_MultiVector> & y,
00059                     const Teuchos::RefCountPtr<Epetra_MultiVector> & u,
00060                     const Teuchos::RefCountPtr<Epetra_MultiVector> & p,
00061                     double tol );
00062 
00063   Teuchos::RefCountPtr<const Epetra_Comm> getCommPtr();
00064 
00065   Teuchos::RefCountPtr<Epetra_FECrsMatrix> getA();
00066   Teuchos::RefCountPtr<Epetra_FECrsMatrix> getB();
00067   Teuchos::RefCountPtr<Epetra_FECrsMatrix> getH();
00068   Teuchos::RefCountPtr<Epetra_FECrsMatrix> getR();
00069   Teuchos::RefCountPtr<Epetra_CrsMatrix> getAugmat();
00070   Teuchos::RefCountPtr<Epetra_FECrsMatrix> getNpy();
00071 
00072   Teuchos::RefCountPtr<Epetra_FEVector> getb();
00073   Teuchos::RefCountPtr<Epetra_FEVector> getq();
00074   Teuchos::RefCountPtr<Epetra_FEVector> getNy();
00075 
00077   void computeNy(const Teuchos::RefCountPtr<const Epetra_MultiVector> & y);
00078 
00080   void computeNpy(const Teuchos::RefCountPtr<const Epetra_MultiVector> & y);
00081 
00083   void computeAugmat();
00084   
00085   Teuchos::RefCountPtr<const Epetra_SerialDenseMatrix> getipcoords();
00086   Teuchos::RefCountPtr<const Epetra_IntSerialDenseVector> getipindx();
00087   Teuchos::RefCountPtr<const Epetra_SerialDenseMatrix> getpcoords();
00088   Teuchos::RefCountPtr<const Epetra_IntSerialDenseVector> getpindx();
00089   Teuchos::RefCountPtr<const Epetra_IntSerialDenseMatrix> gett();
00090   Teuchos::RefCountPtr<const Epetra_IntSerialDenseMatrix> gete();
00091 
00092   double getbeta();
00093   
00095   void PrintVec( const Teuchos::RefCountPtr<const Epetra_Vector> & x );
00096 
00097 private:
00098 
00099   Teuchos::RefCountPtr<const Epetra_Comm> commptr_;
00100           
00102   Teuchos::RefCountPtr<Epetra_SerialDenseMatrix> ipcoords_;
00104   Teuchos::RefCountPtr<Epetra_IntSerialDenseVector> ipindx_;
00106   Teuchos::RefCountPtr<Epetra_SerialDenseMatrix> pcoords_;
00108   Teuchos::RefCountPtr<Epetra_IntSerialDenseVector> pindx_;
00110   Teuchos::RefCountPtr<Epetra_IntSerialDenseMatrix> t_;
00112   Teuchos::RefCountPtr<Epetra_IntSerialDenseMatrix> e_;
00113 
00115   Teuchos::RefCountPtr<Epetra_FECrsMatrix> A_;
00117   Teuchos::RefCountPtr<Epetra_FECrsMatrix> B_;
00119   Teuchos::RefCountPtr<Epetra_FECrsMatrix> H_;
00121   Teuchos::RefCountPtr<Epetra_FECrsMatrix> R_;
00122 
00124   Teuchos::RefCountPtr<Epetra_MultiVector> B_bar_;
00125 
00130   Teuchos::RefCountPtr<Epetra_CrsMatrix> Augmat_;
00131 
00133   Teuchos::RefCountPtr<Epetra_FECrsMatrix> Npy_;
00134 
00136   Teuchos::RefCountPtr<Epetra_FEVector> b_;
00138   Teuchos::RefCountPtr<Epetra_FEVector> q_;
00139 
00140   Teuchos::RefCountPtr<Epetra_FEVector> Ny_;
00141 
00143   double beta_;
00144 
00145 };
00146 
00147 class Usr_Par {
00148 public:
00149 
00150   Usr_Par();
00151 
00152   Epetra_SerialDenseMatrix Nodes;
00153   Epetra_SerialDenseVector Weights;
00154 
00155   Epetra_SerialDenseMatrix N;
00156 
00157   Epetra_SerialDenseMatrix Nx1;
00158 
00159   Epetra_SerialDenseMatrix Nx2;
00160 
00161   Epetra_SerialDenseMatrix S1;
00162   Epetra_SerialDenseMatrix S2;
00163   Epetra_SerialDenseMatrix S3;
00164 
00165   Epetra_SerialDenseVector Nw;
00166 
00167   Epetra_SerialDenseMatrix NNw;
00168 
00169   Epetra_SerialDenseMatrix * NNNw;
00170 
00171   Epetra_SerialDenseMatrix * NdNdx1Nw;
00172 
00173   Epetra_SerialDenseMatrix * NdNdx2Nw;
00174 
00175   ~Usr_Par() {
00176     delete [] NNNw;
00177     delete [] NdNdx1Nw;
00178     delete [] NdNdx2Nw;
00179   }
00180   
00181   void Print(ostream& os) const;
00182 };
00183 
00184 } // namespace GLpApp
00185 
00186 #endif // GLPAPP_GLPYUEPETRADATAPOOL_H

Generated on Wed May 12 21:24:46 2010 for EpetraExt by  doxygen 1.4.7