Actual source code: matimpl.h

  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 (*setuppreallocation)(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 (*setblocksize)(Mat,PetscInt);
 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 (*usescaledform)(Mat,PetscBool );
 95:   /*64*/
 96:   PetscErrorCode (*scalesystem)(Mat,Vec,Vec);
 97:   PetscErrorCode (*unscalesystem)(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 (*matmulttranspose)(Mat,Mat,MatReuse,PetscReal,Mat*);
134:   PetscErrorCode (*matmulttransposesymbolic)(Mat,Mat,PetscReal,Mat*);
135:   PetscErrorCode (*matmulttransposenumeric)(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 (*setsizes)(Mat,PetscInt,PetscInt,PetscInt,PetscInt);
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,Mat*);
167:   /*124*/
168:   PetscErrorCode (*findnonzerorows)(Mat,IS*);
169:   PetscErrorCode (*getcolumnnorms)(Mat,NormType,PetscReal*);
170:   PetscErrorCode (*invertblockdiagonal)(Mat,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: };
176: /*
177:     If you add MatOps entries above also add them to the MATOP enum
178:     in include/petscmat.h and include/finclude/petscmat.h
179: */

181: typedef struct _p_MatBaseName* MatBaseName;
182: struct _p_MatBaseName {
183:   char        *bname,*sname,*mname;
184:   MatBaseName next;
185: };


189: /*
190:    Utility private matrix routines
191: */


202: /* 
203:   The stash is used to temporarily store inserted matrix values that 
204:   belong to another processor. During the assembly phase the stashed 
205:   values are moved to the correct processor and 
206: */

208: typedef struct _MatStashSpace *PetscMatStashSpace;

210: struct _MatStashSpace {
211:   PetscMatStashSpace next;
212:   PetscScalar        *space_head,*val;
213:   PetscInt           *idx,*idy;
214:   PetscInt           total_space_size;
215:   PetscInt           local_used;
216:   PetscInt           local_remaining;
217: };


223: typedef struct {
224:   PetscInt      nmax;                   /* maximum stash size */
225:   PetscInt      umax;                   /* user specified max-size */
226:   PetscInt      oldnmax;                /* the nmax value used previously */
227:   PetscInt      n;                      /* stash size */
228:   PetscInt      bs;                     /* block size of the stash */
229:   PetscInt      reallocs;               /* preserve the no of mallocs invoked */
230:   PetscMatStashSpace space_head,space;  /* linked list to hold stashed global row/column numbers and matrix values */
231:   /* The following variables are used for communication */
232:   MPI_Comm      comm;
233:   PetscMPIInt   size,rank;
234:   PetscMPIInt   tag1,tag2;
235:   MPI_Request   *send_waits;            /* array of send requests */
236:   MPI_Request   *recv_waits;            /* array of receive requests */
237:   MPI_Status    *send_status;           /* array of send status */
238:   PetscInt      nsends,nrecvs;          /* numbers of sends and receives */
239:   PetscScalar   *svalues;               /* sending data */
240:   PetscInt      *sindices;
241:   PetscScalar   **rvalues;              /* receiving data (values) */
242:   PetscInt      **rindices;             /* receiving data (indices) */
243:   PetscInt      nprocessed;             /* number of messages already processed */
244:   PetscMPIInt   *flg_v;                 /* indicates what messages have arrived so far and from whom */
245:   PetscBool     reproduce;
246:   PetscInt      reproduce_count;
247: } MatStash;


261: typedef struct {
262:   PetscInt   dim;
263:   PetscInt   dims[4];
264:   PetscInt   starts[4];
265:   PetscBool  noc;        /* this is a single component problem, hence user will not set MatStencil.c */
266: } MatStencilInfo;

268: /* Info about using compressed row format */
269: typedef struct {
270:   PetscBool  check;                         /* indicates that at MatAssembly() it should check if compressed rows will be efficient */
271:   PetscBool  use;                           /* indicates compressed rows have been checked and will be used */
272:   PetscInt   nrows;                         /* number of non-zero rows */
273:   PetscInt   *i;                            /* compressed row pointer  */
274:   PetscInt   *rindex;                       /* compressed row index               */
275: } Mat_CompressedRow;

278: struct _p_Mat {
279:   PETSCHEADER(struct _MatOps);
280:   PetscLayout            rmap,cmap;
281:   void                   *data;            /* implementation-specific data */
282:   MatFactorType          factortype;       /* MAT_FACTOR_LU, ILU, CHOLESKY or ICC */
283:   PetscBool              assembled;        /* is the matrix assembled? */
284:   PetscBool              was_assembled;    /* new values inserted into assembled mat */
285:   PetscInt               num_ass;          /* number of times matrix has been assembled */
286:   PetscBool              same_nonzero;     /* matrix has same nonzero pattern as previous */
287:   MatInfo                info;             /* matrix information */
288:   InsertMode             insertmode;       /* have values been inserted in matrix or added? */
289:   MatStash               stash,bstash;     /* used for assembling off-proc mat emements */
290:   MatNullSpace           nullsp;
291:   PetscBool              preallocated;
292:   MatStencilInfo         stencil;          /* information for structured grid */
293:   PetscBool              symmetric,hermitian,structurally_symmetric,spd;
294:   PetscBool              symmetric_set,hermitian_set,structurally_symmetric_set,spd_set; /* if true, then corresponding flag is correct*/
295:   PetscBool              symmetric_eternal;
296:   PetscBool              nooffprocentries,nooffproczerorows;
297: #if defined(PETSC_HAVE_CUSP)
298:   PetscCUSPFlag          valid_GPU_matrix; /* flag pointing to the matrix on the gpu*/
299: #endif
300:   void                   *spptr;          /* pointer for special library like SuperLU */
301:   MatSolverPackage       solvertype;
302:   };

304: #define MatPreallocated(A)  ((!(A)->preallocated) ? MatSetUpPreallocation(A) : 0)
307: /*
308:     Object for partitioning graphs
309: */

311: typedef struct _MatPartitioningOps *MatPartitioningOps;
312: struct _MatPartitioningOps {
313:   PetscErrorCode (*apply)(MatPartitioning,IS*);
314:   PetscErrorCode (*setfromoptions)(MatPartitioning);
315:   PetscErrorCode (*destroy)(MatPartitioning);
316:   PetscErrorCode (*view)(MatPartitioning,PetscViewer);
317: };

319: struct _p_MatPartitioning {
320:   PETSCHEADER(struct _MatPartitioningOps);
321:   Mat         adj;
322:   PetscInt    *vertex_weights;
323:   PetscReal   *part_weights;
324:   PetscInt    n;                                 /* number of partitions */
325:   void        *data;
326:   PetscInt    setupcalled;
327: };

329: /*
330:     MatFDColoring is used to compute Jacobian matrices efficiently
331:   via coloring. The data structure is explained below in an example.

333:    Color =   0    1     0    2   |   2      3       0 
334:    ---------------------------------------------------
335:             00   01              |          05
336:             10   11              |   14     15               Processor  0
337:                        22    23  |          25
338:                        32    33  | 
339:    ===================================================
340:                                  |   44     45     46
341:             50                   |          55               Processor 1
342:                                  |   64            66
343:    ---------------------------------------------------

345:     ncolors = 4;

347:     ncolumns      = {2,1,1,0}
348:     columns       = {{0,2},{1},{3},{}}
349:     nrows         = {4,2,3,3}
350:     rows          = {{0,1,2,3},{0,1},{1,2,3},{0,1,2}}
351:     columnsforrow = {{0,0,2,2},{1,1},{4,3,3},{5,5,5}}
352:     vscaleforrow  = {{,,,},{,},{,,},{,,}}
353:     vwscale       = {dx(0),dx(1),dx(2),dx(3)}               MPI Vec
354:     vscale        = {dx(0),dx(1),dx(2),dx(3),dx(4),dx(5)}   Seq Vec

356:     ncolumns      = {1,0,1,1}
357:     columns       = {{6},{},{4},{5}}
358:     nrows         = {3,0,2,2}
359:     rows          = {{0,1,2},{},{1,2},{1,2}}
360:     columnsforrow = {{6,0,6},{},{4,4},{5,5}}
361:     vscaleforrow =  {{,,},{},{,},{,}}
362:     vwscale       = {dx(4),dx(5),dx(6)}              MPI Vec
363:     vscale        = {dx(0),dx(4),dx(5),dx(6)}        Seq Vec

365:     See the routine MatFDColoringApply() for how this data is used
366:     to compute the Jacobian.

368: */

370: struct  _p_MatFDColoring{
371:   PETSCHEADER(int);
372:   PetscInt       M,N,m;            /* total rows, columns; local rows */
373:   PetscInt       rstart;           /* first row owned by local processor */
374:   PetscInt       ncolors;          /* number of colors */
375:   PetscInt       *ncolumns;        /* number of local columns for a color */
376:   PetscInt       **columns;        /* lists the local columns of each color (using global column numbering) */
377:   PetscInt       *nrows;           /* number of local rows for each color */
378:   PetscInt       **rows;           /* lists the local rows for each color (using the local row numbering) */
379:   PetscInt       **columnsforrow;  /* lists the corresponding columns for those rows (using the global column) */
380:   PetscReal      error_rel;        /* square root of relative error in computing function */
381:   PetscReal      umin;             /* minimum allowable u'dx value */
382:   Vec            w1,w2,w3;         /* work vectors used in computing Jacobian */
383:   PetscErrorCode (*f)(void);       /* function that defines Jacobian */
384:   void           *fctx;            /* optional user-defined context for use by the function f */
385:   PetscInt       **vscaleforrow;   /* location in vscale for each columnsforrow[] entry */
386:   Vec            vscale;           /* holds FD scaling, i.e. 1/dx for each perturbed column */
387:   Vec            F;                /* current value of user provided function; can set with MatFDColoringSetF() */
388:   PetscInt       currentcolor;     /* color for which function evaluation is being done now */
389:   const char     *htype;            /* "wp" or "ds" */
390:   ISColoringType ctype;            /* IS_COLORING_GLOBAL or IS_COLORING_GHOSTED */

392:   void           *ftn_func_pointer,*ftn_func_cntx; /* serve the same purpose as *fortran_func_pointers in PETSc objects */
393: };

395: /*
396:    Null space context for preconditioner/operators
397: */
398: struct _p_MatNullSpace {
399:   PETSCHEADER(int);
400:   PetscBool      has_cnst;
401:   PetscInt       n;
402:   Vec*           vecs;
403:   PetscScalar*   alpha;                 /* for projections */
404:   Vec            vec;                   /* for out of place removals */
405:   PetscErrorCode (*remove)(MatNullSpace,Vec,void*);  /* for user provided removal function */
406:   void*          rmctx;                 /* context for remove() function */
407: };

409: /* 
410:    Checking zero pivot for LU, ILU preconditioners.
411: */
412: typedef struct {
413:   PetscInt       nshift,nshift_max;
414:   PetscReal      shift_amount,shift_lo,shift_hi,shift_top,shift_fraction;
415:   PetscBool      newshift;
416:   PetscReal      rs;  /* active row sum of abs(offdiagonals) */
417:   PetscScalar    pv;  /* pivot of the active row */
418: } FactorShiftCtx;


424: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_nz(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
425: {
426:   PetscReal _rs   = sctx->rs;
427:   PetscReal _zero = info->zeropivot*_rs;

430:   if (PetscAbsScalar(sctx->pv) <= _zero){
431:     /* force |diag| > zeropivot*rs */
432:     if (!sctx->nshift) sctx->shift_amount = info->shiftamount;
433:     else sctx->shift_amount *= 2.0;
434:     sctx->newshift = PETSC_TRUE;
435:     (sctx->nshift)++;
436:   } else {
437:     sctx->newshift = PETSC_FALSE;
438:   }
439:   return(0);
440: }

444: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_pd(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
445: {
446:   PetscReal _rs   = sctx->rs;
447:   PetscReal _zero = info->zeropivot*_rs;

450:   if (PetscRealPart(sctx->pv) <= _zero){
451:     /* force matfactor to be diagonally dominant */
452:     if (sctx->nshift == sctx->nshift_max) {
453:       sctx->shift_fraction = sctx->shift_hi;
454:     } else {
455:       sctx->shift_lo = sctx->shift_fraction;
456:       sctx->shift_fraction = (sctx->shift_hi+sctx->shift_lo)/2.;
457:     }
458:     sctx->shift_amount = sctx->shift_fraction * sctx->shift_top;
459:     sctx->nshift++;
460:     sctx->newshift = PETSC_TRUE;
461:   } else {
462:     sctx->newshift = PETSC_FALSE;
463:   }
464:   return(0);
465: }

469: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_inblocks(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
470: {
471:   PetscReal _zero = info->zeropivot;

474:   if (PetscAbsScalar(sctx->pv) <= _zero){
475:     sctx->pv          += info->shiftamount;
476:     sctx->shift_amount = 0.0;
477:     sctx->nshift++;
478:   }
479:   sctx->newshift = PETSC_FALSE;
480:   return(0);
481: }

485: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck_none(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
486: {
487:   PetscReal _zero = info->zeropivot;

490:   sctx->newshift = PETSC_FALSE;
491:   if (PetscAbsScalar(sctx->pv) <= _zero) {
493:     PetscBool      flg = PETSC_FALSE;
494: 
495:     PetscOptionsGetBool(PETSC_NULL,"-mat_dump",&flg,PETSC_NULL);
496:     if (flg) {
497:       MatView(mat,PETSC_VIEWER_BINARY_(((PetscObject)mat)->comm));
498:     }
499:     SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_MAT_LU_ZRPVT,"Zero pivot row %D value %G tolerance %G",row,PetscAbsScalar(sctx->pv),_zero);
500:   }
501:   return(0);
502: }

506: PETSC_STATIC_INLINE PetscErrorCode MatPivotCheck(Mat mat,const MatFactorInfo *info,FactorShiftCtx *sctx,PetscInt row)
507: {

511:   if (info->shifttype == (PetscReal) MAT_SHIFT_NONZERO){
512:     MatPivotCheck_nz(mat,info,sctx,row);
513:   } else if (info->shifttype == (PetscReal) MAT_SHIFT_POSITIVE_DEFINITE){
514:     MatPivotCheck_pd(mat,info,sctx,row);
515:   } else if (info->shifttype == (PetscReal) MAT_SHIFT_INBLOCKS){
516:     MatPivotCheck_inblocks(mat,info,sctx,row);
517:   } else {
518:     MatPivotCheck_none(mat,info,sctx,row);
519:   }
520:   return(0);
521: }

523: /* 
524:   Create and initialize a linked list 
525:   Input Parameters:
526:     idx_start - starting index of the list
527:     lnk_max   - max value of lnk indicating the end of the list
528:     nlnk      - max length of the list
529:   Output Parameters:
530:     lnk       - list initialized
531:     bt        - PetscBT (bitarray) with all bits set to false
532: */
533: #define PetscLLCreate(idx_start,lnk_max,nlnk,lnk,bt) \
534:   (PetscMalloc(nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,0))

536: /*
537:   Add an index set into a sorted linked list
538:   Input Parameters:
539:     nidx      - number of input indices
540:     indices   - interger array
541:     idx_start - starting index of the list
542:     lnk       - linked list(an integer array) that is created
543:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
544:   output Parameters:
545:     nlnk      - number of newly added indices
546:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
547:     bt        - updated PetscBT (bitarray) 
548: */
549: #define PetscLLAdd(nidx,indices,idx_start,nlnk,lnk,bt) 0;\
550: {\
551:   PetscInt _k,_entry,_location,_lnkdata;\
552:   nlnk     = 0;\
553:   _lnkdata = idx_start;\
554:   for (_k=0; _k<nidx; _k++){\
555:     _entry = indices[_k];\
556:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
557:       /* search for insertion location */\
558:       /* start from the beginning if _entry < previous _entry */\
559:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
560:       do {\
561:         _location = _lnkdata;\
562:         _lnkdata  = lnk[_location];\
563:       } while (_entry > _lnkdata);\
564:       /* insertion location is found, add entry into lnk */\
565:       lnk[_location] = _entry;\
566:       lnk[_entry]    = _lnkdata;\
567:       nlnk++;\
568:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
569:     }\
570:   }\
571: }

573: /*
574:   Add a permuted index set into a sorted linked list
575:   Input Parameters:
576:     nidx      - number of input indices
577:     indices   - interger array
578:     perm      - permutation of indices
579:     idx_start - starting index of the list
580:     lnk       - linked list(an integer array) that is created
581:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
582:   output Parameters:
583:     nlnk      - number of newly added indices
584:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
585:     bt        - updated PetscBT (bitarray) 
586: */
587: #define PetscLLAddPerm(nidx,indices,perm,idx_start,nlnk,lnk,bt) 0;\
588: {\
589:   PetscInt _k,_entry,_location,_lnkdata;\
590:   nlnk     = 0;\
591:   _lnkdata = idx_start;\
592:   for (_k=0; _k<nidx; _k++){\
593:     _entry = perm[indices[_k]];\
594:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
595:       /* search for insertion location */\
596:       /* start from the beginning if _entry < previous _entry */\
597:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
598:       do {\
599:         _location = _lnkdata;\
600:         _lnkdata  = lnk[_location];\
601:       } while (_entry > _lnkdata);\
602:       /* insertion location is found, add entry into lnk */\
603:       lnk[_location] = _entry;\
604:       lnk[_entry]    = _lnkdata;\
605:       nlnk++;\
606:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
607:     }\
608:   }\
609: }

611: /*
612:   Add a SORTED index set into a sorted linked list
613:   Input Parameters:
614:     nidx      - number of input indices
615:     indices   - sorted interger array 
616:     idx_start - starting index of the list
617:     lnk       - linked list(an integer array) that is created
618:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
619:   output Parameters:
620:     nlnk      - number of newly added indices
621:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
622:     bt        - updated PetscBT (bitarray) 
623: */
624: #define PetscLLAddSorted(nidx,indices,idx_start,nlnk,lnk,bt) 0;\
625: {\
626:   PetscInt _k,_entry,_location,_lnkdata;\
627:   nlnk      = 0;\
628:   _lnkdata  = idx_start;\
629:   for (_k=0; _k<nidx; _k++){\
630:     _entry = indices[_k];\
631:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
632:       /* search for insertion location */\
633:       do {\
634:         _location = _lnkdata;\
635:         _lnkdata  = lnk[_location];\
636:       } while (_entry > _lnkdata);\
637:       /* insertion location is found, add entry into lnk */\
638:       lnk[_location] = _entry;\
639:       lnk[_entry]    = _lnkdata;\
640:       nlnk++;\
641:       _lnkdata = _entry; /* next search starts from here */\
642:     }\
643:   }\
644: }

646: /*
647:   Add a SORTED index set into a sorted linked list used for LUFactorSymbolic()
648:   Same as PetscLLAddSorted() with an additional operation:
649:        count the number of input indices that are no larger than 'diag' 
650:   Input Parameters:
651:     indices   - sorted interger array 
652:     idx_start - starting index of the list, index of pivot row
653:     lnk       - linked list(an integer array) that is created
654:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
655:     diag      - index of the active row in LUFactorSymbolic
656:     nzbd      - number of input indices with indices <= idx_start
657:     im        - im[idx_start] is initialized as num of nonzero entries in row=idx_start
658:   output Parameters:
659:     nlnk      - number of newly added indices
660:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
661:     bt        - updated PetscBT (bitarray) 
662:     im        - im[idx_start]: unchanged if diag is not an entry 
663:                              : num of entries with indices <= diag if diag is an entry
664: */
665: #define PetscLLAddSortedLU(indices,idx_start,nlnk,lnk,bt,diag,nzbd,im) 0;\
666: {\
667:   PetscInt _k,_entry,_location,_lnkdata,_nidx;\
668:   nlnk     = 0;\
669:   _lnkdata = idx_start;\
670:   _nidx = im[idx_start] - nzbd; /* num of entries with idx_start < index <= diag */\
671:   for (_k=0; _k<_nidx; _k++){\
672:     _entry = indices[_k];\
673:     nzbd++;\
674:     if ( _entry== diag) im[idx_start] = nzbd;\
675:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
676:       /* search for insertion location */\
677:       do {\
678:         _location = _lnkdata;\
679:         _lnkdata  = lnk[_location];\
680:       } while (_entry > _lnkdata);\
681:       /* insertion location is found, add entry into lnk */\
682:       lnk[_location] = _entry;\
683:       lnk[_entry]    = _lnkdata;\
684:       nlnk++;\
685:       _lnkdata = _entry; /* next search starts from here */\
686:     }\
687:   }\
688: }

690: /*
691:   Copy data on the list into an array, then initialize the list 
692:   Input Parameters:
693:     idx_start - starting index of the list 
694:     lnk_max   - max value of lnk indicating the end of the list 
695:     nlnk      - number of data on the list to be copied
696:     lnk       - linked list
697:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
698:   output Parameters:
699:     indices   - array that contains the copied data
700:     lnk       - linked list that is cleaned and initialize
701:     bt        - PetscBT (bitarray) with all bits set to false
702: */
703: #define PetscLLClean(idx_start,lnk_max,nlnk,lnk,indices,bt) 0;\
704: {\
705:   PetscInt _j,_idx=idx_start;\
706:   for (_j=0; _j<nlnk; _j++){\
707:     _idx = lnk[_idx];\
708:     *(indices+_j) = _idx;\
709:     PetscBTClear(bt,_idx);\
710:   }\
711:   lnk[idx_start] = lnk_max;\
712: }
713: /*
714:   Free memories used by the list
715: */
716: #define PetscLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(bt))

718: /* Routines below are used for incomplete matrix factorization */
719: /* 
720:   Create and initialize a linked list and its levels
721:   Input Parameters:
722:     idx_start - starting index of the list
723:     lnk_max   - max value of lnk indicating the end of the list
724:     nlnk      - max length of the list
725:   Output Parameters:
726:     lnk       - list initialized
727:     lnk_lvl   - array of size nlnk for storing levels of lnk
728:     bt        - PetscBT (bitarray) with all bits set to false
729: */
730: #define PetscIncompleteLLCreate(idx_start,lnk_max,nlnk,lnk,lnk_lvl,bt)\
731:   (PetscMalloc(2*nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,lnk_lvl = lnk + nlnk,0))

733: /*
734:   Initialize a sorted linked list used for ILU and ICC
735:   Input Parameters:
736:     nidx      - number of input idx
737:     idx       - interger array used for storing column indices
738:     idx_start - starting index of the list
739:     perm      - indices of an IS
740:     lnk       - linked list(an integer array) that is created
741:     lnklvl    - levels of lnk
742:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
743:   output Parameters:
744:     nlnk     - number of newly added idx
745:     lnk      - the sorted(increasing order) linked list containing new and non-redundate entries from idx
746:     lnklvl   - levels of lnk
747:     bt       - updated PetscBT (bitarray) 
748: */
749: #define PetscIncompleteLLInit(nidx,idx,idx_start,perm,nlnk,lnk,lnklvl,bt) 0;\
750: {\
751:   PetscInt _k,_entry,_location,_lnkdata;\
752:   nlnk     = 0;\
753:   _lnkdata = idx_start;\
754:   for (_k=0; _k<nidx; _k++){\
755:     _entry = perm[idx[_k]];\
756:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
757:       /* search for insertion location */\
758:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
759:       do {\
760:         _location = _lnkdata;\
761:         _lnkdata  = lnk[_location];\
762:       } while (_entry > _lnkdata);\
763:       /* insertion location is found, add entry into lnk */\
764:       lnk[_location]  = _entry;\
765:       lnk[_entry]     = _lnkdata;\
766:       lnklvl[_entry] = 0;\
767:       nlnk++;\
768:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
769:     }\
770:   }\
771: }

773: /*
774:   Add a SORTED index set into a sorted linked list for ILU
775:   Input Parameters:
776:     nidx      - number of input indices
777:     idx       - sorted interger array used for storing column indices
778:     level     - level of fill, e.g., ICC(level)
779:     idxlvl    - level of idx 
780:     idx_start - starting index of the list
781:     lnk       - linked list(an integer array) that is created
782:     lnklvl    - levels of lnk
783:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
784:     prow      - the row number of idx
785:   output Parameters:
786:     nlnk     - number of newly added idx
787:     lnk      - the sorted(increasing order) linked list containing new and non-redundate entries from idx
788:     lnklvl   - levels of lnk
789:     bt       - updated PetscBT (bitarray) 

791:   Note: the level of factor(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(i,prow)+lvl(prow,j)+1)
792:         where idx = non-zero columns of U(prow,prow+1:n-1), prow<i
793: */
794: #define PetscILULLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,lnklvl_prow) 0;\
795: {\
796:   PetscInt _k,_entry,_location,_lnkdata,_incrlev,_lnklvl_prow=lnklvl[prow];\
797:   nlnk     = 0;\
798:   _lnkdata = idx_start;\
799:   for (_k=0; _k<nidx; _k++){\
800:     _incrlev = idxlvl[_k] + _lnklvl_prow + 1;\
801:     if (_incrlev > level) continue;\
802:     _entry = idx[_k];\
803:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
804:       /* search for insertion location */\
805:       do {\
806:         _location = _lnkdata;\
807:         _lnkdata  = lnk[_location];\
808:       } while (_entry > _lnkdata);\
809:       /* insertion location is found, add entry into lnk */\
810:       lnk[_location]  = _entry;\
811:       lnk[_entry]     = _lnkdata;\
812:       lnklvl[_entry] = _incrlev;\
813:       nlnk++;\
814:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
815:     } else { /* existing entry: update lnklvl */\
816:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
817:     }\
818:   }\
819: }

821: /*
822:   Add a index set into a sorted linked list
823:   Input Parameters:
824:     nidx      - number of input idx
825:     idx   - interger array used for storing column indices
826:     level     - level of fill, e.g., ICC(level)
827:     idxlvl - level of idx 
828:     idx_start - starting index of the list
829:     lnk       - linked list(an integer array) that is created
830:     lnklvl   - levels of lnk
831:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
832:   output Parameters:
833:     nlnk      - number of newly added idx
834:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from idx
835:     lnklvl   - levels of lnk
836:     bt        - updated PetscBT (bitarray) 
837: */
838: #define PetscIncompleteLLAdd(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\
839: {\
840:   PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
841:   nlnk     = 0;\
842:   _lnkdata = idx_start;\
843:   for (_k=0; _k<nidx; _k++){\
844:     _incrlev = idxlvl[_k] + 1;\
845:     if (_incrlev > level) continue;\
846:     _entry = idx[_k];\
847:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
848:       /* search for insertion location */\
849:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
850:       do {\
851:         _location = _lnkdata;\
852:         _lnkdata  = lnk[_location];\
853:       } while (_entry > _lnkdata);\
854:       /* insertion location is found, add entry into lnk */\
855:       lnk[_location]  = _entry;\
856:       lnk[_entry]     = _lnkdata;\
857:       lnklvl[_entry] = _incrlev;\
858:       nlnk++;\
859:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
860:     } else { /* existing entry: update lnklvl */\
861:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
862:     }\
863:   }\
864: }

