Actual source code: petsc-matimpl.h
petsc-3.3-p7 2013-05-11
2: #ifndef __MATIMPL_H
5: #include <petscmat.h>
7: /*
8: This file defines the parts of the matrix data structure that are
9: shared by all matrix types.
10: */
12: /*
13: If you add entries here also add them to the MATOP enum
14: in include/petscmat.h and include/finclude/petscmat.h
15: */
16: typedef struct _MatOps *MatOps;
17: struct _MatOps {
18: /* 0*/
19: PetscErrorCode (*setvalues)(Mat,PetscInt,const PetscInt[],PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
20: PetscErrorCode (*getrow)(Mat,PetscInt,PetscInt *,PetscInt*[],PetscScalar*[]);
21: PetscErrorCode (*restorerow)(Mat,PetscInt,PetscInt *,PetscInt *[],PetscScalar *[]);
22: PetscErrorCode (*mult)(Mat,Vec,Vec);
23: PetscErrorCode (*multadd)(Mat,Vec,Vec,Vec);
24: /* 5*/
25: PetscErrorCode (*multtranspose)(Mat,Vec,Vec);
26: PetscErrorCode (*multtransposeadd)(Mat,Vec,Vec,Vec);
27: PetscErrorCode (*solve)(Mat,Vec,Vec);
28: PetscErrorCode (*solveadd)(Mat,Vec,Vec,Vec);
29: PetscErrorCode (*solvetranspose)(Mat,Vec,Vec);
30: /*10*/
31: PetscErrorCode (*solvetransposeadd)(Mat,Vec,Vec,Vec);
32: PetscErrorCode (*lufactor)(Mat,IS,IS,const MatFactorInfo*);
33: PetscErrorCode (*choleskyfactor)(Mat,IS,const MatFactorInfo*);
34: PetscErrorCode (*sor)(Mat,Vec,PetscReal,MatSORType,PetscReal,PetscInt,PetscInt,Vec);
35: PetscErrorCode (*transpose)(Mat,MatReuse,Mat *);
36: /*15*/
37: PetscErrorCode (*getinfo)(Mat,MatInfoType,MatInfo*);
38: PetscErrorCode (*equal)(Mat,Mat,PetscBool *);
39: PetscErrorCode (*getdiagonal)(Mat,Vec);
40: PetscErrorCode (*diagonalscale)(Mat,Vec,Vec);
41: PetscErrorCode (*norm)(Mat,NormType,PetscReal*);
42: /*20*/
43: PetscErrorCode (*assemblybegin)(Mat,MatAssemblyType);
44: PetscErrorCode (*assemblyend)(Mat,MatAssemblyType);
45: PetscErrorCode (*setoption)(Mat,MatOption,PetscBool );
46: PetscErrorCode (*zeroentries)(Mat);
47: /*24*/
48: PetscErrorCode (*zerorows)(Mat,PetscInt,const PetscInt[],PetscScalar,Vec,Vec);
49: PetscErrorCode (*lufactorsymbolic)(Mat,Mat,IS,IS,const MatFactorInfo*);
50: PetscErrorCode (*lufactornumeric)(Mat,Mat,const MatFactorInfo*);
51: PetscErrorCode (*choleskyfactorsymbolic)(Mat,Mat,IS,const MatFactorInfo*);
52: PetscErrorCode (*choleskyfactornumeric)(Mat,Mat,const MatFactorInfo*);
53: /*29*/
54: PetscErrorCode (*setup)(Mat);
55: PetscErrorCode (*ilufactorsymbolic)(Mat,Mat,IS,IS,const MatFactorInfo*);
56: PetscErrorCode (*iccfactorsymbolic)(Mat,Mat,IS,const MatFactorInfo*);
57: PetscErrorCode (*getarray)(Mat,PetscScalar**);
58: PetscErrorCode (*restorearray)(Mat,PetscScalar**);
59: /*34*/
60: PetscErrorCode (*duplicate)(Mat,MatDuplicateOption,Mat*);
61: PetscErrorCode (*forwardsolve)(Mat,Vec,Vec);
62: PetscErrorCode (*backwardsolve)(Mat,Vec,Vec);
63: PetscErrorCode (*ilufactor)(Mat,IS,IS,const MatFactorInfo*);
64: PetscErrorCode (*iccfactor)(Mat,IS,const MatFactorInfo*);
65: /*39*/
66: PetscErrorCode (*axpy)(Mat,PetscScalar,Mat,MatStructure);
67: PetscErrorCode (*getsubmatrices)(Mat,PetscInt,const IS[],const IS[],MatReuse,Mat *[]);
68: PetscErrorCode (*increaseoverlap)(Mat,PetscInt,IS[],PetscInt);
69: PetscErrorCode (*getvalues)(Mat,PetscInt,const PetscInt[],PetscInt,const PetscInt[],PetscScalar []);
70: PetscErrorCode (*copy)(Mat,Mat,MatStructure);
71: /*44*/
72: PetscErrorCode (*getrowmax)(Mat,Vec,PetscInt[]);
73: PetscErrorCode (*scale)(Mat,PetscScalar);
74: PetscErrorCode (*shift)(Mat,PetscScalar);
75: PetscErrorCode (*diagonalset)(Mat,Vec,InsertMode);
76: PetscErrorCode (*zerorowscolumns)(Mat,PetscInt,const PetscInt[],PetscScalar,Vec,Vec);
77: /*49*/
78: PetscErrorCode (*dummy6)(void);
79: PetscErrorCode (*getrowij)(Mat,PetscInt,PetscBool ,PetscBool ,PetscInt*,PetscInt *[],PetscInt *[],PetscBool *);
80: PetscErrorCode (*restorerowij)(Mat,PetscInt,PetscBool ,PetscBool ,PetscInt *,PetscInt *[],PetscInt *[],PetscBool *);
81: PetscErrorCode (*getcolumnij)(Mat,PetscInt,PetscBool ,PetscBool ,PetscInt*,PetscInt *[],PetscInt *[],PetscBool *);
82: PetscErrorCode (*restorecolumnij)(Mat,PetscInt,PetscBool ,PetscBool ,PetscInt*,PetscInt *[],PetscInt *[],PetscBool *);
83: /*54*/
84: PetscErrorCode (*fdcoloringcreate)(Mat,ISColoring,MatFDColoring);
85: PetscErrorCode (*coloringpatch)(Mat,PetscInt,PetscInt,ISColoringValue[],ISColoring*);
86: PetscErrorCode (*setunfactored)(Mat);
87: PetscErrorCode (*permute)(Mat,IS,IS,Mat*);
88: PetscErrorCode (*setvaluesblocked)(Mat,PetscInt,const PetscInt[],PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
89: /*59*/
90: PetscErrorCode (*getsubmatrix)(Mat,IS,IS,MatReuse,Mat*);
91: PetscErrorCode (*destroy)(Mat);
92: PetscErrorCode (*view)(Mat,PetscViewer);
93: PetscErrorCode (*convertfrom)(Mat, const MatType,MatReuse,Mat*);
94: PetscErrorCode (*dummy10)(Mat,PetscBool );
95: /*64*/
96: PetscErrorCode (*dummy8)(Mat,Vec,Vec);
97: PetscErrorCode (*dummy9)(Mat,Vec,Vec);
98: PetscErrorCode (*setlocaltoglobalmapping)(Mat,ISLocalToGlobalMapping,ISLocalToGlobalMapping);
99: PetscErrorCode (*setvalueslocal)(Mat,PetscInt,const PetscInt[],PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
100: PetscErrorCode (*zerorowslocal)(Mat,PetscInt,const PetscInt[],PetscScalar,Vec,Vec);
101: /*69*/
102: PetscErrorCode (*getrowmaxabs)(Mat,Vec,PetscInt[]);
103: PetscErrorCode (*getrowminabs)(Mat,Vec,PetscInt[]);
104: PetscErrorCode (*convert)(Mat, const MatType,MatReuse,Mat*);
105: PetscErrorCode (*setcoloring)(Mat,ISColoring);
106: PetscErrorCode (*setvaluesadic)(Mat,void*);
107: /*74*/
108: PetscErrorCode (*setvaluesadifor)(Mat,PetscInt,void*);
109: PetscErrorCode (*fdcoloringapply)(Mat,MatFDColoring,Vec,MatStructure*,void*);
110: PetscErrorCode (*setfromoptions)(Mat);
111: PetscErrorCode (*multconstrained)(Mat,Vec,Vec);
112: PetscErrorCode (*multtransposeconstrained)(Mat,Vec,Vec);
113: /*79*/
114: PetscErrorCode (*findzerodiagonals)(Mat,IS*);
115: PetscErrorCode (*mults)(Mat, Vecs, Vecs);
116: PetscErrorCode (*solves)(Mat, Vecs, Vecs);
117: PetscErrorCode (*getinertia)(Mat,PetscInt*,PetscInt*,PetscInt*);
118: PetscErrorCode (*load)(Mat, PetscViewer);
119: /*84*/
120: PetscErrorCode (*issymmetric)(Mat,PetscReal,PetscBool *);
121: PetscErrorCode (*ishermitian)(Mat,PetscReal,PetscBool *);
122: PetscErrorCode (*isstructurallysymmetric)(Mat,PetscBool *);
123: PetscErrorCode (*setvaluesblockedlocal)(Mat,PetscInt,const PetscInt[],PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
124: PetscErrorCode (*getvecs)(Mat,Vec*,Vec*);
125: /*89*/
126: PetscErrorCode (*matmult)(Mat,Mat,MatReuse,PetscReal,Mat*);
127: PetscErrorCode (*matmultsymbolic)(Mat,Mat,PetscReal,Mat*);
128: PetscErrorCode (*matmultnumeric)(Mat,Mat,Mat);
129: PetscErrorCode (*ptap)(Mat,Mat,MatReuse,PetscReal,Mat*);
130: PetscErrorCode (*ptapsymbolic)(Mat,Mat,PetscReal,Mat*); /* double dispatch wrapper routine */
131: /*94*/
132: PetscErrorCode (*ptapnumeric)(Mat,Mat,Mat); /* double dispatch wrapper routine */
133: PetscErrorCode (*mattransposemult)(Mat,Mat,MatReuse,PetscReal,Mat*);
134: PetscErrorCode (*mattransposemultsymbolic)(Mat,Mat,PetscReal,Mat*);
135: PetscErrorCode (*mattransposemultnumeric)(Mat,Mat,Mat);
136: PetscErrorCode (*ptapsymbolic_seqaij)(Mat,Mat,PetscReal,Mat*); /* actual implememtation, A=seqaij */
137: /*99*/
138: PetscErrorCode (*ptapnumeric_seqaij)(Mat,Mat,Mat); /* actual implememtation, A=seqaij */
139: PetscErrorCode (*ptapsymbolic_mpiaij)(Mat,Mat,PetscReal,Mat*); /* actual implememtation, A=mpiaij */
140: PetscErrorCode (*ptapnumeric_mpiaij)(Mat,Mat,Mat); /* actual implememtation, A=mpiaij */
141: PetscErrorCode (*conjugate)(Mat); /* complex conjugate */
142: PetscErrorCode (*dummy5)(void);
143: /*104*/
144: PetscErrorCode (*setvaluesrow)(Mat,PetscInt,const PetscScalar[]);
145: PetscErrorCode (*realpart)(Mat);
146: PetscErrorCode (*imaginarypart)(Mat);
147: PetscErrorCode (*getrowuppertriangular)(Mat);
148: PetscErrorCode (*restorerowuppertriangular)(Mat);
149: /*109*/
150: PetscErrorCode (*matsolve)(Mat,Mat,Mat);
151: PetscErrorCode (*getredundantmatrix)(Mat,PetscInt,MPI_Comm,PetscInt,MatReuse,Mat*);
152: PetscErrorCode (*getrowmin)(Mat,Vec,PetscInt[]);
153: PetscErrorCode (*getcolumnvector)(Mat,Vec,PetscInt);
154: PetscErrorCode (*missingdiagonal)(Mat,PetscBool *,PetscInt*);
155: /*114*/
156: PetscErrorCode (*getseqnonzerostructure)(Mat,Mat *);
157: PetscErrorCode (*create)(Mat);
158: PetscErrorCode (*getghosts)(Mat,PetscInt*,const PetscInt *[]);
159: PetscErrorCode (*getlocalsubmatrix)(Mat,IS,IS,Mat*);
160: PetscErrorCode (*restorelocalsubmatrix)(Mat,IS,IS,Mat*);
161: /*119*/
162: PetscErrorCode (*multdiagonalblock)(Mat,Vec,Vec);
163: PetscErrorCode (*hermitiantranspose)(Mat,MatReuse,Mat*);
164: PetscErrorCode (*multhermitiantranspose)(Mat,Vec,Vec);
165: PetscErrorCode (*multhermitiantransposeadd)(Mat,Vec,Vec,Vec);
166: PetscErrorCode (*getmultiprocblock)(Mat,MPI_Comm,MatReuse,Mat*);
167: /*124*/
168: PetscErrorCode (*findnonzerorows)(Mat,IS*);
169: PetscErrorCode (*getcolumnnorms)(Mat,NormType,PetscReal*);
170: PetscErrorCode (*invertblockdiagonal)(Mat,const PetscScalar**);
171: PetscErrorCode (*dummy4)(Mat,Vec,Vec,Vec);
172: PetscErrorCode (*getsubmatricesparallel)(Mat,PetscInt,const IS[], const IS[], MatReuse, Mat**);
173: /*129*/
174: PetscErrorCode (*setvaluesbatch)(Mat,PetscInt,PetscInt,PetscInt*,const PetscScalar*);
175: PetscErrorCode (*transposematmult)(Mat,Mat,MatReuse,PetscReal,Mat*);
176: PetscErrorCode (*transposematmultsymbolic)(Mat,Mat,PetscReal,Mat*);
177: PetscErrorCode (*transposematmultnumeric)(Mat,Mat,Mat);
178: PetscErrorCode (*transposecoloringcreate)(Mat,ISColoring,MatTransposeColoring);
179: /*134*/
180: PetscErrorCode (*transcoloringapplysptoden)(MatTransposeColoring,Mat,Mat);
181: PetscErrorCode (*transcoloringapplydentosp)(MatTransposeColoring,Mat,Mat);
182: PetscErrorCode (*rart)(Mat,Mat,MatReuse,PetscReal,Mat*);
183: PetscErrorCode (*rartsymbolic)(Mat,Mat,PetscReal,Mat*); /* double dispatch wrapper routine */
184: PetscErrorCode (*rartnumeric)(Mat,Mat,Mat); /* double dispatch wrapper routine */
185: };
186: /*
187: If you add MatOps entries above also add them to the MATOP enum
188: in include/petscmat.h and include/finclude/petscmat.h
189: */
191: #include <petscsys.h>
192: PETSC_EXTERN PetscErrorCode MatRegisterOp(MPI_Comm, const char[], PetscVoidFunction, const char[], PetscInt, ...);
193: PETSC_EXTERN PetscErrorCode MatQueryOp(MPI_Comm, PetscVoidFunction*, const char[], PetscInt, ...);
195: typedef struct _p_MatBaseName* MatBaseName;
196: struct _p_MatBaseName {
197: char *bname,*sname,*mname;
198: MatBaseName next;
199: };
201: PETSC_EXTERN MatBaseName MatBaseNameList;
203: /*
204: Utility private matrix routines
205: */
206: PETSC_EXTERN PetscErrorCode MatConvert_Basic(Mat, const MatType,MatReuse,Mat*);
207: PETSC_EXTERN PetscErrorCode MatCopy_Basic(Mat,Mat,MatStructure);
208: PETSC_EXTERN PetscErrorCode MatView_Private(Mat);
210: PETSC_EXTERN PetscErrorCode MatHeaderMerge(Mat,Mat);
211: PETSC_EXTERN PetscErrorCode MatHeaderReplace(Mat,Mat);
212: PETSC_EXTERN PetscErrorCode MatAXPYGetxtoy_Private(PetscInt,PetscInt*,PetscInt*,PetscInt*, PetscInt*,PetscInt*,PetscInt*, PetscInt**);
213: PETSC_EXTERN PetscErrorCode MatPtAP_Basic(Mat,Mat,MatReuse,PetscReal,Mat*);
214: PETSC_EXTERN PetscErrorCode MatDiagonalSet_Default(Mat,Vec,InsertMode);
216: #if defined(PETSC_USE_DEBUG)
217: # define MatCheckPreallocated(A,arg) do { \
218: if (PetscUnlikely(!(A)->preallocated)) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call MatXXXSetPreallocation() or MatSetUp() on argument %D \"%s\" before %s()",(arg),#A,PETSC_FUNCTION_NAME); \
219: } while (0)
220: #else
221: # define MatCheckPreallocated(A,arg) do {} while (0)
222: #endif
224: /*
225: The stash is used to temporarily store inserted matrix values that
226: belong to another processor. During the assembly phase the stashed
227: values are moved to the correct processor and
228: */
230: typedef struct _MatStashSpace *PetscMatStashSpace;
232: struct _MatStashSpace {
233: PetscMatStashSpace next;
234: PetscScalar *space_head,*val;
235: PetscInt *idx,*idy;
236: PetscInt total_space_size;
237: PetscInt local_used;
238: PetscInt local_remaining;
239: };
241: PETSC_EXTERN PetscErrorCode PetscMatStashSpaceGet(PetscInt,PetscInt,PetscMatStashSpace *);
242: PETSC_EXTERN PetscErrorCode PetscMatStashSpaceContiguous(PetscInt,PetscMatStashSpace *,PetscScalar *,PetscInt *,PetscInt *);
243: PETSC_EXTERN PetscErrorCode PetscMatStashSpaceDestroy(PetscMatStashSpace*);
245: typedef struct {
246: PetscInt nmax; /* maximum stash size */
247: PetscInt umax; /* user specified max-size */
248: PetscInt oldnmax; /* the nmax value used previously */
249: PetscInt n; /* stash size */
250: PetscInt bs; /* block size of the stash */
251: PetscInt reallocs; /* preserve the no of mallocs invoked */
252: PetscMatStashSpace space_head,space; /* linked list to hold stashed global row/column numbers and matrix values */
253: /* The following variables are used for communication */
254: MPI_Comm comm;
255: PetscMPIInt size,rank;
256: PetscMPIInt tag1,tag2;
257: MPI_Request *send_waits; /* array of send requests */
258: MPI_Request *recv_waits; /* array of receive requests */
259: MPI_Status *send_status; /* array of send status */
260: PetscInt nsends,nrecvs; /* numbers of sends and receives */
261: PetscScalar *svalues; /* sending data */
262: PetscInt *sindices;
263: PetscScalar **rvalues; /* receiving data (values) */
264: PetscInt **rindices; /* receiving data (indices) */
265: PetscInt nprocessed; /* number of messages already processed */
266: PetscMPIInt *flg_v; /* indicates what messages have arrived so far and from whom */
267: PetscBool reproduce;
268: PetscInt reproduce_count;
269: } MatStash;
271: PETSC_EXTERN PetscErrorCode MatStashCreate_Private(MPI_Comm,PetscInt,MatStash*);
272: PETSC_EXTERN PetscErrorCode MatStashDestroy_Private(MatStash*);
273: PETSC_EXTERN PetscErrorCode MatStashScatterEnd_Private(MatStash*);
274: PETSC_EXTERN PetscErrorCode MatStashSetInitialSize_Private(MatStash*,PetscInt);
275: PETSC_EXTERN PetscErrorCode MatStashGetInfo_Private(MatStash*,PetscInt*,PetscInt*);
276: PETSC_EXTERN PetscErrorCode MatStashValuesRow_Private(MatStash*,PetscInt,PetscInt,const PetscInt[],const PetscScalar[],PetscBool );
277: PETSC_EXTERN PetscErrorCode MatStashValuesCol_Private(MatStash*,PetscInt,PetscInt,const PetscInt[],const PetscScalar[],PetscInt,PetscBool );
278: PETSC_EXTERN PetscErrorCode MatStashValuesRowBlocked_Private(MatStash*,PetscInt,PetscInt,const PetscInt[],const PetscScalar[],PetscInt,PetscInt,PetscInt);
279: PETSC_EXTERN PetscErrorCode MatStashValuesColBlocked_Private(MatStash*,PetscInt,PetscInt,const PetscInt[],const PetscScalar[],PetscInt,PetscInt,PetscInt);
280: PETSC_EXTERN PetscErrorCode MatStashScatterBegin_Private(Mat,MatStash*,PetscInt*);
281: PETSC_EXTERN PetscErrorCode MatStashScatterGetMesg_Private(MatStash*,PetscMPIInt*,PetscInt**,PetscInt**,PetscScalar**,PetscInt*);
283: typedef struct {
284: PetscInt dim;
285: PetscInt dims[4];
286: PetscInt starts[4];
287: PetscBool noc; /* this is a single component problem, hence user will not set MatStencil.c */
288: } MatStencilInfo;
290: /* Info about using compressed row format */
291: typedef struct {
292: PetscBool check; /* indicates that at MatAssembly() it should check if compressed rows will be efficient */
293: PetscBool use; /* indicates compressed rows have been checked and will be used */
294: PetscInt nrows; /* number of non-zero rows */
295: PetscInt *i; /* compressed row pointer */
296: PetscInt *rindex; /* compressed row index */
297: } Mat_CompressedRow;
298: PETSC_EXTERN PetscErrorCode MatCheckCompressedRow(Mat,Mat_CompressedRow*,PetscInt*,PetscInt,PetscReal);
300: struct _p_Mat {
301: PETSCHEADER(struct _MatOps);
302: PetscLayout rmap,cmap;
303: void *data; /* implementation-specific data */
304: MatFactorType factortype; /* MAT_FACTOR_LU, ILU, CHOLESKY or ICC */
305: PetscBool assembled; /* is the matrix assembled? */
306: PetscBool was_assembled; /* new values inserted into assembled mat */
307: PetscInt num_ass; /* number of times matrix has been assembled */
308: PetscBool same_nonzero; /* matrix has same nonzero pattern as previous */
309: MatInfo info; /* matrix information */
310: InsertMode insertmode; /* have values been inserted in matrix or added? */
311: MatStash stash,bstash; /* used for assembling off-proc mat emements */
312: MatNullSpace nullsp; /* null space (operator is singular) */
313: MatNullSpace nearnullsp; /* near null space to be used by multigrid methods */
314: PetscBool preallocated;
315: MatStencilInfo stencil; /* information for structured grid */
316: PetscBool symmetric,hermitian,structurally_symmetric,spd;
317: PetscBool symmetric_set,hermitian_set,structurally_symmetric_set,spd_set; /* if true, then corresponding flag is correct*/
318: PetscBool symmetric_eternal;
319: PetscBool nooffprocentries,nooffproczerorows;
320: #if defined(PETSC_HAVE_CUSP)
321: PetscCUSPFlag valid_GPU_matrix; /* flag pointing to the matrix on the gpu*/
322: #endif
323: void *spptr; /* pointer for special library like SuperLU */
324: MatSolverPackage solvertype;
325: };
327: PETSC_EXTERN PetscErrorCode MatAXPY_Basic(Mat,PetscScalar,Mat,MatStructure);
328: PETSC_EXTERN PetscErrorCode MatAXPY_BasicWithPreallocation(Mat,Mat,PetscScalar,Mat,MatStructure);
329: /*
330: Object for partitioning graphs
331: */
333: typedef struct _MatPartitioningOps *MatPartitioningOps;
334: struct _MatPartitioningOps {
335: PetscErrorCode (*apply)(MatPartitioning,IS*);
336: PetscErrorCode (*setfromoptions)(MatPartitioning);
337: PetscErrorCode (*destroy)(MatPartitioning);
338: PetscErrorCode (*view)(MatPartitioning,PetscViewer);
339: };
341: struct _p_MatPartitioning {
342: PETSCHEADER(struct _MatPartitioningOps);
343: Mat adj;
344: PetscInt *vertex_weights;
345: PetscReal *part_weights;
346: PetscInt n; /* number of partitions */
347: void *data;
348: PetscInt setupcalled;
349: };
351: /*
352: Object for coarsen graphs
353: */
354: typedef struct _MatCoarsenOps *MatCoarsenOps;
355: struct _MatCoarsenOps {
356: PetscErrorCode (*apply)(MatCoarsen);
357: PetscErrorCode (*setfromoptions)(MatCoarsen);
358: PetscErrorCode (*destroy)(MatCoarsen);
359: PetscErrorCode (*view)(MatCoarsen,PetscViewer);
360: };
362: struct _p_MatCoarsen {
363: PETSCHEADER(struct _MatCoarsenOps);
364: Mat graph;
365: PetscInt verbose;
366: PetscInt setupcalled;
367: void *subctx;
368: /* */
369: PetscBool strict_aggs;
370: IS perm;
371: PetscCoarsenData *agg_lists;
372: };
374: PetscErrorCode PetscCDCreate( PetscInt chsz, PetscCoarsenData **ail );
375: PetscErrorCode PetscCDDestroy( PetscCoarsenData *ail );
376: PetscErrorCode PetscLLNSetID( PetscCDIntNd *a_this, PetscInt a_gid );
377: PetscErrorCode PetscLLNGetID( const PetscCDIntNd *a_this, PetscInt *a_gid );
378: PetscErrorCode PetscCDAppendID( PetscCoarsenData *ail, PetscInt a_idx, PetscInt a_gid );
379: PetscErrorCode PetscCDAppendRemove( PetscCoarsenData *ail, PetscInt a_destidx, PetscInt a_srcidx );
380: PetscErrorCode PetscCDAppendNode( PetscCoarsenData *ail, PetscInt a_idx, PetscCDIntNd *a_n );
381: PetscErrorCode PetscCDRemoveNextNode( PetscCoarsenData *ail, PetscInt a_idx, PetscCDIntNd *a_last );
382: PetscErrorCode PetscCDRemoveAllAt( PetscCoarsenData *ail, PetscInt a_idx );
383: PetscErrorCode PetscCDSizeAt( const PetscCoarsenData *ail, PetscInt a_idx, PetscInt *a_sz );
384: PetscErrorCode PetscCDEmptyAt( const PetscCoarsenData *ail, PetscInt a_idx, PetscBool *a_empty );
385: PetscErrorCode PetscCDSetChuckSize( PetscCoarsenData *ail, PetscInt a_sz );
386: PetscErrorCode PetscCDPrint( const PetscCoarsenData *ail, MPI_Comm comm );
387: PetscErrorCode PetscCDGetMIS( PetscCoarsenData *ail, IS * );
388: PetscErrorCode PetscCDGetMat( const PetscCoarsenData *ail, Mat * );
389: PetscErrorCode PetscCDSetMat( PetscCoarsenData *ail, Mat );
390: typedef PetscCDIntNd* PetscCDPos;
391: PetscErrorCode PetscCDGetHeadPos( const PetscCoarsenData *ail, PetscInt idx, PetscCDPos *cpos );
392: PetscErrorCode PetscCDGetNextPos( const PetscCoarsenData *ail, PetscInt idx, PetscCDPos *cpos );
393: PetscErrorCode PetscCDGetASMBlocks( const PetscCoarsenData *ail, const PetscInt, PetscInt *, IS** );
394: PetscErrorCode PetscCDSetRemovedIS( PetscCoarsenData *ail, MPI_Comm, const PetscInt, PetscInt[] );
395: PetscErrorCode PetscCDGetRemovedIS( PetscCoarsenData *ail, IS * );
397: /*
398: MatFDColoring is used to compute Jacobian matrices efficiently
399: via coloring. The data structure is explained below in an example.
401: Color = 0 1 0 2 | 2 3 0
402: ---------------------------------------------------
403: 00 01 | 05
404: 10 11 | 14 15 Processor 0
405: 22 23 | 25
406: 32 33 |
407: ===================================================
408: | 44 45 46
409: 50 | 55 Processor 1
410: | 64 66
411: ---------------------------------------------------
413: ncolors = 4;
415: ncolumns = {2,1,1,0}
416: columns = {{0,2},{1},{3},{}}
417: nrows = {4,2,3,3}
418: rows = {{0,1,2,3},{0,1},{1,2,3},{0,1,2}}
419: columnsforrow = {{0,0,2,2},{1,1},{4,3,3},{5,5,5}}
420: vscaleforrow = {{,,,},{,},{,,},{,,}}
421: vwscale = {dx(0),dx(1),dx(2),dx(3)} MPI Vec
422: vscale = {dx(0),dx(1),dx(2),dx(3),dx(4),dx(5)} Seq Vec
424: ncolumns = {1,0,1,1}
425: columns = {{6},{},{4},{5}}
426: nrows = {3,0,2,2}
427: rows = {{0,1,2},{},{1,2},{1,2}}
428: columnsforrow = {{6,0,6},{},{4,4},{5,5}}
429: vscaleforrow = {{,,},{},{,},{,}}
430: vwscale = {dx(4),dx(5),dx(6)} MPI Vec
431: vscale = {dx(0),dx(4),dx(5),dx(6)} Seq Vec
433: See the routine MatFDColoringApply() for how this data is used
434: to compute the Jacobian.
436: */
438: struct _p_MatFDColoring{
439: PETSCHEADER(int);
440: PetscInt M,N,m; /* total rows, columns; local rows */
441: PetscInt rstart; /* first row owned by local processor */
442: PetscInt ncolors; /* number of colors */
443: PetscInt *ncolumns; /* number of local columns for a color */
444: PetscInt **columns; /* lists the local columns of each color (using global column numbering) */
445: PetscInt *nrows; /* number of local rows for each color */
446: PetscInt **rows; /* lists the local rows for each color (using the local row numbering) */
447: PetscInt **columnsforrow; /* lists the corresponding columns for those rows (using the global column) */
448: PetscReal error_rel; /* square root of relative error in computing function */
449: PetscReal umin; /* minimum allowable u'dx value */
450: Vec w1,w2,w3; /* work vectors used in computing Jacobian */
451: PetscErrorCode (*f)(void); /* function that defines Jacobian */
452: void *fctx; /* optional user-defined context for use by the function f */
453: PetscInt **vscaleforrow; /* location in vscale for each columnsforrow[] entry */
454: Vec vscale; /* holds FD scaling, i.e. 1/dx for each perturbed column */
455: Vec F; /* current value of user provided function; can set with MatFDColoringSetF() */
456: PetscInt currentcolor; /* color for which function evaluation is being done now */
457: const char *htype; /* "wp" or "ds" */
458: ISColoringType ctype; /* IS_COLORING_GLOBAL or IS_COLORING_GHOSTED */
460: void *ftn_func_pointer,*ftn_func_cntx; /* serve the same purpose as *fortran_func_pointers in PETSc objects */
461: };
463: struct _p_MatTransposeColoring{
464: PETSCHEADER(int);
465: PetscInt M,N,m; /* total rows, columns; local rows */
466: PetscInt rstart; /* first row owned by local processor */
467: PetscInt ncolors; /* number of colors */
468: PetscInt *ncolumns; /* number of local columns for a color */
469: PetscInt *nrows; /* number of local rows for each color */
470: PetscInt currentcolor; /* color for which function evaluation is being done now */
471: ISColoringType ctype; /* IS_COLORING_GLOBAL or IS_COLORING_GHOSTED */
472:
473: PetscInt *colorforrow,*colorforcol; /* pointer to rows and columns */
474: PetscInt *rows; /* lists the local rows for each color (using the local row numbering) */
475: PetscInt *columnsforspidx; /* maps (row,color) in the dense matrix to index of sparse matrix arrays a->j and a->a */
476: PetscInt *columns; /* lists the local columns of each color (using global column numbering) */
477: };
479: /*
480: Null space context for preconditioner/operators
481: */
482: struct _p_MatNullSpace {
483: PETSCHEADER(int);
484: PetscBool has_cnst;
485: PetscInt n;
486: Vec* vecs;
487: PetscScalar* alpha; /* for projections */
488: Vec vec; /* for out of place removals */
489: PetscErrorCode (*remove)(MatNullSpace,Vec,void*); /* for user provided removal function */
490: void* rmctx; /* context for remove() function */
491: };
493: /*
494: Checking zero pivot for LU, ILU preconditioners.
495: */
496: typedef struct {
497: PetscInt nshift,nshift_max;
498: PetscReal shift_amount,shift_lo,shift_hi,shift_top,shift_fraction;
499: PetscBool newshift;
500: PetscReal rs; /* active row sum of abs(offdiagonals) */
501: PetscScalar pv; /* pivot of the active row */
502: } FactorShiftCtx;
504: PETSC_EXTERN PetscErrorCode MatFactorDumpMatrix(Mat);
508: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_nz(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
509: {
510: PetscReal _rs = sctx->rs;
511: PetscReal _zero = info->zeropivot*_rs;
514: if (PetscAbsScalar(sctx->pv) <= _zero){
515: /* force |diag| > zeropivot*rs */
516: if (!sctx->nshift) sctx->shift_amount = info->shiftamount;
517: else sctx->shift_amount *= 2.0;
518: sctx->newshift = PETSC_TRUE;
519: (sctx->nshift)++;
520: } else {
521: sctx->newshift = PETSC_FALSE;
522: }
523: return(0);
524: }
528: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_pd(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
529: {
530: PetscReal _rs = sctx->rs;
531: PetscReal _zero = info->zeropivot*_rs;
534: if (PetscRealPart(sctx->pv) <= _zero){
535: /* force matfactor to be diagonally dominant */
536: if (sctx->nshift == sctx->nshift_max) {
537: sctx->shift_fraction = sctx->shift_hi;
538: } else {
539: sctx->shift_lo = sctx->shift_fraction;
540: sctx->shift_fraction = (sctx->shift_hi+sctx->shift_lo)/2.;
541: }
542: sctx->shift_amount = sctx->shift_fraction * sctx->shift_top;
543: sctx->nshift++;
544: sctx->newshift = PETSC_TRUE;
545: } else {
546: sctx->newshift = PETSC_FALSE;
547: }
548: return(0);
549: }
553: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_inblocks(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
554: {
555: PetscReal _zero = info->zeropivot;
558: if (PetscAbsScalar(sctx->pv) <= _zero){
559: sctx->pv += info->shiftamount;
560: sctx->shift_amount = 0.0;
561: sctx->nshift++;
562: }
563: sctx->newshift = PETSC_FALSE;
564: return(0);
565: }
569: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_none(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
570: {
571: PetscReal _zero = info->zeropivot;
574: sctx->newshift = PETSC_FALSE;
575: if (PetscAbsScalar(sctx->pv) <= _zero) {
577: PetscBool flg = PETSC_FALSE;
578:
579: PetscOptionsGetBool(PETSC_NULL,"-mat_dump",&flg,PETSC_NULL);
580: if (flg) {
581: MatView(mat,PETSC_VIEWER_BINARY_(((PetscObject)mat)->comm));
582: }
583: SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_MAT_LU_ZRPVT,"Zero pivot row %D value %G tolerance %G",row,PetscAbsScalar(sctx->pv),_zero);
584: }
585: return(0);
586: }
590: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
591: {
595: if (info->shifttype == (PetscReal) MAT_SHIFT_NONZERO){
596: MatPivotCheck_nz(mat,info,sctx,row);
597: } else if (info->shifttype == (PetscReal) MAT_SHIFT_POSITIVE_DEFINITE){
598: MatPivotCheck_pd(mat,info,sctx,row);
599: } else if (info->shifttype == (PetscReal) MAT_SHIFT_INBLOCKS){
600: MatPivotCheck_inblocks(mat,info,sctx,row);
601: } else {
602: MatPivotCheck_none(mat,info,sctx,row);
603: }
604: return(0);
605: }
607: /*
608: Create and initialize a linked list
609: Input Parameters:
610: idx_start - starting index of the list
611: lnk_max - max value of lnk indicating the end of the list
612: nlnk - max length of the list
613: Output Parameters:
614: lnk - list initialized
615: bt - PetscBT (bitarray) with all bits set to false
616: lnk_empty - flg indicating the list is empty
617: */
618: #define PetscLLCreate(idx_start,lnk_max,nlnk,lnk,bt) \
619: (PetscMalloc(nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,&(bt)) || (lnk[idx_start] = lnk_max,0))
621: #define PetscLLCreate_new(idx_start,lnk_max,nlnk,lnk,bt,lnk_empty)\
622: (PetscMalloc(nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,&(bt)) || (lnk_empty = PETSC_TRUE,0) ||(lnk[idx_start] = lnk_max,0))
624: /*
625: Add an index set into a sorted linked list
626: Input Parameters:
627: nidx - number of input indices
628: indices - interger array
629: idx_start - starting index of the list
630: lnk - linked list(an integer array) that is created
631: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
632: output Parameters:
633: nlnk - number of newly added indices
634: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices
635: bt - updated PetscBT (bitarray)
636: */
637: #define PetscLLAdd(nidx,indices,idx_start,nlnk,lnk,bt) 0;\
638: {\
639: PetscInt _k,_entry,_location,_lnkdata;\
640: nlnk = 0;\
641: _lnkdata = idx_start;\
642: for (_k=0; _k<nidx; _k++){\
643: _entry = indices[_k];\
644: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
645: /* search for insertion location */\
646: /* start from the beginning if _entry < previous _entry */\
647: if (_k && _entry < _lnkdata) _lnkdata = idx_start;\
648: do {\
649: _location = _lnkdata;\
650: _lnkdata = lnk[_location];\
651: } while (_entry > _lnkdata);\
652: /* insertion location is found, add entry into lnk */\
653: lnk[_location] = _entry;\
654: lnk[_entry] = _lnkdata;\
655: nlnk++;\
656: _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
657: }\
658: }\
659: }
661: /*
662: Add a permuted index set into a sorted linked list
663: Input Parameters:
664: nidx - number of input indices
665: indices - interger array
666: perm - permutation of indices
667: idx_start - starting index of the list
668: lnk - linked list(an integer array) that is created
669: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
670: output Parameters:
671: nlnk - number of newly added indices
672: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices
673: bt - updated PetscBT (bitarray)
674: */
675: #define PetscLLAddPerm(nidx,indices,perm,idx_start,nlnk,lnk,bt) 0;\
676: {\
677: PetscInt _k,_entry,_location,_lnkdata;\
678: nlnk = 0;\
679: _lnkdata = idx_start;\
680: for (_k=0; _k<nidx; _k++){\
681: _entry = perm[indices[_k]];\
682: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
683: /* search for insertion location */\
684: /* start from the beginning if _entry < previous _entry */\
685: if (_k && _entry < _lnkdata) _lnkdata = idx_start;\
686: do {\
687: _location = _lnkdata;\
688: _lnkdata = lnk[_location];\
689: } while (_entry > _lnkdata);\
690: /* insertion location is found, add entry into lnk */\
691: lnk[_location] = _entry;\
692: lnk[_entry] = _lnkdata;\
693: nlnk++;\
694: _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
695: }\
696: }\
697: }
699: /*
700: Add a SORTED ascending index set into a sorted linked list - same as PetscLLAdd() bus skip 'if (_k && _entry < _lnkdata) _lnkdata = idx_start;'
701: Input Parameters:
702: nidx - number of input indices
703: indices - sorted interger array
704: idx_start - starting index of the list
705: lnk - linked list(an integer array) that is created
706: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
707: output Parameters:
708: nlnk - number of newly added indices
709: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices
710: bt - updated PetscBT (bitarray)
711: */
712: #define PetscLLAddSorted(nidx,indices,idx_start,nlnk,lnk,bt) 0;\
713: {\
714: PetscInt _k,_entry,_location,_lnkdata;\
715: nlnk = 0;\
716: _lnkdata = idx_start;\
717: for (_k=0; _k<nidx; _k++){\
718: _entry = indices[_k];\
719: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
720: /* search for insertion location */\
721: do {\
722: _location = _lnkdata;\
723: _lnkdata = lnk[_location];\
724: } while (_entry > _lnkdata);\
725: /* insertion location is found, add entry into lnk */\
726: lnk[_location] = _entry;\
727: lnk[_entry] = _lnkdata;\
728: nlnk++;\
729: _lnkdata = _entry; /* next search starts from here */\
730: }\
731: }\
732: }
734: #define PetscLLAddSorted_new(nidx,indices,idx_start,lnk_empty,nlnk,lnk,bt) 0; \
735: {\
736: PetscInt _k,_entry,_location,_lnkdata;\
737: if (lnk_empty){\
738: _lnkdata = idx_start; \
739: for (_k=0; _k<nidx; _k++){ \
740: _entry = indices[_k]; \
741: PetscBTSet(bt,_entry); /* mark the new entry */ \
742: _location = _lnkdata; \
743: _lnkdata = lnk[_location]; \
744: /* insertion location is found, add entry into lnk */ \
745: lnk[_location] = _entry; \
746: lnk[_entry] = _lnkdata; \
747: _lnkdata = _entry; /* next search starts from here */ \
748: } \
749: /*\
750: lnk[indices[nidx-1]] = lnk[idx_start];\
751: lnk[idx_start] = indices[0];\
752: PetscBTSet(bt,indices[0]); \
753: for (_k=1; _k<nidx; _k++){ \
754: PetscBTSet(bt,indices[_k]); \
755: lnk[indices[_k-1]] = indices[_k]; \
756: } \
757: */\
758: nlnk = nidx;\
759: lnk_empty = PETSC_FALSE;\
760: } else {\
761: nlnk = 0; \
762: _lnkdata = idx_start; \
763: for (_k=0; _k<nidx; _k++){ \
764: _entry = indices[_k]; \
765: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */ \
766: /* search for insertion location */ \
767: do { \
768: _location = _lnkdata; \
769: _lnkdata = lnk[_location]; \
770: } while (_entry > _lnkdata); \
771: /* insertion location is found, add entry into lnk */ \
772: lnk[_location] = _entry; \
773: lnk[_entry] = _lnkdata; \
774: nlnk++; \
775: _lnkdata = _entry; /* next search starts from here */ \
776: } \
777: } \
778: } \
779: }
781: /*
782: Add a SORTED index set into a sorted linked list used for LUFactorSymbolic()
783: Same as PetscLLAddSorted() with an additional operation:
784: count the number of input indices that are no larger than 'diag'
785: Input Parameters:
786: indices - sorted interger array
787: idx_start - starting index of the list, index of pivot row
788: lnk - linked list(an integer array) that is created
789: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
790: diag - index of the active row in LUFactorSymbolic
791: nzbd - number of input indices with indices <= idx_start
792: im - im[idx_start] is initialized as num of nonzero entries in row=idx_start
793: output Parameters:
794: nlnk - number of newly added indices
795: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices
796: bt - updated PetscBT (bitarray)
797: im - im[idx_start]: unchanged if diag is not an entry
798: : num of entries with indices <= diag if diag is an entry
799: */
800: #define PetscLLAddSortedLU(indices,idx_start,nlnk,lnk,bt,diag,nzbd,im) 0;\
801: {\
802: PetscInt _k,_entry,_location,_lnkdata,_nidx;\
803: nlnk = 0;\
804: _lnkdata = idx_start;\
805: _nidx = im[idx_start] - nzbd; /* num of entries with idx_start < index <= diag */\
806: for (_k=0; _k<_nidx; _k++){\
807: _entry = indices[_k];\
808: nzbd++;\
809: if ( _entry== diag) im[idx_start] = nzbd;\
810: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
811: /* search for insertion location */\
812: do {\
813: _location = _lnkdata;\
814: _lnkdata = lnk[_location];\
815: } while (_entry > _lnkdata);\
816: /* insertion location is found, add entry into lnk */\
817: lnk[_location] = _entry;\
818: lnk[_entry] = _lnkdata;\
819: nlnk++;\
820: _lnkdata = _entry; /* next search starts from here */\
821: }\
822: }\
823: }
825: /*
826: Copy data on the list into an array, then initialize the list
827: Input Parameters:
828: idx_start - starting index of the list
829: lnk_max - max value of lnk indicating the end of the list
830: nlnk - number of data on the list to be copied
831: lnk - linked list
832: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
833: output Parameters:
834: indices - array that contains the copied data
835: lnk - linked list that is cleaned and initialize
836: bt - PetscBT (bitarray) with all bits set to false
837: */
838: #define PetscLLClean(idx_start,lnk_max,nlnk,lnk,indices,bt) 0;\
839: {\
840: PetscInt _j,_idx=idx_start;\
841: for (_j=0; _j<nlnk; _j++){\
842: _idx = lnk[_idx];\
843: indices[_j] = _idx;\
844: PetscBTClear(bt,_idx);\
845: }\
846: lnk[idx_start] = lnk_max;\
847: }
848: /*
849: Free memories used by the list
850: */
851: #define PetscLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(&(bt)))
853: /* Routines below are used for incomplete matrix factorization */
854: /*
855: Create and initialize a linked list and its levels
856: Input Parameters:
857: idx_start - starting index of the list
858: lnk_max - max value of lnk indicating the end of the list
859: nlnk - max length of the list
860: Output Parameters:
861: lnk - list initialized
862: lnk_lvl - array of size nlnk for storing levels of lnk
863: bt - PetscBT (bitarray) with all bits set to false
864: */
865: #define PetscIncompleteLLCreate(idx_start,lnk_max,nlnk,lnk,lnk_lvl,bt)\
866: (PetscMalloc(2*nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,&(bt)) || (lnk[idx_start] = lnk_max,lnk_lvl = lnk + nlnk,0))
868: /*
869: Initialize a sorted linked list used for ILU and ICC
870: Input Parameters:
871: nidx - number of input idx
872: idx - interger array used for storing column indices
873: idx_start - starting index of the list
874: perm - indices of an IS
875: lnk - linked list(an integer array) that is created
876: lnklvl - levels of lnk
877: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
878: output Parameters:
879: nlnk - number of newly added idx
880: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx
881: lnklvl - levels of lnk
882: bt - updated PetscBT (bitarray)
883: */
884: #define PetscIncompleteLLInit(nidx,idx,idx_start,perm,nlnk,lnk,lnklvl,bt) 0;\
885: {\
886: PetscInt _k,_entry,_location,_lnkdata;\
887: nlnk = 0;\
888: _lnkdata = idx_start;\
889: for (_k=0; _k<nidx; _k++){\
890: _entry = perm[idx[_k]];\
891: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
892: /* search for insertion location */\
893: if (_k && _entry < _lnkdata) _lnkdata = idx_start;\
894: do {\
895: _location = _lnkdata;\
896: _lnkdata = lnk[_location];\
897: } while (_entry > _lnkdata);\
898: /* insertion location is found, add entry into lnk */\
899: lnk[_location] = _entry;\
900: lnk[_entry] = _lnkdata;\
901: lnklvl[_entry] = 0;\
902: nlnk++;\
903: _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
904: }\
905: }\
906: }
908: /*
909: Add a SORTED index set into a sorted linked list for ILU
910: Input Parameters:
911: nidx - number of input indices
912: idx - sorted interger array used for storing column indices
913: level - level of fill, e.g., ICC(level)
914: idxlvl - level of idx
915: idx_start - starting index of the list
916: lnk - linked list(an integer array) that is created
917: lnklvl - levels of lnk
918: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
919: prow - the row number of idx
920: output Parameters:
921: nlnk - number of newly added idx
922: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx
923: lnklvl - levels of lnk
924: bt - updated PetscBT (bitarray)
926: Note: the level of factor(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(i,prow)+lvl(prow,j)+1)
927: where idx = non-zero columns of U(prow,prow+1:n-1), prow<i
928: */
929: #define PetscILULLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,lnklvl_prow) 0;\
930: {\
931: PetscInt _k,_entry,_location,_lnkdata,_incrlev,_lnklvl_prow=lnklvl[prow];\
932: nlnk = 0;\
933: _lnkdata = idx_start;\
934: for (_k=0; _k<nidx; _k++){\
935: _incrlev = idxlvl[_k] + _lnklvl_prow + 1;\
936: if (_incrlev > level) continue;\
937: _entry = idx[_k];\
938: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
939: /* search for insertion location */\
940: do {\
941: _location = _lnkdata;\
942: _lnkdata = lnk[_location];\
943: } while (_entry > _lnkdata);\
944: /* insertion location is found, add entry into lnk */\
945: lnk[_location] = _entry;\
946: lnk[_entry] = _lnkdata;\
947: lnklvl[_entry] = _incrlev;\
948: nlnk++;\
949: _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
950: } else { /* existing entry: update lnklvl */\
951: if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
952: }\
953: }\
954: }
956: /*
957: Add a index set into a sorted linked list
958: Input Parameters:
959: nidx - number of input idx
960: idx - interger array used for storing column indices
961: level - level of fill, e.g., ICC(level)
962: idxlvl - level of idx
963: idx_start - starting index of the list
964: lnk - linked list(an integer array) that is created
965: lnklvl - levels of lnk
966: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
967: output Parameters:
968: nlnk - number of newly added idx
969: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx
970: lnklvl - levels of lnk
971: bt - updated PetscBT (bitarray)
972: */
973: #define PetscIncompleteLLAdd(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\
974: {\
975: PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
976: nlnk = 0;\
977: _lnkdata = idx_start;\
978: for (_k=0; _k<nidx; _k++){\
979: _incrlev = idxlvl[_k] + 1;\
980: if (_incrlev > level) continue;\
981: _entry = idx[_k];\
982: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
983: /* search for insertion location */\
984: if (_k && _entry < _lnkdata) _lnkdata = idx_start;\
985: do {\
986: _location = _lnkdata;\
987: _lnkdata = lnk[_location];\
988: } while (_entry > _lnkdata);\
989: /* insertion location is found, add entry into lnk */\
990: lnk[_location] = _entry;\
991: lnk[_entry] = _lnkdata;\
992: lnklvl[_entry] = _incrlev;\
993: nlnk++;\
994: _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
995: } else { /* existing entry: update lnklvl */\
996: if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
997: }\
998: }\
999: }
1001: /*
1002: Add a SORTED index set into a sorted linked list
1003: Input Parameters:
1004: nidx - number of input indices
1005: idx - sorted interger array used for storing column indices
1006: level - level of fill, e.g., ICC(level)
1007: idxlvl - level of idx
1008: idx_start - starting index of the list
1009: lnk - linked list(an integer array) that is created
1010: lnklvl - levels of lnk
1011: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
1012: output Parameters:
1013: nlnk - number of newly added idx
1014: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx
1015: lnklvl - levels of lnk
1016: bt - updated PetscBT (bitarray)
1017: */
1018: #define PetscIncompleteLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\
1019: {\
1020: PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
1021: nlnk = 0;\
1022: _lnkdata = idx_start;\
1023: for (_k=0; _k<nidx; _k++){\
1024: _incrlev = idxlvl[_k] + 1;\
1025: if (_incrlev > level) continue;\
1026: _entry = idx[_k];\
1027: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
1028: /* search for insertion location */\
1029: do {\
1030: _location = _lnkdata;\
1031: _lnkdata = lnk[_location];\
1032: } while (_entry > _lnkdata);\
1033: /* insertion location is found, add entry into lnk */\
1034: lnk[_location] = _entry;\
1035: lnk[_entry] = _lnkdata;\
1036: lnklvl[_entry] = _incrlev;\
1037: nlnk++;\
1038: _lnkdata = _entry; /* next search starts from here */\
1039: } else { /* existing entry: update lnklvl */\
1040: if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
1041: }\
1042: }\
1043: }
1045: /*
1046: Add a SORTED index set into a sorted linked list for ICC
1047: Input Parameters:
1048: nidx - number of input indices
1049: idx - sorted interger array used for storing column indices
1050: level - level of fill, e.g., ICC(level)
1051: idxlvl - level of idx
1052: idx_start - starting index of the list
1053: lnk - linked list(an integer array) that is created
1054: lnklvl - levels of lnk
1055: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
1056: idxlvl_prow - idxlvl[prow], where prow is the row number of the idx
1057: output Parameters:
1058: nlnk - number of newly added indices
1059: lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx
1060: lnklvl - levels of lnk
1061: bt - updated PetscBT (bitarray)
1062: Note: the level of U(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(prow,i)+lvl(prow,j)+1)
1063: where idx = non-zero columns of U(prow,prow+1:n-1), prow<i
1064: */
1065: #define PetscICCLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,idxlvl_prow) 0;\
1066: {\
1067: PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
1068: nlnk = 0;\
1069: _lnkdata = idx_start;\
1070: for (_k=0; _k<nidx; _k++){\
1071: _incrlev = idxlvl[_k] + idxlvl_prow + 1;\
1072: if (_incrlev > level) continue;\
1073: _entry = idx[_k];\
1074: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\
1075: /* search for insertion location */\
1076: do {\
1077: _location = _lnkdata;\
1078: _lnkdata = lnk[_location];\
1079: } while (_entry > _lnkdata);\
1080: /* insertion location is found, add entry into lnk */\
1081: lnk[_location] = _entry;\
1082: lnk[_entry] = _lnkdata;\
1083: lnklvl[_entry] = _incrlev;\
1084: nlnk++;\
1085: _lnkdata = _entry; /* next search starts from here */\
1086: } else { /* existing entry: update lnklvl */\
1087: if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
1088: }\
1089: }\
1090: }
1092: /*
1093: Copy data on the list into an array, then initialize the list
1094: Input Parameters:
1095: idx_start - starting index of the list
1096: lnk_max - max value of lnk indicating the end of the list
1097: nlnk - number of data on the list to be copied
1098: lnk - linked list
1099: lnklvl - level of lnk
1100: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
1101: output Parameters:
1102: indices - array that contains the copied data
1103: lnk - linked list that is cleaned and initialize
1104: lnklvl - level of lnk that is reinitialized
1105: bt - PetscBT (bitarray) with all bits set to false
1106: */
1107: #define PetscIncompleteLLClean(idx_start,lnk_max,nlnk,lnk,lnklvl,indices,indiceslvl,bt) 0;\
1108: {\
1109: PetscInt _j,_idx=idx_start;\
1110: for (_j=0; _j<nlnk; _j++){\
1111: _idx = lnk[_idx];\
1112: *(indices+_j) = _idx;\
1113: *(indiceslvl+_j) = lnklvl[_idx];\
1114: lnklvl[_idx] = -1;\
1115: PetscBTClear(bt,_idx);\
1116: }\
1117: lnk[idx_start] = lnk_max;\
1118: }
1119: /*
1120: Free memories used by the list
1121: */
1122: #define PetscIncompleteLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(&(bt)))
1124: /* -------------------------------------------------------------------------------------------------------*/
1125: #include <petscbt.h>
1128: /*
1129: Create and initialize a condensed linked list -
1130: same as PetscLLCreate(), but uses a scalable array 'lnk' with size of max number of entries, not O(N).
1131: Barry suggested this approach (Dec. 6, 2011):
1132: I've thought of an alternative way of representing a linked list that is efficient but doesn't have the O(N) scaling issue
1133: (it may be faster than the O(N) even sequentially due to less crazy memory access).
1135: Instead of having some like a 2 -> 4 -> 11 -> 22 list that uses slot 2 4 11 and 22 in a big array use a small array with two slots
1136: for each entry for example [ 2 1 | 4 3 | 22 -1 | 11 2] so the first number (of the pair) is the value while the second tells you where
1137: in the list the next entry is. Inserting a new link means just append another pair at the end. For example say we want to insert 13 into the
1138: list it would then become [2 1 | 4 3 | 22 -1 | 11 4 | 13 2 ] you just add a pair at the end and fix the point for the one that points to it.
1139: That is 11 use to point to the 2 slot, after the change 11 points to the 4th slot which has the value 13. Note that values are always next
1140: to each other so memory access is much better than using the big array.
1142: Example:
1143: nlnk_max=5, lnk_max=36:
1144: Initial list: [0, 0 | 36, 2 | 0, 0 | 0, 0 | 0, 0 | 0, 0 | 0, 0]
1145: here, head_node has index 2 with value lnk[2]=lnk_max=36,
1146: 0-th entry is used to store the number of entries in the list,
1147: The initial lnk represents head -> tail(marked by 36) with number of entries = lnk[0]=0.
1148:
1149: Now adding a sorted set {2,4}, the list becomes
1150: [2, 0 | 36, 4 |2, 6 | 4, 2 | 0, 0 | 0, 0 | 0, 0 ]
1151: represents head -> 2 -> 4 -> tail with number of entries = lnk[0]=2.
1153: Then adding a sorted set {0,3,35}, the list
1154: [5, 0 | 36, 8 | 2, 10 | 4, 12 | 0, 4 | 3, 6 | 35, 2 ]
1155: represents head -> 0 -> 2 -> 3 -> 4 -> 35 -> tail with number of entries = lnk[0]=5.
1156:
1157: Input Parameters:
1158: nlnk_max - max length of the list
1159: lnk_max - max value of the entries
1160: Output Parameters:
1161: lnk - list created and initialized
1162: bt - PetscBT (bitarray) with all bits set to false. Note: bt has size lnk_max, not nln_max!
1163: */
1164: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedCreate(PetscInt nlnk_max,PetscInt lnk_max,PetscInt **lnk,PetscBT *bt)
1165: {
1167: PetscInt *llnk;
1170: PetscMalloc(2*(nlnk_max+2)*sizeof(PetscInt),lnk);
1171: PetscBTCreate(lnk_max,bt);
1172: llnk = *lnk;
1173: llnk[0] = 0; /* number of entries on the list */
1174: llnk[2] = lnk_max; /* value in the head node */
1175: llnk[3] = 2; /* next for the head node */
1176: return(0);
1177: }
1181: /*
1182: Add a SORTED ascending index set into a sorted linked list. See PetscLLCondensedCreate() for detailed description.
1183: Input Parameters:
1184: nidx - number of input indices
1185: indices - sorted interger array
1186: lnk - condensed linked list(an integer array) that is created
1187: bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
1188: output Parameters:
1189: lnk - the sorted(increasing order) linked list containing previous and newly added non-redundate indices
1190: bt - updated PetscBT (bitarray)
1191: */
1192: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedAddSorted(PetscInt nidx,const PetscInt indices[],PetscInt lnk[],PetscBT bt)
1193: {
1194: PetscInt _k,_entry,_location,_next,_lnkdata,_nlnk,_newnode;
1197: _nlnk = lnk[0]; /* num of entries on the input lnk */
1198: _location = 2; /* head */
1199: for (_k=0; _k<nidx; _k++){
1200: _entry = indices[_k];
1201: if (!PetscBTLookupSet(bt,_entry)){ /* new entry */
1202: /* search for insertion location */
1203: do {
1204: _next = _location + 1; /* link from previous node to next node */
1205: _location = lnk[_next]; /* idx of next node */
1206: _lnkdata = lnk[_location];/* value of next node */
1207: } while (_entry > _lnkdata);
1208: /* insertion location is found, add entry into lnk */
1209: _newnode = 2*(_nlnk+2); /* index for this new node */
1210: lnk[_next] = _newnode; /* connect previous node to the new node */
1211: lnk[_newnode] = _entry; /* set value of the new node */
1212: lnk[_newnode+1] = _location; /* connect new node to next node */
1213: _location = _newnode; /* next search starts from the new node */
1214: _nlnk++;
1215: } \
1216: }\
1217: lnk[0] = _nlnk; /* number of entries in the list */
1218: return(0);
1219: }
1223: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedClean(PetscInt lnk_max,PetscInt nidx,PetscInt *indices,PetscInt lnk[],PetscBT bt)
1224: {
1226: PetscInt _k,_next,_nlnk;
1229: _next = lnk[3]; /* head node */
1230: _nlnk = lnk[0]; /* num of entries on the list */
1231: for (_k=0; _k<_nlnk; _k++){
1232: indices[_k] = lnk[_next];
1233: _next = lnk[_next + 1];
1234: PetscBTClear(bt,indices[_k]);
1235: }
1236: lnk[0] = 0; /* num of entries on the list */
1237: lnk[2] = lnk_max; /* initialize head node */
1238: lnk[3] = 2; /* head node */
1239: return(0);
1240: }
1244: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedView(PetscInt *lnk)
1245: {
1247: PetscInt k;
1250: PetscPrintf(PETSC_COMM_SELF,"LLCondensed of size %d, (val, next)\n",lnk[0]);
1251: for (k=2; k< lnk[0]+2; k++){
1252: PetscPrintf(PETSC_COMM_SELF," %D: (%D, %D)\n",2*k,lnk[2*k],lnk[2*k+1]);
1253: }
1254: return(0);
1255: }
1259: /*
1260: Free memories used by the list
1261: */
1262: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedDestroy(PetscInt *lnk,PetscBT bt)
1263: {
1267: PetscFree(lnk);
1268: PetscBTDestroy(&bt);
1269: return(0);
1270: }
1272: /* -------------------------------------------------------------------------------------------------------*/
1275: /*
1276: Same as PetscLLCondensedCreate(), but does not use O(lnk_max) bitarray
1277: */
1278: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedCreate_Scalable(PetscInt lnk_max,PetscInt **lnk)
1279: {
1281: PetscInt *llnk;
1284: PetscMalloc(2*(lnk_max+2)*sizeof(PetscInt),lnk);
1285: llnk = *lnk;
1286: llnk[0] = 0; /* number of entries on the list */
1287: llnk[2] = PETSC_MAX_INT; /* value in the head node */
1288: llnk[3] = 2; /* next for the head node */
1289: return(0);
1290: }
1294: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedAddSorted_Scalable(PetscInt nidx,const PetscInt indices[],PetscInt lnk[])
1295: {
1296: PetscInt _k,_entry,_location,_next,_lnkdata,_nlnk,_newnode;
1297: _nlnk = lnk[0]; /* num of entries on the input lnk */
1298: _location = 2; /* head */ \
1299: for (_k=0; _k<nidx; _k++){
1300: _entry = indices[_k];
1301: /* search for insertion location */
1302: do {
1303: _next = _location + 1; /* link from previous node to next node */
1304: _location = lnk[_next]; /* idx of next node */
1305: _lnkdata = lnk[_location];/* value of next node */
1306: } while (_entry > _lnkdata);
1307: if (_entry < _lnkdata) {
1308: /* insertion location is found, add entry into lnk */
1309: _newnode = 2*(_nlnk+2); /* index for this new node */
1310: lnk[_next] = _newnode; /* connect previous node to the new node */
1311: lnk[_newnode] = _entry; /* set value of the new node */
1312: lnk[_newnode+1] = _location; /* connect new node to next node */
1313: _location = _newnode; /* next search starts from the new node */
1314: _nlnk++;
1315: }
1316: }
1317: lnk[0] = _nlnk; /* number of entries in the list */
1318: return 0;
1319: }
1323: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedClean_Scalable(PetscInt nidx,PetscInt *indices,PetscInt *lnk)
1324: {
1325: PetscInt _k,_next,_nlnk;
1326: _next = lnk[3]; /* head node */
1327: _nlnk = lnk[0];
1328: for (_k=0; _k<_nlnk; _k++){
1329: indices[_k] = lnk[_next];
1330: _next = lnk[_next + 1];
1331: }
1332: lnk[0] = 0; /* num of entries on the list */
1333: lnk[3] = 2; /* head node */
1334: return 0;
1335: }
1339: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedDestroy_Scalable(PetscInt *lnk)
1340: {
1341: return PetscFree(lnk);
1342: }
1344: /* -------------------------------------------------------------------------------------------------------*/
1345: /*
1346: lnk[0] number of links
1347: lnk[1] number of entries
1348: lnk[3n] value
1349: lnk[3n+1] len
1350: lnk[3n+2] link to next value
1352: The next three are always the first link
1354: lnk[3] PETSC_MIN_INT+1
1355: lnk[4] 1
1356: lnk[5] link to first real entry
1358: The next three are always the last link
1360: lnk[6] PETSC_MAX_INT - 1
1361: lnk[7] 1
1362: lnk[8] next valid link (this is the same as lnk[0] but without the decreases)
1363: */
1367: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedCreate_fast(PetscInt lnk_max,PetscInt **lnk)
1368: {
1370: PetscInt *llnk;
1373: PetscMalloc(3*(lnk_max+3)*sizeof(PetscInt),lnk);
1374: llnk = *lnk;
1375: llnk[0] = 0; /* nlnk: number of entries on the list */
1376: llnk[1] = 0; /* number of integer entries represented in list */
1377: llnk[3] = PETSC_MIN_INT+1; /* value in the first node */
1378: llnk[4] = 1; /* count for the first node */
1379: llnk[5] = 6; /* next for the first node */
1380: llnk[6] = PETSC_MAX_INT-1; /* value in the last node */
1381: llnk[7] = 1; /* count for the last node */
1382: llnk[8] = 0; /* next valid node to be used */
1383: return(0);
1384: }
1386: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedAddSorted_fast(PetscInt nidx,const PetscInt indices[],PetscInt lnk[])
1387: {
1388: PetscInt k,entry,prev,next;
1389: prev = 3; /* first value */
1390: next = lnk[prev+2];
1391: for (k=0; k<nidx; k++){
1392: entry = indices[k];
1393: /* search for insertion location */
1394: while (entry >= lnk[next]) {
1395: prev = next;
1396: next = lnk[next+2];
1397: }
1398: /* entry is in range of previous list */
1399: if (entry < lnk[prev]+lnk[prev+1]) continue;
1400: lnk[1]++;
1401: /* entry is right after previous list */
1402: if (entry == lnk[prev]+lnk[prev+1]) {
1403: lnk[prev+1]++;
1404: if (lnk[next] == entry+1) { /* combine two contiquous strings */
1405: lnk[prev+1] += lnk[next+1];
1406: lnk[prev+2] = lnk[next+2];
1407: next = lnk[next+2];
1408: lnk[0]--;
1409: }
1410: continue;
1411: }
1412: /* entry is right before next list */
1413: if (entry == lnk[next]-1) {
1414: lnk[next]--;
1415: lnk[next+1]++;
1416: prev = next;
1417: next = lnk[prev+2];
1418: continue;
1419: }
1420: /* add entry into lnk */
1421: lnk[prev+2] = 3*((lnk[8]++)+3); /* connect previous node to the new node */
1422: prev = lnk[prev+2];
1423: lnk[prev] = entry; /* set value of the new node */
1424: lnk[prev+1] = 1; /* number of values in contiquous string is one to start */
1425: lnk[prev+2] = next; /* connect new node to next node */
1426: lnk[0]++;
1427: }
1428: return 0;
1429: }
1431: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedClean_fast(PetscInt nidx,PetscInt *indices,PetscInt *lnk)
1432: {
1433: PetscInt _k,_next,_nlnk,cnt,j;
1434: _next = lnk[5]; /* first node */
1435: _nlnk = lnk[0];
1436: cnt = 0;
1437: for (_k=0; _k<_nlnk; _k++){
1438: for (j=0; j<lnk[_next+1]; j++) {
1439: indices[cnt++] = lnk[_next] + j;
1440: }
1441: _next = lnk[_next + 2];
1442: }
1443: lnk[0] = 0; /* nlnk: number of links */
1444: lnk[1] = 0; /* number of integer entries represented in list */
1445: lnk[3] = PETSC_MIN_INT+1; /* value in the first node */
1446: lnk[4] = 1; /* count for the first node */
1447: lnk[5] = 6; /* next for the first node */
1448: lnk[6] = PETSC_MAX_INT-1; /* value in the last node */
1449: lnk[7] = 1; /* count for the last node */
1450: lnk[8] = 0; /* next valid location to make link */
1451: return 0;
1452: }
1454: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedView_fast(PetscInt *lnk)
1455: {
1456: PetscInt k,next,nlnk;
1457: next = lnk[5]; /* first node */
1458: nlnk = lnk[0];
1459: for (k=0; k<nlnk; k++){
1460: #if 0 /* Debugging code */
1461: printf("%d value %d len %d next %d\n",next,lnk[next],lnk[next+1],lnk[next+2]);
1462: #endif
1463: next = lnk[next + 2];
1464: }
1465: return 0;
1466: }
1468: PETSC_STATIC_INLINE PetscErrorCode PetscLLCondensedDestroy_fast(PetscInt *lnk)
1469: {
1470: return PetscFree(lnk);
1471: }
1473: PETSC_EXTERN PetscLogEvent MAT_Mult, MAT_MultMatrixFree, MAT_Mults, MAT_MultConstrained, MAT_MultAdd, MAT_MultTranspose;
1474: PETSC_EXTERN PetscLogEvent MAT_MultTransposeConstrained, MAT_MultTransposeAdd, MAT_Solve, MAT_Solves, MAT_SolveAdd, MAT_SolveTranspose;
1475: PETSC_EXTERN PetscLogEvent MAT_SolveTransposeAdd, MAT_SOR, MAT_ForwardSolve, MAT_BackwardSolve, MAT_LUFactor, MAT_LUFactorSymbolic;
1476: PETSC_EXTERN PetscLogEvent MAT_LUFactorNumeric, MAT_CholeskyFactor, MAT_CholeskyFactorSymbolic, MAT_CholeskyFactorNumeric, MAT_ILUFactor;
1477: PETSC_EXTERN PetscLogEvent MAT_ILUFactorSymbolic, MAT_ICCFactorSymbolic, MAT_Copy, MAT_Convert, MAT_Scale, MAT_AssemblyBegin;
1478: PETSC_EXTERN PetscLogEvent MAT_AssemblyEnd, MAT_SetValues, MAT_GetValues, MAT_GetRow, MAT_GetRowIJ, MAT_GetSubMatrices, MAT_GetColoring, MAT_GetOrdering, MAT_GetRedundantMatrix;
1479: PETSC_EXTERN PetscLogEvent MAT_IncreaseOverlap, MAT_Partitioning, MAT_Coarsen, MAT_ZeroEntries, MAT_Load, MAT_View, MAT_AXPY, MAT_FDColoringCreate, MAT_TransposeColoringCreate;
1480: PETSC_EXTERN PetscLogEvent MAT_FDColoringApply, MAT_Transpose, MAT_FDColoringFunction;
1481: PETSC_EXTERN PetscLogEvent MAT_MatMult, MAT_MatSolve,MAT_MatMultSymbolic, MAT_MatMultNumeric,MAT_Getlocalmatcondensed,MAT_GetBrowsOfAcols,MAT_GetBrowsOfAocols;
1482: PETSC_EXTERN PetscLogEvent MAT_PtAP, MAT_PtAPSymbolic, MAT_PtAPNumeric,MAT_Seqstompinum,MAT_Seqstompisym,MAT_Seqstompi,MAT_Getlocalmat;
1483: PETSC_EXTERN PetscLogEvent MAT_RARt, MAT_RARtSymbolic, MAT_RARtNumeric;
1484: PETSC_EXTERN PetscLogEvent MAT_MatTransposeMult, MAT_MatTransposeMultSymbolic, MAT_MatTransposeMultNumeric;
1485: PETSC_EXTERN PetscLogEvent MAT_TransposeMatMult, MAT_TransposeMatMultSymbolic, MAT_TransposeMatMultNumeric;
1486: PETSC_EXTERN PetscLogEvent MAT_Applypapt, MAT_Applypapt_symbolic, MAT_Applypapt_numeric;
1487: PETSC_EXTERN PetscLogEvent MAT_Getsymtranspose, MAT_Transpose_SeqAIJ, MAT_Getsymtransreduced,MAT_GetSequentialNonzeroStructure;
1489: PETSC_EXTERN PetscLogEvent MATMFFD_Mult;
1490: PETSC_EXTERN PetscLogEvent MAT_GetMultiProcBlock;
1491: PETSC_EXTERN PetscLogEvent MAT_CUSPCopyToGPU, MAT_SetValuesBatch, MAT_SetValuesBatchI, MAT_SetValuesBatchII, MAT_SetValuesBatchIII, MAT_SetValuesBatchIV;
1492: PETSC_EXTERN PetscLogEvent MAT_Merge;
1494: #endif