Actual source code: petscmat.h

  1: /* $Id: petscmat.h,v 1.219 2001/04/10 22:35:01 balay Exp $ */
  2: /*
  3:      Include file for the matrix component of PETSc
  4: */
  5: #ifndef __PETSCMAT_H
  7: #include "petscvec.h"

  9: #define MAT_COOKIE         PETSC_COOKIE+5

 11: /*S
 12:      Mat - Abstract PETSc matrix object

 14:    Level: beginner

 16:   Concepts: matrix; linear operator

 18: .seealso:  MatCreate(), MatType, MatSetType()
 19: S*/
 20: typedef struct _p_Mat*           Mat;

 22: /*E
 23:     MatType - String with the name of a PETSc matrix or the creation function
 24:        with an optional dynamic library name, for example
 25:        http://www.mcs.anl.gov/petsc/lib.a:mymatcreate()

 27:    Level: beginner

 29: .seealso: MatSetType(), Mat
 30: E*/
 31: #define MATSAME     "same"
 32: #define MATSEQMAIJ  "seqmaij"
 33: #define MATMPIMAIJ  "mpimaij"
 34: #define MATIS       "is"
 35: #define MATMPIROWBS "mpirowbs"
 36: #define MATSEQDENSE "seqdense"
 37: #define MATSEQAIJ   "seqaij"
 38: #define MATMPIAIJ   "mpiaij"
 39: #define MATSHELL    "shell"
 40: #define MATSEQBDIAG "seqbdiag"
 41: #define MATMPIBDIAG "mpibdiag"
 42: #define MATMPIDENSE "mpidense"
 43: #define MATSEQBAIJ  "seqbaij"
 44: #define MATMPIBAIJ  "mpibaij"
 45: #define MATMPIADJ   "mpiadj"
 46: #define MATSEQSBAIJ "seqsbaij"
 47: #define MATMPISBAIJ "mpisbaij"
 48: typedef char* MatType;

 50: EXTERN int MatCreate(MPI_Comm,int,int,int,int,Mat*);
 51: EXTERN int MatSetType(Mat,MatType);
 52: EXTERN int MatSetFromOptions(Mat);
 53: EXTERN int MatSetUpPreallocation(Mat);
 54: EXTERN int MatRegisterAll(char*);
 55: EXTERN int MatRegister(char*,char*,char*,int(*)(Mat));
 56: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
 57: #define MatRegisterDynamic(a,b,c,d) MatRegister(a,b,c,0)
 58: #else
 59: #define MatRegisterDynamic(a,b,c,d) MatRegister(a,b,c,d)
 60: #endif
 61: extern PetscTruth MatRegisterAllCalled;
 62: extern PetscFList MatList;

 64: EXTERN int MatCreate(MPI_Comm,int,int,int,int,Mat*);
 65: EXTERN int MatCreateSeqDense(MPI_Comm,int,int,Scalar*,Mat*);
 66: EXTERN int MatCreateMPIDense(MPI_Comm,int,int,int,int,Scalar*,Mat*);
 67: EXTERN int MatCreateSeqAIJ(MPI_Comm,int,int,int,int*,Mat*);
 68: EXTERN int MatCreateMPIAIJ(MPI_Comm,int,int,int,int,int,int*,int,int*,Mat*);
 69: EXTERN int MatCreateMPIRowbs(MPI_Comm,int,int,int,int*,Mat*);
 70: EXTERN int MatCreateSeqBDiag(MPI_Comm,int,int,int,int,int*,Scalar**,Mat*);
 71: EXTERN int MatCreateMPIBDiag(MPI_Comm,int,int,int,int,int,int*,Scalar**,Mat*);
 72: EXTERN int MatCreateSeqBAIJ(MPI_Comm,int,int,int,int,int*,Mat*);
 73: EXTERN int MatCreateMPIBAIJ(MPI_Comm,int,int,int,int,int,int,int*,int,int*,Mat*);
 74: EXTERN int MatCreateMPIAdj(MPI_Comm,int,int,int*,int*,int *,Mat*);
 75: EXTERN int MatCreateSeqSBAIJ(MPI_Comm,int,int,int,int,int*,Mat*);
 76: EXTERN int MatCreateMPISBAIJ(MPI_Comm,int,int,int,int,int,int,int*,int,int*,Mat*);
 77: EXTERN int MatCreateShell(MPI_Comm,int,int,int,int,void *,Mat*);
 78: EXTERN int MatDestroy(Mat);

 80: EXTERN int MatPrintHelp(Mat);
 81: EXTERN int MatGetMaps(Mat,Map*,Map*);

 83: /* ------------------------------------------------------------*/
 84: EXTERN int MatSetValues(Mat,int,int*,int,int*,Scalar*,InsertMode);
 85: EXTERN int MatSetValuesBlocked(Mat,int,int*,int,int*,Scalar*,InsertMode);

 87: /*S
 88:      MatStencil - Data structure (C struct) for storing information about a single row or
 89:         column of a matrix as index on an associated grid.

 91:    Level: beginner

 93:   Concepts: matrix; linear operator

 95: .seealso:  MatSetValuesStencil(), MatSetStencil()
 96: S*/
 97: typedef struct {
 98:   int k,j,i,c;
 99: } MatStencil;

101: EXTERN int MatSetValuesStencil(Mat,int,MatStencil*,int,MatStencil*,Scalar*,InsertMode);
102: EXTERN int MatSetValuesBlockedStencil(Mat,int,MatStencil*,int,MatStencil*,Scalar*,InsertMode);
103: EXTERN int MatSetStencil(Mat,int,int*,int*,int);

105: /*E
106:     MatAssemblyType - Indicates if the matrix is now to be used, or if you plan 
107:      to continue to add values to it

109:     Level: beginner

111: .seealso: MatAssemblyBegin(), MatAssemblyEnd()
112: E*/
113: typedef enum {MAT_FLUSH_ASSEMBLY=1,MAT_FINAL_ASSEMBLY=0} MatAssemblyType;
114: EXTERN int MatAssemblyBegin(Mat,MatAssemblyType);
115: EXTERN int MatAssemblyEnd(Mat,MatAssemblyType);
116: EXTERN int MatAssembled(Mat,PetscTruth*);

118: #define MatSetValue(v,i,j,va,mode) 
119: 0; {int _ierr,_row = i,_col = j; Scalar _va = va; 
120:   _MatSetValues(v,1,&_row,1,&_col,&_va,mode);CHKERRQ(_ierr); 
121: }
122: #define MatGetValue(v,i,j,va) 
123: 0; {int _ierr,_row = i,_col = j; 
124:   _MatGetValues(v,1,&_row,1,&_col,&va);CHKERRQ(_ierr); 
125: }
126: #define MatSetValueLocal(v,i,j,va,mode) 
127: 0; {int _ierr,_row = i,_col = j; Scalar _va = va; 
128:   _MatSetValuesLocal(v,1,&_row,1,&_col,&_va,mode);CHKERRQ(_ierr); 
129: }
130: /*E
131:     MatOption - Options that may be set for a matrix and its behavior or storage

133:     Level: beginner

135:    Any additions/changes here MUST also be made in include/finclude/petscmat.h

137: .seealso: MatSetOption()
138: E*/
139: typedef enum {MAT_ROW_ORIENTED=1,MAT_COLUMN_ORIENTED=2,MAT_ROWS_SORTED=4,
140:               MAT_COLUMNS_SORTED=8,MAT_NO_NEW_NONZERO_LOCATIONS=16,
141:               MAT_YES_NEW_NONZERO_LOCATIONS=32,MAT_SYMMETRIC=64,
142:               MAT_STRUCTURALLY_SYMMETRIC=65,MAT_NO_NEW_DIAGONALS=66,
143:               MAT_YES_NEW_DIAGONALS=67,MAT_INODE_LIMIT_1=68,MAT_INODE_LIMIT_2=69,
144:               MAT_INODE_LIMIT_3=70,MAT_INODE_LIMIT_4=71,MAT_INODE_LIMIT_5=72,
145:               MAT_IGNORE_OFF_PROC_ENTRIES=73,MAT_ROWS_UNSORTED=74,
146:               MAT_COLUMNS_UNSORTED=75,MAT_NEW_NONZERO_LOCATION_ERR=76,
147:               MAT_NEW_NONZERO_ALLOCATION_ERR=77,MAT_USE_HASH_TABLE=78,
148:               MAT_KEEP_ZEROED_ROWS=79,MAT_IGNORE_ZERO_ENTRIES=80,MAT_USE_INODES=81,
149:               MAT_DO_NOT_USE_INODES=82} MatOption;
150: EXTERN int MatSetOption(Mat,MatOption);
151: EXTERN int MatGetType(Mat,MatType*);

153: EXTERN int MatGetValues(Mat,int,int*,int,int*,Scalar*);
154: EXTERN int MatGetRow(Mat,int,int *,int **,Scalar**);
155: EXTERN int MatRestoreRow(Mat,int,int *,int **,Scalar**);
156: EXTERN int MatGetColumn(Mat,int,int *,int **,Scalar**);
157: EXTERN int MatRestoreColumn(Mat,int,int *,int **,Scalar**);
158: EXTERN int MatGetColumnVector(Mat,Vec,int);
159: EXTERN int MatGetArray(Mat,Scalar **);
160: EXTERN int MatRestoreArray(Mat,Scalar **);
161: EXTERN int MatGetBlockSize(Mat,int *);

163: EXTERN int MatMult(Mat,Vec,Vec);
164: EXTERN int MatMultAdd(Mat,Vec,Vec,Vec);
165: EXTERN int MatMultTranspose(Mat,Vec,Vec);
166: EXTERN int MatMultTransposeAdd(Mat,Vec,Vec,Vec);

168: /*E
169:     MatDuplicateOption - Indicates if a duplicated sparse matrix should have
170:   its numerical values copied over or just its nonzero structure.

172:     Level: beginner

174:    Any additions/changes here MUST also be made in include/finclude/petscmat.h

176: .seealso: MatDuplicate()
177: E*/
178: typedef enum {MAT_DO_NOT_COPY_VALUES,MAT_COPY_VALUES} MatDuplicateOption;

180: EXTERN int MatConvertRegister(char*,char*,char*,int (*)(Mat,MatType,Mat*));
181: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
182: #define MatConvertRegisterDynamic(a,b,c,d) MatConvertRegister(a,b,c,0)
183: #else
184: #define MatConvertRegisterDynamic(a,b,c,d) MatConvertRegister(a,b,c,d)
185: #endif
186: EXTERN int        MatConvertRegisterAll(char*);
187: EXTERN int        MatConvertRegisterDestroy(void);
188: extern PetscTruth MatConvertRegisterAllCalled;
189: extern PetscFList MatConvertList;
190: EXTERN int        MatConvert(Mat,MatType,Mat*);
191: EXTERN int        MatDuplicate(Mat,MatDuplicateOption,Mat*);

193: /*E
194:     MatStructure - Indicates if the matrix has the same nonzero structure

196:     Level: beginner

198:    Any additions/changes here MUST also be made in include/finclude/petscmat.h

200: .seealso: MatCopy(), SLESSetOperators(), PCSetOperators()
201: E*/
202: typedef enum {SAME_NONZERO_PATTERN,DIFFERENT_NONZERO_PATTERN,SAME_PRECONDITIONER} MatStructure;

204: EXTERN int MatCopy(Mat,Mat,MatStructure);
205: EXTERN int MatView(Mat,PetscViewer);

207: EXTERN int MatLoadRegister(char*,char*,char*,int (*)(PetscViewer,MatType,Mat*));
208: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
209: #define MatLoadRegisterDynamic(a,b,c,d) MatLoadRegister(a,b,c,0)
210: #else
211: #define MatLoadRegisterDynamic(a,b,c,d) MatLoadRegister(a,b,c,d)
212: #endif
213: EXTERN int        MatLoadRegisterAll(char*);
214: EXTERN int        MatLoadRegisterDestroy(void);
215: extern PetscTruth MatLoadRegisterAllCalled;
216: extern PetscFList MatLoadList;
217: EXTERN int        MatLoad(PetscViewer,MatType,Mat*);

219: EXTERN int MatGetRowIJ(Mat,int,PetscTruth,int*,int **,int **,PetscTruth *);
220: EXTERN int MatRestoreRowIJ(Mat,int,PetscTruth,int *,int **,int **,PetscTruth *);
221: EXTERN int MatGetColumnIJ(Mat,int,PetscTruth,int*,int **,int **,PetscTruth *);
222: EXTERN int MatRestoreColumnIJ(Mat,int,PetscTruth,int *,int **,int **,PetscTruth *);

224: /*S
225:      MatInfo - Context of matrix information, used with MatGetInfo()

227:    In Fortran this is simply a double precision array of dimension MAT_INFO_SIZE

229:    Level: intermediate

231:   Concepts: matrix^nonzero information

233: .seealso:  MatGetInfo(), MatInfoType
234: S*/
235: typedef struct {
236:   PetscLogDouble rows_global,columns_global;         /* number of global rows and columns */
237:   PetscLogDouble rows_local,columns_local;           /* number of local rows and columns */
238:   PetscLogDouble block_size;                         /* block size */
239:   PetscLogDouble nz_allocated,nz_used,nz_unneeded;   /* number of nonzeros */
240:   PetscLogDouble memory;                             /* memory allocated */
241:   PetscLogDouble assemblies;                         /* number of matrix assemblies called */
242:   PetscLogDouble mallocs;                            /* number of mallocs during MatSetValues() */
243:   PetscLogDouble fill_ratio_given,fill_ratio_needed; /* fill ratio for LU/ILU */
244:   PetscLogDouble factor_mallocs;                     /* number of mallocs during factorization */
245: } MatInfo;

247: /*E
248:     MatInfoType - Indicates if you want information about the local part of the matrix,
249:      the entire parallel matrix or the maximum over all the local parts.

251:     Level: beginner

253:    Any additions/changes here MUST also be made in include/finclude/petscmat.h

255: .seealso: MatGetInfo(), MatInfo
256: E*/
257: typedef enum {MAT_LOCAL=1,MAT_GLOBAL_MAX=2,MAT_GLOBAL_SUM=3} MatInfoType;
258: EXTERN int MatGetInfo(Mat,MatInfoType,MatInfo*);
259: EXTERN int MatValid(Mat,PetscTruth*);
260: EXTERN int MatGetDiagonal(Mat,Vec);
261: EXTERN int MatGetRowMax(Mat,Vec);
262: EXTERN int MatTranspose(Mat,Mat*);
263: EXTERN int MatPermute(Mat,IS,IS,Mat *);
264: EXTERN int MatDiagonalScale(Mat,Vec,Vec);
265: EXTERN int MatDiagonalSet(Mat,Vec,InsertMode);
266: EXTERN int MatEqual(Mat,Mat,PetscTruth*);

268: EXTERN int MatNorm(Mat,NormType,double *);
269: EXTERN int MatZeroEntries(Mat);
270: EXTERN int MatZeroRows(Mat,IS,Scalar*);
271: EXTERN int MatZeroColumns(Mat,IS,Scalar*);

273: EXTERN int MatUseScaledForm(Mat,PetscTruth);
274: EXTERN int MatScaleSystem(Mat,Vec,Vec);
275: EXTERN int MatUnScaleSystem(Mat,Vec,Vec);

277: EXTERN int MatGetSize(Mat,int*,int*);
278: EXTERN int MatGetLocalSize(Mat,int*,int*);
279: EXTERN int MatGetOwnershipRange(Mat,int*,int*);