866: /*
867:   Add a SORTED index set into a sorted linked list
868:   Input Parameters:
869:     nidx      - number of input indices
870:     idx   - sorted interger array used for storing column indices
871:     level     - level of fill, e.g., ICC(level)
872:     idxlvl - level of idx 
873:     idx_start - starting index of the list
874:     lnk       - linked list(an integer array) that is created
875:     lnklvl    - levels of lnk
876:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
877:   output Parameters:
878:     nlnk      - number of newly added idx
879:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from idx
880:     lnklvl    - levels of lnk
881:     bt        - updated PetscBT (bitarray) 
882: */
883: #define PetscIncompleteLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\
884: {\
885:   PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
886:   nlnk = 0;\
887:   _lnkdata = idx_start;\
888:   for (_k=0; _k<nidx; _k++){\
889:     _incrlev = idxlvl[_k] + 1;\
890:     if (_incrlev > level) continue;\
891:     _entry = idx[_k];\
892:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
893:       /* search for insertion location */\
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] = _incrlev;\
902:       nlnk++;\
903:       _lnkdata = _entry; /* next search starts from here */\
904:     } else { /* existing entry: update lnklvl */\
905:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
906:     }\
907:   }\
908: }

910: /*
911:   Add a SORTED index set into a sorted linked list for ICC
912:   Input Parameters:
913:     nidx      - number of input indices
914:     idx       - sorted interger array used for storing column indices
915:     level     - level of fill, e.g., ICC(level)
916:     idxlvl    - level of idx 
917:     idx_start - starting index of the list
918:     lnk       - linked list(an integer array) that is created
919:     lnklvl    - levels of lnk
920:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
921:     idxlvl_prow - idxlvl[prow], where prow is the row number of the idx
922:   output Parameters:
923:     nlnk   - number of newly added indices
924:     lnk    - the sorted(increasing order) linked list containing new and non-redundate entries from idx
925:     lnklvl - levels of lnk
926:     bt     - updated PetscBT (bitarray) 
927:   Note: the level of U(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(prow,i)+lvl(prow,j)+1)
928:         where idx = non-zero columns of U(prow,prow+1:n-1), prow<i
929: */
930: #define PetscICCLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,idxlvl_prow) 0;\
931: {\
932:   PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
933:   nlnk = 0;\
934:   _lnkdata = idx_start;\
935:   for (_k=0; _k<nidx; _k++){\
936:     _incrlev = idxlvl[_k] + idxlvl_prow + 1;\
937:     if (_incrlev > level) continue;\
938:     _entry = idx[_k];\
939:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
940:       /* search for insertion location */\
941:       do {\
942:         _location = _lnkdata;\
943:         _lnkdata  = lnk[_location];\
944:       } while (_entry > _lnkdata);\
945:       /* insertion location is found, add entry into lnk */\
946:       lnk[_location] = _entry;\
947:       lnk[_entry]    = _lnkdata;\
948:       lnklvl[_entry] = _incrlev;\
949:       nlnk++;\
950:       _lnkdata = _entry; /* next search starts from here */\
951:     } else { /* existing entry: update lnklvl */\
952:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
953:     }\
954:   }\
955: }

