FEI Version of the Day
fei_Aztec_LinSysCore.hpp
00001 #ifndef _fei_Aztec_LinSysCore_hpp_
00002 #define _fei_Aztec_LinSysCore_hpp_
00003 
00004 /*--------------------------------------------------------------------*/
00005 /*    Copyright 2005 Sandia Corporation.                              */
00006 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00007 /*    non-exclusive license for use of this work by or on behalf      */
00008 /*    of the U.S. Government.  Export of this program may require     */
00009 /*    a license from the United States Government.                    */
00010 /*--------------------------------------------------------------------*/
00011 
00012 #include <fei_macros.hpp>
00013 #include <fei_defs.h>
00014 #include <fei_Data.hpp>
00015 #include <fei_LinearSystemCore.hpp>
00016 #include <fei_SharedPtr.hpp>
00017 
00018 #include <string>
00019 #include <map>
00020 
00021 #include <az_aztec.h>
00022 
00023 //
00024 //This is the Aztec implementation of LinSysCore.
00025 //
00026 namespace fei_trilinos {
00027 
00028 class Aztec_Map;
00029 class Aztec_BlockMap;
00030 class Aztec_LSVector;
00031 class AztecDMSR_Matrix;
00032 class AztecDVBR_Matrix;
00033 
00034 class Aztec_LinSysCore: public LinearSystemCore {
00035  public:
00036    Aztec_LinSysCore(MPI_Comm comm);
00037    virtual ~Aztec_LinSysCore();
00038 
00039    //for creating another instance of LinearSystemCore without knowing
00040    //the run-time type of 'this' one.
00041    LinearSystemCore* clone();
00042 
00043    //int parameters:
00044    //for setting generic argc/argv style parameters.
00045 
00046    int parameters(int numParams, const char*const * params);
00047 
00048    int setLookup(Lookup& lookup);
00049 
00050    int setGlobalOffsets(int len, int* nodeOffsets, int* eqnOffsets,
00051                                  int* blkEqnOffsets);
00052 
00053    int setConnectivities(GlobalID elemBlock,
00054                           int numElements,
00055                           int numNodesPerElem,
00056                           const GlobalID* elemIDs,
00057                           const int* const* connNodes) ;
00058 
00059    int setStiffnessMatrices(GlobalID,
00060           int,
00061           const GlobalID*,
00062           const double *const *const *,
00063           int,
00064           const int *const *)
00065      { return(0); }
00066 
00067    int setLoadVectors(GlobalID,
00068           int,
00069           const GlobalID*,
00070           const double *const *,
00071           int,
00072           const int *const *)
00073      { return(0); }
00074 
00075    int setMatrixStructure(int** ptColIndices,
00076                            int* ptRowLengths,
00077                            int** blkColIndices,
00078                            int* blkRowLengths,
00079                            int* ptRowsPerBlkRow) ;
00080 
00081    int setMultCREqns(int,
00082          int, int,
00083          int**, int**,
00084          int*,
00085          int*)
00086      { return(0); }
00087 
00088    int setPenCREqns(int,
00089          int, int,
00090          int**, int**,
00091          int*)
00092      { return(0); }
00093 
00094    //int resetMatrixAndVector:
00095    //don't destroy the structure of the matrix, but set the value 's'
00096    //throughout the matrix and vectors.
00097 
00098    int resetMatrixAndVector(double s);
00099    int resetMatrix(double s);
00100    int resetRHSVector(double s);
00101 
00102    //int sumIntoSystemMatrix:
00103    //this is the primary assembly function. The coefficients 'values'
00104    //are to be accumumlated into (added to any values already in place)
00105    //global (0-based) equations in 'ptRows' of the matrix.
00106 
00107    int sumIntoSystemMatrix(int numPtRows, const int* ptRows,
00108                             int numPtCols, const int* ptCols,
00109                             int numBlkRows, const int* blkRows,
00110                             int numBlkCols, const int* blkCols,
00111                             const double* const* values);
00112 
00113    int sumIntoSystemMatrix(int numPtRows, const int* ptRows,
00114                             int numPtCols, const int* ptCols,
00115                             const double* const* values);
00116    int putIntoSystemMatrix(int numPtRows, const int* ptRows,
00117                             int numPtCols, const int* ptCols,
00118                             const double* const* values);
00119 
00120    int getMatrixRowLength(int row, int& length);
00121 
00122    int getMatrixRow(int row, double* coefs, int* indices,
00123                             int len, int& rowLength);
00124 
00125    //int sumIntoRHSVector:
00126    //this is the rhs vector equivalent to sumIntoSystemMatrix above.
00127 
00128    int sumIntoRHSVector(int num,
00129                           const double* values,
00130                           const int* indices);
00131    int putIntoRHSVector(int num,
00132                           const double* values,
00133                           const int* indices);
00134    int getFromRHSVector(int num,
00135                           double* values,
00136                           const int* indices);
00137 
00138    //int matrixLoadComplete:
00139    //do any internal synchronization/communication.
00140 
00141    int matrixLoadComplete();
00142    
00143    //functions for enforcing boundary conditions.
00144    int enforceEssentialBC(int* globalEqn,
00145                            double* alpha,
00146                            double* gamma, int len);
00147 
00148    int enforceBlkEssentialBC(int* blkEqn, int* blkOffset,
00149                               double* alpha, double* gamma,
00150                               int len);
00151 
00152    int enforceRemoteEssBCs(int numEqns, int* globalEqns,
00153                             int** colIndices, int* colIndLen,
00154                             double** coefs);
00155 
00156    int enforceBlkRemoteEssBCs(int numEqns, int* blkEqns,
00157                                int** blkColInds, int** blkColOffsets,
00158                                int* blkColLens, double** remEssBCCoefs);
00159 
00160    //functions for getting/setting matrix or vector pointers.
00161 
00162    //getMatrixPtr:
00163    //obtain a pointer to the 'A' matrix. This should be considered a
00164    //constant pointer -- i.e., this class remains responsible for the
00165    //matrix (e.g., de-allocation upon destruction). 
00166    int getMatrixPtr(Data& data);
00167 
00168    //copyInMatrix:
00169    //replaces the internal matrix with a copy of the input argument, scaled
00170    //by the coefficient 'scalar'.
00171 
00172    int copyInMatrix(double scalar, const Data& data);
00173 
00174    //copyOutMatrix:
00175    //passes out a copy of the internal matrix, scaled by the coefficient
00176    //'scalar'.
00177 
00178    int copyOutMatrix(double scalar, Data& data);
00179 
00180    //sumInMatrix:
00181    //accumulate (sum) a copy of the input argument into the internal
00182    //matrix, scaling the input by the coefficient 'scalar'.
00183 
00184    int sumInMatrix(double scalar, const Data& data);
00185 
00186    //get/setRHSVectorPtr:
00187    //the same semantics apply here as for the matrixPtr functions above.
00188 
00189    int getRHSVectorPtr(Data& data);
00190 
00191    //copyInRHSVector/copyOutRHSVector/sumInRHSVector:
00192    //the same semantics apply here as for the matrix functions above.
00193 
00194    int copyInRHSVector(double scalar, const Data& data);
00195    int copyOutRHSVector(double scalar, Data& data);
00196    int sumInRHSVector(double scalar, const Data& data);
00197 
00198    //destroyMatrixData/destroyVectorData:
00199    //Utility function for destroying the matrix (or vector) in Data
00200 
00201    int destroyMatrixData(Data& data);
00202    int destroyVectorData(Data& data);
00203 
00204    //functions for managing multiple rhs vectors
00205    int setNumRHSVectors(int numRHSs, const int* rhsIDs);
00206 
00207    //int setRHSID:
00208    //set the 'current' rhs context, assuming multiple rhs vectors.
00209    int setRHSID(int rhsID);
00210 
00211    //int putInitialGuess:
00212    //function for setting (a subset of) the initial-guess
00213    //solution values (i.e., in the 'x' vector).
00214 
00215    int putInitialGuess(const int* eqnNumbers, const double* values,
00216                         int len);
00217 
00218    //function for getting all of the answers ('x' vector).
00219    int getSolution(double* answers, int len);
00220 
00221    //function for getting the (single) entry at equation
00222    //number 'eqnNumber'.
00223    int getSolnEntry(int eqnNumber, double& answer);
00224 
00225    //function for obtaining the residual vector (using current solution or
00226    //initial guess)
00227    int formResidual(double* values, int len);
00228 
00229    //function for launching the linear solver
00230    int launchSolver(int& solveStatus, int& iterations);
00231 
00232    int putNodalFieldData(int, int, int*, int, const double*)
00233      { return(0); }
00234 
00235    int writeSystem(const char* name);
00236 
00237  private:        //functions
00238 
00239    int createMiscStuff();
00240 
00241    int allocateMatrix(int** ptColIndices, int* ptRowLengths,
00242                        int** blkColIndices, int* blkRowLengths,
00243                        int* ptRowsPerBlkRow);
00244 
00245    int VBRmatPlusScaledMat(AztecDVBR_Matrix* A, double scalar,
00246                             AztecDVBR_Matrix* source);
00247 
00248    int MSRmatPlusScaledMat(AztecDMSR_Matrix* A, double scalar,
00249                             AztecDMSR_Matrix* source);
00250 
00251    int createBlockMatrix(int** blkColIndices,
00252                           int* blkRowLengths,
00253                           int* ptRowsPerBlkRow);
00254 
00255    int sumIntoBlockRow(int numBlkRows, const int* blkRows,
00256                         int numBlkCols, const int* blkCols,
00257                         const double* const* values,
00258                         int numPtCols,
00259            bool overwriteInsteadOfAccumulate);
00260 
00261    int copyBlockRow(int i, const int* blkRows,
00262                      int numBlkCols, const int* blkCols,
00263                      const double* const* values,
00264                      double* coefs);
00265 
00266    int modifyRHSforBCs();
00267 
00268    int explicitlySetDirichletBCs();
00269 
00270    int blockRowToPointRow(int blkRow);
00271 
00272    int getBlockRow(int blkRow, double*& val, int& valLen,
00273                    int*& blkColInds, int& blkColIndLen,
00274                    int& numNzBlks, int& numNNZ);
00275 
00276    int getBlkEqnsAndOffsets(int* ptEqns, int* blkEqns, int* blkOffsets,
00277                              int numEqns);
00278 
00279    int getBlockSize(int blkInd);
00280 
00281    int sumIntoPointRow(int numPtRows, const int* ptRows,
00282            int numPtCols, const int* ptColIndices,
00283            const double* const* values,
00284            bool overwriteInsteadOfAccumulate);
00285 
00286    int sumPointIntoBlockRow(int blkRow, int rowOffset,
00287           int blkCol, int colOffset, double value);
00288 
00289    int setMatrixType(const char* name);
00290    int selectSolver(const char* name);
00291    int selectPreconditioner(const char* name);
00292    void setSubdomainSolve(const char* name);
00293    void setScalingOption(const char* param);
00294    void setConvTest(const char* param);
00295    void setPreCalc(const char* param);
00296    void setTypeOverlap(const char* param);
00297    void setOverlap(const char* param);
00298    void setOrthog(const char* param);
00299    void setAuxVec(const char* param);
00300    void setAZ_output(const char* param);
00301 
00302    void recordUserParams();
00303 
00304    void checkForParam(const char* paramName, int numParams_,
00305                       char** paramStrings,
00306                       double& param);
00307 
00308    void recordUserOptions();
00309 
00310    void checkForOption(const char* paramName, int numParams_,
00311                        char** paramStrings,
00312                        int& param);
00313 
00314    int blkRowEssBCMod(int blkEqn, int blkOffset, double* val, int* blkCols,
00315                         int numCols, int numPtNNZ, double alpha, double gamma);
00316 
00317    int blkColEssBCMod(int blkRow, int blkEqn, int blkOffset, double* val,
00318                       int* blkCols,
00319                       int numCols, int numPtNNZ, double alpha, double gamma);
00320 
00321    void setDebugOutput(const char* path, const char* name);
00322 
00323    void debugOutput(const char* msg) const;
00324 
00325    int writeA(const char* name);
00326    int writeVec(Aztec_LSVector* v, const char* name);
00327 
00328    int messageAbort(const char* msg) const;
00329 
00330  private:            //variables
00331 
00332    MPI_Comm comm_;
00333 
00334    Lookup* lookup_;
00335    bool haveLookup_;
00336 
00337    int numProcs_;
00338    int thisProc_;
00339    int masterProc_;
00340 
00341    int* update_;
00342    fei::SharedPtr<Aztec_Map> map_;
00343    AztecDMSR_Matrix *A_;
00344    AztecDMSR_Matrix *A_ptr_;
00345    Aztec_LSVector *x_, **b_, *bc_;
00346    int* essBCindices_;
00347    int numEssBCs_;
00348    bool bcsLoaded_;
00349    bool explicitDirichletBCs_;
00350    bool BCenforcement_no_column_mod_;
00351    Aztec_LSVector *b_ptr_;
00352    bool matrixAllocated_;
00353    bool vectorsAllocated_;
00354    bool blkMatrixAllocated_;
00355    bool matrixLoaded_;
00356    bool rhsLoaded_;
00357    bool needNewPreconditioner_;
00358 
00359    bool tooLateToChooseBlock_;
00360    bool blockMatrix_;
00361    fei::SharedPtr<Aztec_BlockMap> blkMap_;
00362    AztecDVBR_Matrix *blkA_;
00363    AztecDVBR_Matrix *blkA_ptr_;
00364    int* blkUpdate_;
00365 
00366    AZ_MATRIX *azA_;
00367    AZ_PRECOND *azP_;
00368    bool precondCreated_;
00369    AZ_SCALING *azS_;
00370    bool scalingCreated_;
00371 
00372    int *aztec_options_;
00373    double *aztec_params_;
00374    double *aztec_status_;
00375 
00376    double* tmp_x_;
00377    bool tmp_x_touched_;
00378    double** tmp_b_;
00379    double* tmp_bc_;
00380    bool tmp_b_allocated_;
00381 
00382    bool ML_Vanek_;
00383    int numLevels_;
00384 
00385    int* rhsIDs_;
00386    int numRHSs_;
00387 
00388    int currentRHS_;
00389 
00390    int numGlobalEqns_;
00391    int localOffset_;
00392    int numLocalEqns_;
00393 
00394    int numGlobalEqnBlks_;
00395    int localBlkOffset_;
00396    int numLocalEqnBlks_;
00397    int* localBlockSizes_;
00398 
00399    int numNonzeroBlocks_;
00400 
00401    int outputLevel_;
00402    int numParams_;
00403    char** paramStrings_;
00404 
00405    std::string name_;
00406    int debugOutput_;
00407    int debugFileCounter_;
00408    char* debugPath_;
00409    char* debugFileName_;
00410    FILE* debugFile_;
00411 
00412    std::map<std::string,unsigned>& named_solve_counter_;
00413 };
00414 
00415 }//namespace fei_trilinos
00416 
00417 #endif
00418 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends