• Main Page
  • Related Pages
  • Namespaces
  • Classes
  • Files
  • File List
  • File Members

ml_epetra_utils.h

Go to the documentation of this file.
00001 
00009 /* ******************************************************************** */
00010 /* See the file COPYRIGHT for a complete copyright notice, contact      */
00011 /* person and disclaimer.                                               */        
00012 /* ******************************************************************** */
00013 
00014 
00015 #ifndef _ML_EPETRA_UTILS_H_
00016 #define _ML_EPETRA_UTILS_H_
00017 
00018 class Epetra_Comm;
00019 class Epetra_BlockMap;
00020 class Epetra_MultiVector;
00021 class Epetra_RowMatrix;
00022 class Epetra_Map;
00023 class Epetra_Vector;
00024 class Epetra_IntVector;
00025 class Epetra_Import;
00026 class Epetra_Object;
00027 class Epetra_CrsGraph;
00028 class Epetra_CrsMatrix;
00029 class Epetra_RowMatrix;
00030 class Epetra_LinearProblem;
00031 class Epetra_SerialDenseMatrix;
00032 namespace Teuchos {
00033   class ParameterList;
00034 }
00035 
00036 #include "ml_common.h"
00037 
00038 #ifdef ML_MPI
00039 #ifndef EPETRA_MPI
00040 #define EPETRA_MPI
00041 #endif
00042 #include "mpi.h"
00043 #endif
00044 #include "ml_include.h"
00045 #include <iostream>
00046 
00047 #ifdef HAVE_ML_EPETRAEXT
00048 #include "EpetraExt_SolverMap_CrsMatrix.h"
00049 #endif
00050 
00051 
00052 #ifndef ML_CPP
00053 extern "C" {
00054 #endif
00055 
00056 #include "Epetra_DataAccess.h"
00057   
00058 // ====================================================================== 
00060 
00062 int ML_Epetra_matvec(ML_Operator *data, int in, double *p,
00063                                         int out, double *ap);
00064 int ML_Epetra_RowMatrix_matvec(ML_Operator *data, int in, double *p,
00065                                                   int out, double *ap);
00066 int ML_Epetra_CrsMatrix_matvec(ML_Operator *data, int in, double *p,
00067                                                   int out, double *ap);
00068 int ML_Epetra_VbrMatrix_matvec(ML_Operator *data, int in, double *p,
00069                                                   int out, double *ap);
00070 int Epetra_ML_GetCrsDataptrs(ML_Operator *data, double **values, int **cols, int **rowptr);
00071 
00072 #ifdef WKC
00073 int ML_Epetra_matvec_WKC(ML_Operator *data, int in, double *p, int out,
00074                  double *ap);
00075 #endif 
00076 
00078 
00105 int ML_Epetra_getrow(ML_Operator *data, int N_requested_rows,
00106                  int requested_rows[], int allocated_space, int columns[],
00107                  double values[], int row_lengths[]);
00108 
00109 int ML_Epetra_RowMatrix_getrow(ML_Operator *data, int N_requested_rows, 
00110                                int requested_rows[], int allocated_space, 
00111                                int columns[], double values[],
00112                                int row_lengths[]);
00113   
00114 int ML_Epetra_CrsMatrix_getrow(ML_Operator *data, int N_requested_rows,
00115             int requested_rows[], 
00116         int allocated_space, int columns[], double values[],
00117                                int row_lengths[]);
00118 int ML_Epetra_CrsMatrix_get_one_row(ML_Operator *data, int N_requested_rows,
00119                                     int requested_rows[], 
00120                                     int allocated_space, int columns[], double values[],
00121                                     int row_lengths[]);
00122   
00123 
00124 int ML_Epetra_VbrMatrix_getrow(ML_Operator *data,
00125             int N_requested_rows, int requested_rows[], 
00126         int allocated_space, int columns[], double values[],
00127         int row_lengths[]);
00128 
00129 void ML_Set_Filter(Teuchos::ParameterList& List);
00130 
00131 int ML_Epetra_matvec_Filter(ML_Operator *data, int in, double *p, 
00132                             int out, double *ap);
00133 
00134 int ML_Epetra_getrow_Filter(ML_Operator *data, int N_requested_rows,
00135                             int requested_rows[], int allocated_space, int columns[],
00136                             double values[], int row_lengths[]);
00138 
00150 int ML_Epetra_comm_wrapper(double vec[], void *data);
00151 int ML_Epetra_CrsMatrix_comm_wrapper(double vec[], void *data);
00152 int ML_Epetra_VbrMatrix_comm_wrapper(double vec[], void *data);
00153 
00154 // wrappers for Epetra_CrsGraph
00155 int ML_Epetra_CrsGraph_comm_wrapper(double vec[], void *data);
00156 int ML_Epetra_CrsGraph_matvec(ML_Operator *data, int in, double *p,
00157                               int out, double *ap);
00158 int ML_Epetra_CrsGraph_getrow(ML_Operator *data, int N_requested_rows,
00159                               int requested_rows[], int allocated_space, 
00160                               int columns[], double values[],
00161                               int row_lengths[]);
00162 int ML_Operator_WrapEpetraCrsGraph(Epetra_CrsGraph* Graph, ML_Operator *newMatrix);
00163 
00164 #ifdef HAVE_ML_TEUCHOS
00165 void ML_CreateSublists(Teuchos::ParameterList &List,
00166                       Teuchos::ParameterList &newList,
00167                       int NumLevels);
00168 #endif
00169 
00170 #ifndef ML_CPP
00171 }
00172 #endif
00173 
00175 
00182 int EpetraMatrix2MLMatrix(ML *ml_handle, int level,
00183                                 Epetra_RowMatrix * Amat);
00184 
00186 
00194 int ML_Operator_WrapEpetraMatrix(Epetra_RowMatrix * A, ML_Operator *Result);
00195 
00196 
00198 
00209 int ML_Operator_WrapEpetraCrsMatrix(Epetra_CrsMatrix * A, ML_Operator *newMatrix, bool verbose=false);
00210 
00212 
00214 void Epetra_CrsMatrix_Wrap_ML_Operator(ML_Operator * A, const Epetra_Comm &Comm, const Epetra_Map &RowMap,Epetra_CrsMatrix **Result,Epetra_DataAccess CV=View,int base=0);
00215 //void Epetra_CrsMatrix_Wrap_ML_Operator(ML_Operator * A, Epetra_CrsMatrix *Result);
00216 
00217 
00219 Epetra_CrsMatrix *Epetra_MatrixMult(Epetra_RowMatrix *B, Epetra_RowMatrix *Bt);
00220 
00222 Epetra_CrsMatrix *Epetra_MatrixAdd(Epetra_RowMatrix *B, Epetra_RowMatrix *Bt, double scalar);
00223 int ML_Epetra_CRSinsert(ML_Operator *, int, int *, double *, int);
00224 
00226 
00232 int ML_Operator2EpetraCrsMatrix(ML_Operator *Ke, Epetra_CrsMatrix * &
00233         CrsMatrix, int & MaxNumNonzeros,
00234         bool CheckNonzeroRow, double &, int base=0,bool verbose=false);
00235 
00236 inline int ML_Operator2EpetraCrsMatrix(ML_Operator *Ke, Epetra_CrsMatrix * &
00237         CrsMatrix)
00238 {
00239   int MaxNumNonzeros;
00240   double CPUTime;
00241   return(ML_Operator2EpetraCrsMatrix(Ke, CrsMatrix, MaxNumNonzeros, false, CPUTime));
00242 }
00243 
00244 Epetra_Map* Epetra_ML_readupdatevector(char* filename, Epetra_Comm& comm);
00245 Epetra_CrsMatrix* Epetra_ML_readaztecmatrix(char* filename,Epetra_Map& map,
00246                                             Epetra_Comm& comm);
00247 
00248 
00249 
00250 namespace ML_Epetra{
00251 
00253   int ML_Epetra_PtAP(const Epetra_CrsMatrix & A, const Epetra_CrsMatrix & P, Epetra_CrsMatrix *&Result,bool verbose=false);
00254   
00256   //treatment
00260   int* FindLocalDiricheltRowsFromOnesAndZeros(const Epetra_CrsMatrix & Matrix, int &numBCRows);
00261 
00263 
00268   void Apply_BCsToMatrixColumns(const int *dirichletRows, int numBCRows, const Epetra_CrsMatrix & Matrix);
00269   
00271   void Apply_BCsToMatrixRows(const int *dirichletRows, int numBCRows, const Epetra_CrsMatrix & Matrix);
00272 
00274 
00279   void Apply_BCsToMatrixColumns(const Epetra_RowMatrix & iBoundaryMatrix, const Epetra_RowMatrix & iMatrix);
00280   void Apply_BCsToMatrixColumns(const Epetra_IntVector &dirichletColumns,const Epetra_CrsMatrix & Matrix);
00281 
00282   
00284   void Apply_BCsToGradient( const Epetra_RowMatrix & EdgeMatrix,
00285                             const Epetra_RowMatrix & T);
00286 
00287 
00289   void Apply_OAZToMatrix(int *dirichletRows, int numBCRows, const Epetra_CrsMatrix & Matrix);
00290 
00292   Epetra_IntVector * LocalRowstoColumns(int *Rows, int numRows,const Epetra_CrsMatrix & Matrix);
00293 
00294 
00296   Epetra_IntVector * FindLocalDirichletColumnsFromRows(const int *dirichletRows, int numBCRows,const Epetra_CrsMatrix & Matrix);
00297 
00299   void Remove_Zeroed_Rows(const Epetra_CrsMatrix & Matrix, double tol=0.0);
00300 
00301 
00303 
00311 #ifdef HAVE_ML_EPETRAEXT
00312   Epetra_RowMatrix* ModifyEpetraMatrixColMap( const Epetra_RowMatrix &A,
00313                                    EpetraExt::CrsMatrix_SolverMap &transform,
00314                                    const char* matrixName=0, bool verbose=false);
00315 #endif
00316 
00317   
00318 }
00319 
00320 
00321 #ifdef FIXME
00322 string ML_toString(const int& x);
00323 string ML_toString(const double& x);
00324 #endif
00325 
00326 
00327 #ifdef __cplusplus
00328 extern "C" {
00329 #endif
00330 int ML_Operator_Destroy_DiscreteLaplacian();
00331 int ML_Operator_DiscreteLaplacian(ML_Operator* Op, int SymmetricPattern,
00332           double* x_coord, double* y_coord,
00333           double* z_coord, double theta,
00334           ML_Operator** NewOp);
00335 bool Epetra_ML_readaztecvector(char* filename, Epetra_MultiVector& Vector, 
00336                                Epetra_Map& map,Epetra_Comm& comm, int ivec);
00337 bool Epetra_ML_readvariableblocks(char* filename, Epetra_Map& map,
00338                                   Epetra_Comm& comm, 
00339                                   int**blocks, int** block_pde);
00340 bool Epetra_ML_writegidviz(char* filename, int label, 
00341                            Epetra_MultiVector& vector, int ivec, 
00342                            Epetra_Map& map, Epetra_Comm& comm);
00343                                   
00344 #ifdef __cplusplus
00345 }
00346 #endif
00347   
00349 
00361 void ML_BreakForDebugger(const Epetra_Comm& comm);
00362 
00363 #endif /* _ML_EPETRA_UTILS_H_ */