957: /*
958:   Copy data on the list into an array, then initialize the list 
959:   Input Parameters:
960:     idx_start - starting index of the list 
961:     lnk_max   - max value of lnk indicating the end of the list 
962:     nlnk      - number of data on the list to be copied
963:     lnk       - linked list
964:     lnklvl    - level of lnk
965:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
966:   output Parameters:
967:     indices - array that contains the copied data
968:     lnk     - linked list that is cleaned and initialize
969:     lnklvl  - level of lnk that is reinitialized 
970:     bt      - PetscBT (bitarray) with all bits set to false
971: */
972: #define PetscIncompleteLLClean(idx_start,lnk_max,nlnk,lnk,lnklvl,indices,indiceslvl,bt) 0;\
973: {\
974:   PetscInt _j,_idx=idx_start;\
975:   for (_j=0; _j<nlnk; _j++){\
976:     _idx = lnk[_idx];\
977:     *(indices+_j) = _idx;\
978:     *(indiceslvl+_j) = lnklvl[_idx];\
979:     lnklvl[_idx] = -1;\
980:     PetscBTClear(bt,_idx);\
981:   }\
982:   lnk[idx_start] = lnk_max;\
983: }
984: /*
985:   Free memories used by the list
986: */
987: #define PetscIncompleteLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(bt))




1008: #endif