281: /*E
282:     MatReuse - Indicates if matrices obtained from a previous call to MatGetSubMatrices()
283:      or MatGetSubMatrix() are to be reused to store the new matrix values.

285:     Level: beginner

287:    Any additions/changes here MUST also be made in include/finclude/petscmat.h

289: .seealso: MatGetSubMatrices(), MatGetSubMatrix(), MatDestroyMatrices()
290: E*/
291: typedef enum {MAT_INITIAL_MATRIX,MAT_REUSE_MATRIX} MatReuse;
292: EXTERN int MatGetSubMatrices(Mat,int,IS *,IS *,MatReuse,Mat **);
293: EXTERN int MatDestroyMatrices(int,Mat **);
294: EXTERN int MatGetSubMatrix(Mat,IS,IS,int,MatReuse,Mat *);

296: EXTERN int MatIncreaseOverlap(Mat,int,IS *,int);

298: EXTERN int MatAXPY(Scalar *,Mat,Mat);
299: EXTERN int MatAYPX(Scalar *,Mat,Mat);
300: EXTERN int MatCompress(Mat);

302: EXTERN int MatScale(Scalar *,Mat);
303: EXTERN int MatShift(Scalar *,Mat);

305: EXTERN int MatSetLocalToGlobalMapping(Mat,ISLocalToGlobalMapping);
306: EXTERN int MatSetLocalToGlobalMappingBlock(Mat,ISLocalToGlobalMapping);
307: EXTERN int MatZeroRowsLocal(Mat,IS,Scalar*);
308: EXTERN int MatSetValuesLocal(Mat,int,int*,int,int*,Scalar*,InsertMode);
309: EXTERN int MatSetValuesBlockedLocal(Mat,int,int*,int,int*,Scalar*,InsertMode);

311: EXTERN int MatSetStashInitialSize(Mat,int,int);

313: EXTERN int MatInterpolateAdd(Mat,Vec,Vec,Vec);
314: EXTERN int MatInterpolate(Mat,Vec,Vec);
315: EXTERN int MatRestrict(Mat,Vec,Vec);

317: /*
318:       These three (or four) macros MUST be used together. The third one closes the open { of the first one
319: */
320: #define MatPreallocateInitialize(comm,nrows,ncols,dnz,onz) 0; 
321: { 
322:   int __ierr,__tmp = (nrows),__ctmp = (ncols),__rstart,__start,__end; 
323:   __PetscMalloc(2*__tmp*sizeof(int),&dnz);CHKERRQ(__ierr);onz = dnz + __tmp;
324:   __PetscMemzero(dnz,2*__tmp*sizeof(int));CHKERRQ(__ierr);
325:   __MPI_Scan(&__ctmp,&__end,1,MPI_INT,MPI_SUM,comm);CHKERRQ(__ierr); __start = __end - __ctmp;
326:   __MPI_Scan(&__tmp,&__rstart,1,MPI_INT,MPI_SUM,comm);CHKERRQ(__ierr); __rstart = __rstart - __tmp;

328: #define MatPreallocateSetLocal(map,nrows,rows,ncols,cols,dnz,onz) 0;
329: {
330:   int __l;
331:   __ISLocalToGlobalMappingApply(map,nrows,rows,rows);CHKERRQ(__ierr);
332:   __ISLocalToGlobalMappingApply(map,ncols,cols,cols);CHKERRQ(__ierr);
333:   for (__l=0;__l<nrows;__l++) {
334:     __MatPreallocateSet((rows)[__l],ncols,cols,dnz,onz);CHKERRQ(__ierr);
335:   }
336: }
337: 
338: #define MatPreallocateSet(row,nc,cols,dnz,onz) 0;
339: { int __i; 
340:   for (__i=0; __i<nc; __i++) {
341:     if (cols[__i] < __start || cols[__i] >= __end) onz[row - __rstart]++; 
342:   }
343:   dnz[row - __rstart] = nc - onz[row - __rstart];
344: }

346: #define MatPreallocateFinalize(dnz,onz) 0;__PetscFree(dnz);CHKERRQ(__ierr);}

348: /* Routines unique to particular data structures */
349: EXTERN int MatShellGetContext(Mat,void **);

351: EXTERN int MatBDiagGetData(Mat,int*,int*,int**,int**,Scalar***);
352: EXTERN int MatSeqAIJSetColumnIndices(Mat,int *);
353: EXTERN int MatSeqBAIJSetColumnIndices(Mat,int *);
354: EXTERN int MatCreateSeqAIJWithArrays(MPI_Comm,int,int,int*,int*,Scalar *,Mat*);

356: EXTERN int MatSeqBAIJSetPreallocation(Mat,int,int,int*);
357: EXTERN int MatSeqSBAIJSetPreallocation(Mat,int,int,int*);
358: EXTERN int MatSeqAIJSetPreallocation(Mat,int,int*);
359: EXTERN int MatSeqDensePreallocation(Mat,Scalar*);
360: EXTERN int MatSeqBDiagSetPreallocation(Mat,int,int,int*,Scalar**);
361: EXTERN int MatSeqDenseSetPreallocation(Mat,Scalar*);

363: EXTERN int MatMPIBAIJSetPreallocation(Mat,int,int,int*,int,int*);
364: EXTERN int MatMPISBAIJSetPreallocation(Mat,int,int,int*,int,int*);
365: EXTERN int MatMPIAIJSetPreallocation(Mat,int,int*,int,int*);
366: EXTERN int MatMPIDensePreallocation(Mat,Scalar*);
367: EXTERN int MatMPIBDiagSetPreallocation(Mat,int,int,int*,Scalar**);
368: EXTERN int MatMPIAdjSetPreallocation(Mat,int*,int*,int*);
369: EXTERN int MatMPIDenseSetPreallocation(Mat,Scalar*);
370: EXTERN int MatMPIRowbsSetPreallocation(Mat,int,int*);
371: EXTERN int MatMPIAIJGetSeqAIJ(Mat,Mat*,Mat*,int**);
372: EXTERN int MatMPIBAIJGetSeqBAIJ(Mat,Mat*,Mat*,int**);

374: EXTERN int MatStoreValues(Mat);
375: EXTERN int MatRetrieveValues(Mat);

377: /* 
378:   These routines are not usually accessed directly, rather solving is 
379:   done through the SLES, KSP and PC interfaces.
380: */

382: /*E
383:     MatOrderingType - String with the name of a PETSc matrix ordering or the creation function
384:        with an optional dynamic library name, for example 
385:        http://www.mcs.anl.gov/petsc/lib.a:orderingcreate()

387:    Level: beginner

389: .seealso: MatGetOrdering()
390: E*/
391: typedef char* MatOrderingType;
392: #define MATORDERING_NATURAL   "natural"
393: #define MATORDERING_ND        "nd"
394: #define MATORDERING_1WD       "1wd"
395: #define MATORDERING_RCM       "rcm"
396: #define MATORDERING_QMD       "qmd"
397: #define MATORDERING_ROWLENGTH "rowlength"
398: #define MATORDERING_DSC_ND    "dsc_nd"
399: #define MATORDERING_DSC_MMD   "dsc_mmd"
400: #define MATORDERING_DSC_MDF   "dsc_mdf"

402: EXTERN int MatGetOrdering(Mat,MatOrderingType,IS*,IS*);
403: EXTERN int MatOrderingRegister(char*,char*,char*,int(*)(Mat,MatOrderingType,IS*,IS*));
404: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
405: #define MatOrderingRegisterDynamic(a,b,c,d) MatOrderingRegister(a,b,c,0)
406: #else
407: #define MatOrderingRegisterDynamic(a,b,c,d) MatOrderingRegister(a,b,c,d)
408: #endif
409: EXTERN int        MatOrderingRegisterDestroy(void);
410: EXTERN int        MatOrderingRegisterAll(char*);
411: extern PetscTruth MatOrderingRegisterAllCalled;
412: extern PetscFList      MatOrderingList;

414: EXTERN int MatReorderForNonzeroDiagonal(Mat,double,IS,IS);

416: EXTERN int MatCholeskyFactor(Mat,IS,double);
417: EXTERN int MatCholeskyFactorSymbolic(Mat,IS,double,Mat*);
418: EXTERN int MatCholeskyFactorNumeric(Mat,Mat*);

420: /*S 
421:    MatILUInfo - Data based into the matrix ILU factorization routines

423:    In Fortran these are simply double precision arrays of size MAT_ILUINFO_SIZE

425:    Notes: These are not usually directly used by users, instead use the PC type of ILU
426:           All entries are double precision.

428:    Level: developer

430: .seealso: MatILUFactorSymbolic(), MatILUFactor(), MatLUInfo, MatCholeskyInfo

432: S*/
433: typedef struct {
434:   double     levels;         /* ILU(levels) */
435:   double     fill;           /* expected fill; nonzeros in factored matrix/nonzeros in original matrix*/
436:   double     diagonal_fill;  /* force diagonal to fill in if initially not filled */
437:   double     dt;             /* drop tolerance */
438:   double     dtcol;          /* tolerance for pivoting */
439:   double     dtcount;        /* maximum nonzeros to be allowed per row */
440:   double     damping;        /* scaling of identity added to matrix to prevent zero pivots */
441:   double     damp;           /* if is 1.0 and factorization fails, damp until successful */
442: } MatILUInfo;

444: /*S 
445:    MatLUInfo - Data based into the matrix LU factorization routines

447:    In Fortran these are simply double precision arrays of size MAT_LUINFO_SIZE

449:    Notes: These are not usually directly used by users, instead use the PC type of LU
450:           All entries are double precision.

452:    Level: developer

454: .seealso: MatLUFactorSymbolic(), MatILUInfo, MatCholeskyInfo

456: S*/
457: typedef struct {
458:   double     fill;    /* expected fill; nonzeros in factored matrix/nonzeros in original matrix */
459:   double     dtcol;   /* tolerance for pivoting; pivot if off_diagonal*dtcol > diagonal */
460:   double     damping; /* scaling of identity added to matrix to prevent zero pivots */
461:   double     damp;    /* if this is 1.0 and factorization fails, damp until successful */
462: } MatLUInfo;

464: /*S 
465:    MatCholeskyInfo - Data based into the matrix Cholesky factorization routines

467:    In Fortran these are simply double precision arrays of size MAT_CHOLESKYINFO_SIZE

469:    Notes: These are not usually directly used by users, instead use the PC type of Cholesky
470:           All entries are double precision.

472:    Level: developer

474: .seealso: MatCholeskyFactorSymbolic(), MatLUInfo, MatILUInfo

476: S*/
477: typedef struct {
478:   double     fill;    /* expected fill; nonzeros in factored matrix/nonzeros in original matrix */
479:   double     damping; /* scaling of identity added to matrix to prevent zero pivots */
480:   double     damp;    /* if this is 1.0 and factorization fails, damp until successful */
481: } MatCholeskyInfo;

483: EXTERN int MatLUFactor(Mat,IS,IS,MatLUInfo*);
484: EXTERN int MatILUFactor(Mat,IS,IS,MatILUInfo*);
485: EXTERN int MatLUFactorSymbolic(Mat,IS,IS,MatLUInfo*,Mat*);
486: EXTERN int MatILUFactorSymbolic(Mat,IS,IS,MatILUInfo*,Mat*);
487: EXTERN int MatICCFactorSymbolic(Mat,IS,double,int,Mat*);
488: EXTERN int MatICCFactor(Mat,IS,double,int);
489: EXTERN int MatLUFactorNumeric(Mat,Mat*);
490: EXTERN int MatILUDTFactor(Mat,MatILUInfo*,IS,IS,Mat *);

492: EXTERN int MatSolve(Mat,Vec,Vec);
493: EXTERN int MatForwardSolve(Mat,Vec,Vec);
494: EXTERN int MatBackwardSolve(Mat,Vec,Vec);
495: EXTERN int MatSolveAdd(Mat,Vec,Vec,Vec);
496: EXTERN int MatSolveTranspose(Mat,Vec,Vec);
497: EXTERN int MatSolveTransposeAdd(Mat,Vec,Vec,Vec);

499: EXTERN int MatSetUnfactored(Mat);

501: /*  MatSORType may be bitwise ORd together, so do not change the numbers */
502: /*E
503:     MatSORType - What type of (S)SOR to perform

505:     Level: beginner

507:    May be bitwise ORd together

509:    Any additions/changes here MUST also be made in include/finclude/petscmat.h

511: .seealso: MatRelax()
512: E*/
513: typedef enum {SOR_FORWARD_SWEEP=1,SOR_BACKWARD_SWEEP=2,SOR_SYMMETRIC_SWEEP=3,
514:               SOR_LOCAL_FORWARD_SWEEP=4,SOR_LOCAL_BACKWARD_SWEEP=8,
515:               SOR_LOCAL_SYMMETRIC_SWEEP=12,SOR_ZERO_INITIAL_GUESS=16,
516:               SOR_EISENSTAT=32,SOR_APPLY_UPPER=64,SOR_APPLY_LOWER=128} MatSORType;
517: EXTERN int MatRelax(Mat,Vec,double,MatSORType,double,int,Vec);

519: /* 
520:     These routines are for efficiently computing Jacobians via finite differences.
521: */

523: /*E
524:     MatColoringType - String with the name of a PETSc matrix coloring or the creation function
525:        with an optional dynamic library name, for example 
526:        http://www.mcs.anl.gov/petsc/lib.a:coloringcreate()

528:    Level: beginner

530: .seealso: MatGetColoring()
531: E*/
532: typedef char* MatColoringType;
533: #define MATCOLORING_NATURAL "natural"
534: #define MATCOLORING_SL      "sl"
535: #define MATCOLORING_LF      "lf"
536: #define MATCOLORING_ID      "id"

538: EXTERN int MatGetColoring(Mat,MatColoringType,ISColoring*);
539: EXTERN int MatColoringRegister(char*,char*,char*,int(*)(Mat,MatColoringType,ISColoring *));
540: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
541: #define MatColoringRegisterDynamic(a,b,c,d) MatColoringRegister(a,b,c,0)
542: #else
543: #define MatColoringRegisterDynamic(a,b,c,d) MatColoringRegister(a,b,c,d)
544: #endif
545: EXTERN int        MatColoringRegisterAll(char *);
546: extern PetscTruth MatColoringRegisterAllCalled;
547: EXTERN int        MatColoringRegisterDestroy(void);
548: EXTERN int        MatColoringPatch(Mat,int,int,int *,ISColoring*);

550: #define MAT_FDCOLORING_COOKIE PETSC_COOKIE + 23
551: /*S
552:      MatFDColoring - Object for computing a sparse Jacobian via finite differences
553:         and coloring

555:    Level: beginner

557:   Concepts: coloring, sparse Jacobian, finite differences

559: .seealso:  MatFDColoringCreate()
560: S*/
561: typedef struct _p_MatFDColoring *MatFDColoring;

563: EXTERN int MatFDColoringCreate(Mat,ISColoring,MatFDColoring *);
564: EXTERN int MatFDColoringDestroy(MatFDColoring);
565: EXTERN int MatFDColoringView(MatFDColoring,PetscViewer);
566: EXTERN int MatFDColoringSetFunction(MatFDColoring,int (*)(void),void*);
567: EXTERN int MatFDColoringSetParameters(MatFDColoring,double,double);
568: EXTERN int MatFDColoringSetFrequency(MatFDColoring,int);
569: EXTERN int MatFDColoringGetFrequency(MatFDColoring,int*);
570: EXTERN int MatFDColoringSetFromOptions(MatFDColoring);
571: EXTERN int MatFDColoringApply(Mat,MatFDColoring,Vec,MatStructure*,void *);
572: EXTERN int MatFDColoringApplyTS(Mat,MatFDColoring,double,Vec,MatStructure*,void *);
573: EXTERN int MatFDColoringSetRecompute(MatFDColoring);

575: /* 
576:     These routines are for partitioning matrices: currently used only 
577:   for adjacency matrix, MatCreateMPIAdj().
578: */
579: #define MATPARTITIONING_COOKIE PETSC_COOKIE + 25

581: /*S
582:      MatPartitioning - Object for managing the partitioning of a matrix or graph

584:    Level: beginner

586:   Concepts: partitioning

588: .seealso:  MatParitioningCreate(), MatPartitioningType
589: S*/
590: typedef struct _p_MatPartitioning *MatPartitioning;

592: /*E
593:     MatPartitioningType - String with the name of a PETSc matrix partitioing or the creation function
594:        with an optional dynamic library name, for example 
595:        http://www.mcs.anl.gov/petsc/lib.a:partitioningcreate()

597:    Level: beginner

599: .seealso: MatPartitioingCreate(), MatPartitioning
600: E*/
601: typedef char* MatPartitioningType;
602: #define MATPARTITIONING_CURRENT  "current"
603: #define MATPARTITIONING_PARMETIS "parmetis"

605: EXTERN int MatPartitioningCreate(MPI_Comm,MatPartitioning*);
606: EXTERN int MatPartitioningSetType(MatPartitioning,MatPartitioningType);
607: EXTERN int MatPartitioningSetAdjacency(MatPartitioning,Mat);
608: EXTERN int MatPartitioningSetVertexWeights(MatPartitioning,int*);
609: EXTERN int MatPartitioningApply(MatPartitioning,IS*);
610: EXTERN int MatPartitioningDestroy(MatPartitioning);

612: EXTERN int MatPartitioningRegister(char*,char*,char*,int(*)(MatPartitioning));
613: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
614: #define MatPartitioningRegisterDynamic(a,b,c,d) MatPartitioningRegister(a,b,c,0)
615: #else
616: #define MatPartitioningRegisterDynamic(a,b,c,d) MatPartitioningRegister(a,b,c,d)
617: #endif

619: EXTERN int        MatPartitioningRegisterAll(char *);
620: extern PetscTruth MatPartitioningRegisterAllCalled;
621: EXTERN int        MatPartitioningRegisterDestroy(void);

623: EXTERN int MatPartitioningView(MatPartitioning,PetscViewer);
624: EXTERN int MatPartitioningSetFromOptions(MatPartitioning);
625: EXTERN int MatPartitioningGetType(MatPartitioning,MatPartitioningType*);

627: EXTERN int MatPartitioningParmetisSetCoarseSequential(MatPartitioning);

629: /*
630:     If you add entries here you must also add them to finclude/petscmat.h
631: */
632: typedef enum { MATOP_SET_VALUES=0,
633:                MATOP_GET_ROW=1,
634:                MATOP_RESTORE_ROW=2,
635:                MATOP_MULT=3,
636:                MATOP_MULT_ADD=4,
637:                MATOP_MULT_TRANSPOSE=5,
638:                MATOP_MULT_TRANSPOSE_ADD=6,
639:                MATOP_SOLVE=7,
640:                MATOP_SOLVE_ADD=8,
641:                MATOP_SOLVE_TRANSPOSE=9,
642:                MATOP_SOLVE_TRANSPOSE_ADD=10,
643:                MATOP_LUFACTOR=11,
644:                MATOP_CHOLESKYFACTOR=12,
645:                MATOP_RELAX=13,
646:                MATOP_TRANSPOSE=14,
647:                MATOP_GETINFO=15,
648:                MATOP_EQUAL=16,
649:                MATOP_GET_DIAGONAL=17,
650:                MATOP_DIAGONAL_SCALE=18,
651:                MATOP_NORM=19,
652:                MATOP_ASSEMBLY_BEGIN=20,
653:                MATOP_ASSEMBLY_END=21,
654:                MATOP_COMPRESS=22,
655:                MATOP_SET_OPTION=23,
656:                MATOP_ZERO_ENTRIES=24,
657:                MATOP_ZERO_ROWS=25,
658:                MATOP_LUFACTOR_SYMBOLIC=26,
659:                MATOP_LUFACTOR_NUMERIC=27,
660:                MATOP_CHOLESKY_FACTOR_SYMBOLIC=28,
661:                MATOP_CHOLESKY_FACTOR_NUMERIC=29,
662:                MATOP_GET_SIZE=30,
663:                MATOP_GET_LOCAL_SIZE=31,
664:                MATOP_GET_OWNERSHIP_RANGE=32,
665:                MATOP_ILUFACTOR_SYMBOLIC=33,
666:                MATOP_ICCFACTOR_SYMBOLIC=34,
667:                MATOP_GET_ARRAY=35,
668:                MATOP_RESTORE_ARRAY=36,

670:                MATOP_CONVERT_SAME_TYPE=37,
671:                MATOP_FORWARD_SOLVE=38,
672:                MATOP_BACKWARD_SOLVE=39,
673:                MATOP_ILUFACTOR=40,
674:                MATOP_ICCFACTOR=41,
675:                MATOP_AXPY=42,
676:                MATOP_GET_SUBMATRICES=43,
677:                MATOP_INCREASE_OVERLAP=44,
678:                MATOP_GET_VALUES=45,
679:                MATOP_COPY=46,
680:                MATOP_PRINT_HELP=47,
681:                MATOP_SCALE=48,
682:                MATOP_SHIFT=49,
683:                MATOP_DIAGONAL_SHIFT=50,
684:                MATOP_ILUDT_FACTOR=51,
685:                MATOP_GET_BLOCK_SIZE=52,
686:                MATOP_GET_ROW_IJ=53,
687:                MATOP_RESTORE_ROW_IJ=54,
688:                MATOP_GET_COLUMN_IJ=55,
689:                MATOP_RESTORE_COLUMN_IJ=56,
690:                MATOP_FDCOLORING_CREATE=57,
691:                MATOP_COLORING_PATCH=58,
692:                MATOP_SET_UNFACTORED=59,
693:                MATOP_PERMUTE=60,
694:                MATOP_SET_VALUES_BLOCKED=61,
695:                MATOP_DESTROY=250,
696:                MATOP_VIEW=251
697:              } MatOperation;
698: EXTERN int MatHasOperation(Mat,MatOperation,PetscTruth*);
699: EXTERN int MatShellSetOperation(Mat,MatOperation,void(*)());
700: EXTERN int MatShellGetOperation(Mat,MatOperation,void(**)());
701: EXTERN int MatShellSetContext(Mat,void*);

703: /*
704:    Codes for matrices stored on disk. By default they are
705:  stored in a universal format. By changing the format with 
706:  PetscViewerSetFormat(viewer,PETSC_VIEWER_BINARY_NATIVE); the matrices will
707:  be stored in a way natural for the matrix, for example dense matrices
708:  would be stored as dense. Matrices stored this way may only be
709:  read into matrices of the same time.
710: */
711: #define MATRIX_BINARY_FORMAT_DENSE -1

713: /*
714:      New matrix classes not yet distributed 
715: */
716: /*
717:     MatAIJIndices is a data structure for storing the nonzero location information
718:   for sparse matrices. Several matrices with identical nonzero structure can share
719:   the same MatAIJIndices.
720: */
721: typedef struct _p_MatAIJIndices* MatAIJIndices;

723: EXTERN int MatCreateAIJIndices(int,int,int*,int*,PetscTruth,MatAIJIndices*);
724: EXTERN int MatCreateAIJIndicesEmpty(int,int,int*,PetscTruth,MatAIJIndices*);
725: EXTERN int MatAttachAIJIndices(MatAIJIndices,MatAIJIndices*);
726: EXTERN int MatDestroyAIJIndices(MatAIJIndices);
727: EXTERN int MatCopyAIJIndices(MatAIJIndices,MatAIJIndices*);
728: EXTERN int MatValidateAIJIndices(int,MatAIJIndices);
729: EXTERN int MatShiftAIJIndices(MatAIJIndices);
730: EXTERN int MatShrinkAIJIndices(MatAIJIndices);
731: EXTERN int MatTransposeAIJIndices(MatAIJIndices,MatAIJIndices*);

733: EXTERN int MatCreateSeqCSN(MPI_Comm,int,int,int*,int,Mat*);
734: EXTERN int MatCreateSeqCSN_Single(MPI_Comm,int,int,int*,int,Mat*);
735: EXTERN int MatCreateSeqCSNWithPrecision(MPI_Comm,int,int,int*,int,ScalarPrecision,Mat*);

737: EXTERN int MatCreateSeqCSNIndices(MPI_Comm,MatAIJIndices,int,Mat *);
738: EXTERN int MatCreateSeqCSNIndices_Single(MPI_Comm,MatAIJIndices,int,Mat *);
739: EXTERN int MatCreateSeqCSNIndicesWithPrecision(MPI_Comm,MatAIJIndices,int,ScalarPrecision,Mat *);

741: EXTERN int MatMPIBAIJSetHashTableFactor(Mat,PetscReal);
742: EXTERN int MatSeqAIJGetInodeSizes(Mat,int *,int *[],int *);
743: EXTERN int MatMPIRowbsGetColor(Mat,ISColoring *);

745: /*S
746:      MatNullSpace - Object that removes a null space from a vector, i.e.
747:          orthogonalizes the vector to a subsapce

749:    Level: beginner

751:   Concepts: matrix; linear operator, null space

753: .seealso:  MatNullSpaceCreate()
754: S*/
755: typedef struct _p_MatNullSpace* MatNullSpace;

757: #define MATNULLSPACE_COOKIE    PETSC_COOKIE+17

759: EXTERN int MatNullSpaceCreate(MPI_Comm,int,int,Vec *,MatNullSpace*);
760: EXTERN int MatNullSpaceDestroy(MatNullSpace);
761: EXTERN int MatNullSpaceRemove(MatNullSpace,Vec,Vec*);
762: EXTERN int MatNullSpaceAttach(Mat,MatNullSpace);
763: EXTERN int MatNullSpaceTest(MatNullSpace,Mat);

765: EXTERN int MatReorderingSeqSBAIJ(Mat A,IS isp);
766: EXTERN int MatMPISBAIJSetHashTableFactor(Mat,PetscReal);
767: EXTERN int MatSeqSBAIJSetColumnIndices(Mat,int *);


770: EXTERN int MatCreateMAIJ(Mat,int,Mat*);
771: EXTERN int MatMAIJRedimension(Mat,int,Mat*);
772: EXTERN int MatMAIJGetAIJ(Mat,Mat*);

774: EXTERN int MatMPIAdjSetValues(Mat,int*,int*,int*);

776: EXTERN int MatComputeExplicitOperator(Mat,Mat*);

778: #endif