Actual source code: petscsys.h

  1: /*
  2:     Provides access to system related and general utility routines.
  3: */
 6:  #include petsc.h

  9: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetArchType(char[],size_t);
 10: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetHostName(char[],size_t);
 11: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetUserName(char[],size_t);
 12: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetProgramName(char[],size_t);
 13: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSetProgramName(const char[]);
 14: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetDate(char[],size_t);

 16: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSortInt(PetscInt,PetscInt[]);
 17: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSortIntWithPermutation(PetscInt,const PetscInt[],PetscInt[]);
 18: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSortStrWithPermutation(PetscInt,const char*[],PetscInt[]);
 19: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSortIntWithArray(PetscInt,PetscInt[],PetscInt[]);
 20: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSortIntWithScalarArray(PetscInt,PetscInt[],PetscScalar[]);
 21: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSortReal(PetscInt,PetscReal[]);
 22: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSortRealWithPermutation(PetscInt,const PetscReal[],PetscInt[]);

 24: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscSetDisplay(void);
 25: EXTERN PetscErrorCode PETSC_DLLEXPORT  PetscGetDisplay(char[],size_t);


 29: typedef enum { RANDOM_DEFAULT,RANDOM_DEFAULT_REAL,
 30:                RANDOM_DEFAULT_IMAGINARY } PetscRandomType;

 32: /*S
 33:      PetscRandom - Abstract PETSc object that manages generating random numbers

 35:    Level: intermediate

 37:   Concepts: random numbers

 39: .seealso:  PetscRandomCreate(), PetscRandomGetValue()
 40: S*/
 41: typedef struct _p_PetscRandom*   PetscRandom;

 43: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscRandomCreate(MPI_Comm,PetscRandomType,PetscRandom*);
 44: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscRandomGetValue(PetscRandom,PetscScalar*);
 45: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscRandomSetInterval(PetscRandom,PetscScalar,PetscScalar);
 46: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscRandomDestroy(PetscRandom);

 48: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetFullPath(const char[],char[],size_t);
 49: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetRelativePath(const char[],char[],size_t);
 50: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetWorkingDirectory(char[],size_t);
 51: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetRealPath(char[],char[]);
 52: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetHomeDirectory(char[],size_t);
 53: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscTestFile(const char[],char,PetscTruth*);
 54: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscTestDirectory(const char[],char,PetscTruth*);
 55: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscBinaryRead(int,void*,PetscInt,PetscDataType);
 56: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSynchronizedBinaryRead(MPI_Comm,int,void*,PetscInt,PetscDataType);
 57: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSynchronizedBinaryWrite(MPI_Comm,int,void*,PetscInt,PetscDataType,PetscTruth);
 58: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscBinaryWrite(int,void*,PetscInt,PetscDataType,PetscTruth);
 59: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscBinaryOpen(const char[],int,int *);
 60: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscBinaryClose(int);
 61: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSharedTmp(MPI_Comm,PetscTruth *);
 62: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSharedWorkingDirectory(MPI_Comm,PetscTruth *);
 63: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGetTmp(MPI_Comm,char *,size_t);
 64: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscFileRetrieve(MPI_Comm,const char *,char *,size_t,PetscTruth*);
 65: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscLs(MPI_Comm,const char[],char*,size_t,PetscTruth*);
 66: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscDLLibraryCCAAppend(MPI_Comm,PetscDLLibraryList*,const char[]);

 68: /*
 69:    In binary files variables are stored using the following lengths,
 70:   regardless of how they are stored in memory on any one particular
 71:   machine. Use these rather then sizeof() in computing sizes for 
 72:   PetscBinarySeek().
 73: */
 74: #define PETSC_BINARY_INT_SIZE    (32/8)
 75: #define PETSC_BINARY_FLOAT_SIZE  (32/8)
 76: #define PETSC_BINARY_CHAR_SIZE    (8/8)
 77: #define PETSC_BINARY_SHORT_SIZE  (16/8)
 78: #define PETSC_BINARY_DOUBLE_SIZE (64/8)
 79: #define PETSC_BINARY_SCALAR_SIZE sizeof(PetscScalar)

 81: /*E
 82:   PetscBinarySeekType - argument to PetscBinarySeek()

 84:   Level: advanced

 86: .seealso: PetscBinarySeek(), PetscSynchronizedBinarySeek()
 87: E*/
 88: typedef enum {PETSC_BINARY_SEEK_SET = 0,PETSC_BINARY_SEEK_CUR = 1,PETSC_BINARY_SEEK_END = 2} PetscBinarySeekType;
 89: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscBinarySeek(int,off_t,PetscBinarySeekType,off_t*);
 90: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSynchronizedBinarySeek(MPI_Comm,int,off_t,PetscBinarySeekType,off_t*);

 92: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSetDebugger(const char[],PetscTruth);
 93: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSetDefaultDebugger(void);
 94: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSetDebuggerFromString(char*);
 95: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscAttachDebugger(void);
 96: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscStopForDebugger(void);

 98: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGatherNumberOfMessages(MPI_Comm,PetscMPIInt*,PetscMPIInt*,PetscMPIInt*);
 99: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGatherMessageLengths(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt**,PetscMPIInt**);
100: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscGatherMessageLengths2(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt*,PetscMPIInt**,PetscMPIInt**,PetscMPIInt**);
101: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscPostIrecvInt(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt*,PetscInt***,MPI_Request**);
102: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscPostIrecvScalar(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt*,PetscScalar***,MPI_Request**);

104: EXTERN PetscErrorCode PETSC_DLLEXPORT PetscSSEIsEnabled(MPI_Comm,PetscTruth *,PetscTruth *);

106: /* Parallel communication routines */
107: /*E
108:   InsertMode - Whether entries are inserted or added into vectors or matrices

110:   Level: beginner

112: .seealso: VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(),
113:           VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(),
114:           MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd()
115: E*/
116: typedef enum {NOT_SET_VALUES, INSERT_VALUES, ADD_VALUES, MAX_VALUES} InsertMode;

118: /*M
119:     INSERT_VALUES - Put a value into a vector or matrix, overwrites any previous value

121:     Level: beginner

123: .seealso: InsertMode, VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(),
124:           VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), ADD_VALUES, INSERT_VALUES,
125:           MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd()

127: M*/

129: /*M
130:     ADD_VALUES - Adds a value into a vector or matrix, if there previously was no value, just puts the
131:                 value into that location

133:     Level: beginner

135: .seealso: InsertMode, VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(),
136:           VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), ADD_VALUES, INSERT_VALUES,
137:           MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd()

139: M*/

141: /*M
142:     MAX_VALUES - Puts the maximum of the scattered/gathered value and the current value into each location

144:     Level: beginner

146: .seealso: InsertMode, VecScatterBegin(), VecScatterEnd(), ADD_VALUES, INSERT_VALUES

148: M*/

150: /*E
151:   ScatterMode - Determines the direction of a scatter

153:   Level: beginner

155: .seealso: VecScatter, VecScatterBegin(), VecScatterEnd()
156: E*/
157: typedef enum {SCATTER_FORWARD=0, SCATTER_REVERSE=1, SCATTER_FORWARD_LOCAL=2, SCATTER_REVERSE_LOCAL=3, SCATTER_LOCAL=2} ScatterMode;

159: /*M
160:     SCATTER_FORWARD - Scatters the values as dictated by the VecScatterCreate() call

162:     Level: beginner

164: .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_REVERSE, SCATTER_FORWARD_LOCAL,
165:           SCATTER_REVERSE_LOCAL

167: M*/

169: /*M
170:     SCATTER_REVERSE - Moves the values in the opposite direction then the directions indicated in
171:          in the VecScatterCreate()

173:     Level: beginner

175: .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_FORWARD, SCATTER_FORWARD_LOCAL,
176:           SCATTER_REVERSE_LOCAL

178: M*/

180: /*M
181:     SCATTER_FORWARD_LOCAL - Scatters the values as dictated by the VecScatterCreate() call except NO parallel communication
182:        is done. Any variables that have be moved between processes are ignored

184:     Level: developer

186: .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_REVERSE, SCATTER_FORWARD,
187:           SCATTER_REVERSE_LOCAL

189: M*/

191: /*M
192:     SCATTER_REVERSE_LOCAL - Moves the values in the opposite direction then the directions indicated in
193:          in the VecScatterCreate()  except NO parallel communication
194:        is done. Any variables that have be moved between processes are ignored

196:     Level: developer

198: .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_FORWARD, SCATTER_FORWARD_LOCAL,
199:           SCATTER_REVERSE

201: M*/

203: /* 
204:   Create and initialize a linked list 
205:   Input Parameters:
206:     idx_start - starting index of the list
207:     lnk_max   - max value of lnk indicating the end of the list
208:     nlnk      - max length of the list
209:   Output Parameters:
210:     lnk       - list initialized
211:     bt        - PetscBT (bitarray) with all bits set to false
212: */
213: #define PetscLLCreate(idx_start,lnk_max,nlnk,lnk,bt) \
214:   (PetscMalloc(nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,0))

216: /*
217:   Add a index set into a sorted linked list
218:   Input Parameters:
219:     nidx      - number of input indices
220:     indices   - interger array
221:     idx_start - starting index of the list
222:     lnk       - linked list(an integer array) that is created
223:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
224:   output Parameters:
225:     nlnk      - number of newly added indices
226:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
227:     bt        - updated PetscBT (bitarray) 
228: */
229: #define PetscLLAdd(nidx,indices,idx_start,nlnk,lnk,bt) 0;\
230: {\
231:   PetscInt _k,_entry,_location,_lnkdata;\
232:   nlnk     = 0;\
233:   _lnkdata = idx_start;\
234:   for (_k=0; _k<nidx; _k++){\
235:     _entry = indices[_k];\
236:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
237:       /* search for insertion location */\
238:       /* start from the beginning if _entry < previous _entry */\
239:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
240:       do {\
241:         _location = _lnkdata;\
242:         _lnkdata  = lnk[_location];\
243:       } while (_entry > _lnkdata);\
244:       /* insertion location is found, add entry into lnk */\
245:       lnk[_location] = _entry;\
246:       lnk[_entry]    = _lnkdata;\
247:       nlnk++;\
248:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
249:     }\
250:   }\
251: }

253: /*
254:   Add a SORTED index set into a sorted linked list
255:   Input Parameters:
256:     nidx      - number of input indices
257:     indices   - sorted interger array 
258:     idx_start - starting index of the list
259:     lnk       - linked list(an integer array) that is created
260:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
261:   output Parameters:
262:     nlnk      - number of newly added indices
263:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
264:     bt        - updated PetscBT (bitarray) 
265: */
266: #define PetscLLAddSorted(nidx,indices,idx_start,nlnk,lnk,bt) 0;\
267: {\
268:   PetscInt _k,_entry,_location,_lnkdata;\
269:   nlnk      = 0;\
270:   _lnkdata  = idx_start;\
271:   for (_k=0; _k<nidx; _k++){\
272:     _entry = indices[_k];\
273:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
274:       /* search for insertion location */\
275:       do {\
276:         _location = _lnkdata;\
277:         _lnkdata  = lnk[_location];\
278:       } while (_entry > _lnkdata);\
279:       /* insertion location is found, add entry into lnk */\
280:       lnk[_location] = _entry;\
281:       lnk[_entry]    = _lnkdata;\
282:       nlnk++;\
283:       _lnkdata = _entry; /* next search starts from here */\
284:     }\
285:   }\
286: }

288: /*
289:   Add a SORTED index set into a sorted linked list used for LUFactorSymbolic()
290:   Same as PetscLLAddSorted() with an additional operation:
291:        count the number of input indices that are no larger than 'diag'
292:   Input Parameters:
293:     nidx      - number of input indices
294:     indices   - sorted interger array 
295:     idx_start - starting index of the list
296:     lnk       - linked list(an integer array) that is created
297:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
298:     diag      - index of the active row in LUFactorSymbolic
299:   output Parameters:
300:     nlnk      - number of newly added indices
301:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from indices
302:     bt        - updated PetscBT (bitarray) 
303:     nzbd      - number of input indices that are no larger than 'diag'
304: */
305: #define PetscLLAddSortedLU(nidx,indices,idx_start,nlnk,lnk,bt,diag,nzbd) 0;\
306: {\
307:   PetscInt _k,_entry,_location,_lnkdata;\
308:   nlnk     = 0;\
309:   _lnkdata = idx_start;\
310:   nzbd     = 0;\
311:   for (_k=0; _k<nidx; _k++){\
312:     _entry = indices[_k];\
313:     if (_entry <= diag) nzbd++;\
314:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
315:       /* search for insertion location */\
316:       do {\
317:         _location = _lnkdata;\
318:         _lnkdata  = lnk[_location];\
319:       } while (_entry > _lnkdata);\
320:       /* insertion location is found, add entry into lnk */\
321:       lnk[_location] = _entry;\
322:       lnk[_entry]    = _lnkdata;\
323:       nlnk++;\
324:       _lnkdata = _entry; /* next search starts from here */\
325:     }\
326:   }\
327: }

329: /*
330:   Copy data on the list into an array, then initialize the list 
331:   Input Parameters:
332:     idx_start - starting index of the list 
333:     lnk_max   - max value of lnk indicating the end of the list 
334:     nlnk      - number of data on the list to be copied
335:     lnk       - linked list
336:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
337:   output Parameters:
338:     indices   - array that contains the copied data
339:     lnk       - linked list that is cleaned and initialize
340:     bt        - PetscBT (bitarray) with all bits set to false
341: */
342: #define PetscLLClean(idx_start,lnk_max,nlnk,lnk,indices,bt) 0;\
343: {\
344:   PetscInt _j,_idx=idx_start;\
345:   for (_j=0; _j<nlnk; _j++){\
346:     _idx = lnk[_idx];\
347:     *(indices+_j) = _idx;\
348:     PetscBTClear(bt,_idx);\
349:   }\
350:   lnk[idx_start] = lnk_max;\
351: }
352: /*
353:   Free memories used by the list
354: */
355: #define PetscLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(bt))

357: /* Routines below are used for incomplete matrix factorization */
358: /* 
359:   Create and initialize a linked list and its levels
360:   Input Parameters:
361:     idx_start - starting index of the list
362:     lnk_max   - max value of lnk indicating the end of the list
363:     nlnk      - max length of the list
364:   Output Parameters:
365:     lnk       - list initialized
366:     lnk_lvl   - array of size nlnk for storing levels of lnk
367:     bt        - PetscBT (bitarray) with all bits set to false
368: */
369: #define PetscIncompleteLLCreate(idx_start,lnk_max,nlnk,lnk,lnk_lvl,bt)\
370:   (PetscMalloc(2*nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,lnk_lvl = lnk + nlnk,0))

372: /*
373:   Initialize a sorted linked list used for ILU and ICC
374:   Input Parameters:
375:     nidx      - number of input idx
376:     idx       - interger array used for storing column indices
377:     idx_start - starting index of the list
378:     perm      - indices of an IS
379:     lnk       - linked list(an integer array) that is created
380:     lnklvl    - levels of lnk
381:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
382:   output Parameters:
383:     nlnk     - number of newly added idx
384:     lnk      - the sorted(increasing order) linked list containing new and non-redundate entries from idx
385:     lnklvl   - levels of lnk
386:     bt       - updated PetscBT (bitarray) 
387: */
388: #define PetscIncompleteLLInit(nidx,idx,idx_start,perm,nlnk,lnk,lnklvl,bt) 0;\
389: {\
390:   PetscInt _k,_entry,_location,_lnkdata;\
391:   nlnk     = 0;\
392:   _lnkdata = idx_start;\
393:   for (_k=0; _k<nidx; _k++){\
394:     _entry = perm[idx[_k]];\
395:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
396:       /* search for insertion location */\
397:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
398:       do {\
399:         _location = _lnkdata;\
400:         _lnkdata  = lnk[_location];\
401:       } while (_entry > _lnkdata);\
402:       /* insertion location is found, add entry into lnk */\
403:       lnk[_location]  = _entry;\
404:       lnk[_entry]     = _lnkdata;\
405:       lnklvl[_entry] = 0;\
406:       nlnk++;\
407:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
408:     }\
409:   }\
410: }

412: /*
413:   Add a SORTED index set into a sorted linked list for ILU
414:   Input Parameters:
415:     nidx      - number of input indices
416:     idx       - sorted interger array used for storing column indices
417:     level     - level of fill, e.g., ICC(level)
418:     idxlvl    - level of idx 
419:     idx_start - starting index of the list
420:     lnk       - linked list(an integer array) that is created
421:     lnklvl    - levels of lnk
422:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
423:     prow      - the row number of idx
424:   output Parameters:
425:     nlnk     - number of newly added idx
426:     lnk      - the sorted(increasing order) linked list containing new and non-redundate entries from idx
427:     lnklvl   - levels of lnk
428:     bt       - updated PetscBT (bitarray) 

430:   Note: the level of factor(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(i,prow)+lvl(prow,j)+1)
431:         where idx = non-zero columns of U(prow,prow+1:n-1), prow<i
432: */
433: #define PetscILULLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,lnklvl_prow) 0;\
434: {\
435:   PetscInt _k,_entry,_location,_lnkdata,_incrlev,_lnklvl_prow=lnklvl[prow];\
436:   nlnk     = 0;\
437:   _lnkdata = idx_start;\
438:   for (_k=0; _k<nidx; _k++){\
439:     _incrlev = idxlvl[_k] + _lnklvl_prow + 1;\
440:     if (_incrlev > level) continue;\
441:     _entry = idx[_k];\
442:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
443:       /* search for insertion location */\
444:       do {\
445:         _location = _lnkdata;\
446:         _lnkdata  = lnk[_location];\
447:       } while (_entry > _lnkdata);\
448:       /* insertion location is found, add entry into lnk */\
449:       lnk[_location]  = _entry;\
450:       lnk[_entry]     = _lnkdata;\
451:       lnklvl[_entry] = _incrlev;\
452:       nlnk++;\
453:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
454:     } else { /* existing entry: update lnklvl */\
455:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
456:     }\
457:   }\
458: }

460: /*
461:   Add a index set into a sorted linked list
462:   Input Parameters:
463:     nidx      - number of input idx
464:     idx   - interger array used for storing column indices
465:     level     - level of fill, e.g., ICC(level)
466:     idxlvl - level of idx 
467:     idx_start - starting index of the list
468:     lnk       - linked list(an integer array) that is created
469:     lnklvl   - levels of lnk
470:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
471:   output Parameters:
472:     nlnk      - number of newly added idx
473:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from idx
474:     lnklvl   - levels of lnk
475:     bt        - updated PetscBT (bitarray) 
476: */
477: #define PetscIncompleteLLAdd(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\
478: {\
479:   PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
480:   nlnk     = 0;\
481:   _lnkdata = idx_start;\
482:   for (_k=0; _k<nidx; _k++){\
483:     _incrlev = idxlvl[_k] + 1;\
484:     if (_incrlev > level) continue;\
485:     _entry = idx[_k];\
486:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
487:       /* search for insertion location */\
488:       if (_k && _entry < _lnkdata) _lnkdata  = idx_start;\
489:       do {\
490:         _location = _lnkdata;\
491:         _lnkdata  = lnk[_location];\
492:       } while (_entry > _lnkdata);\
493:       /* insertion location is found, add entry into lnk */\
494:       lnk[_location]  = _entry;\
495:       lnk[_entry]     = _lnkdata;\
496:       lnklvl[_entry] = _incrlev;\
497:       nlnk++;\
498:       _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\
499:     } else { /* existing entry: update lnklvl */\
500:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
501:     }\
502:   }\
503: }

505: /*
506:   Add a SORTED index set into a sorted linked list
507:   Input Parameters:
508:     nidx      - number of input indices
509:     idx   - sorted interger array used for storing column indices
510:     level     - level of fill, e.g., ICC(level)
511:     idxlvl - level of idx 
512:     idx_start - starting index of the list
513:     lnk       - linked list(an integer array) that is created
514:     lnklvl    - levels of lnk
515:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
516:   output Parameters:
517:     nlnk      - number of newly added idx
518:     lnk       - the sorted(increasing order) linked list containing new and non-redundate entries from idx
519:     lnklvl    - levels of lnk
520:     bt        - updated PetscBT (bitarray) 
521: */
522: #define PetscIncompleteLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\
523: {\
524:   PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
525:   nlnk = 0;\
526:   _lnkdata = idx_start;\
527:   for (_k=0; _k<nidx; _k++){\
528:     _incrlev = idxlvl[_k] + 1;\
529:     if (_incrlev > level) continue;\
530:     _entry = idx[_k];\
531:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
532:       /* search for insertion location */\
533:       do {\
534:         _location = _lnkdata;\
535:         _lnkdata  = lnk[_location];\
536:       } while (_entry > _lnkdata);\
537:       /* insertion location is found, add entry into lnk */\
538:       lnk[_location] = _entry;\
539:       lnk[_entry]    = _lnkdata;\
540:       lnklvl[_entry] = _incrlev;\
541:       nlnk++;\
542:       _lnkdata = _entry; /* next search starts from here */\
543:     } else { /* existing entry: update lnklvl */\
544:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
545:     }\
546:   }\
547: }

549: /*
550:   Add a SORTED index set into a sorted linked list for ICC
551:   Input Parameters:
552:     nidx      - number of input indices
553:     idx       - sorted interger array used for storing column indices
554:     level     - level of fill, e.g., ICC(level)
555:     idxlvl    - level of idx 
556:     idx_start - starting index of the list
557:     lnk       - linked list(an integer array) that is created
558:     lnklvl    - levels of lnk
559:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
560:     idxlvl_prow - idxlvl[prow], where prow is the row number of the idx
561:   output Parameters:
562:     nlnk   - number of newly added indices
563:     lnk    - the sorted(increasing order) linked list containing new and non-redundate entries from idx
564:     lnklvl - levels of lnk
565:     bt     - updated PetscBT (bitarray) 
566:   Note: the level of U(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(prow,i)+lvl(prow,j)+1)
567:         where idx = non-zero columns of U(prow,prow+1:n-1), prow<i
568: */
569: #define PetscICCLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,idxlvl_prow) 0;\
570: {\
571:   PetscInt _k,_entry,_location,_lnkdata,_incrlev;\
572:   nlnk = 0;\
573:   _lnkdata = idx_start;\
574:   for (_k=0; _k<nidx; _k++){\
575:     _incrlev = idxlvl[_k] + idxlvl_prow + 1;\
576:     if (_incrlev > level) continue;\
577:     _entry = idx[_k];\
578:     if (!PetscBTLookupSet(bt,_entry)){  /* new entry */\
579:       /* search for insertion location */\
580:       do {\
581:         _location = _lnkdata;\
582:         _lnkdata  = lnk[_location];\
583:       } while (_entry > _lnkdata);\
584:       /* insertion location is found, add entry into lnk */\
585:       lnk[_location] = _entry;\
586:       lnk[_entry]    = _lnkdata;\
587:       lnklvl[_entry] = _incrlev;\
588:       nlnk++;\
589:       _lnkdata = _entry; /* next search starts from here */\
590:     } else { /* existing entry: update lnklvl */\
591:       if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\
592:     }\
593:   }\
594: }

596: /*
597:   Copy data on the list into an array, then initialize the list 
598:   Input Parameters:
599:     idx_start - starting index of the list 
600:     lnk_max   - max value of lnk indicating the end of the list 
601:     nlnk      - number of data on the list to be copied
602:     lnk       - linked list
603:     lnklvl    - level of lnk
604:     bt        - PetscBT (bitarray), bt[idx]=true marks idx is in lnk
605:   output Parameters:
606:     indices - array that contains the copied data
607:     lnk     - linked list that is cleaned and initialize
608:     lnklvl  - level of lnk that is reinitialized 
609:     bt      - PetscBT (bitarray) with all bits set to false
610: */
611: #define PetscIncompleteLLClean(idx_start,lnk_max,nlnk,lnk,lnklvl,indices,indiceslvl,bt) 0;\
612: {\
613:   PetscInt _j,_idx=idx_start;\
614:   for (_j=0; _j<nlnk; _j++){\
615:     _idx = lnk[_idx];\
616:     *(indices+_j) = _idx;\
617:     *(indiceslvl+_j) = lnklvl[_idx];\
618:     lnklvl[_idx] = -1;\
619:     PetscBTClear(bt,_idx);\
620:   }\
621:   lnk[idx_start] = lnk_max;\
622: }
623: /*
624:   Free memories used by the list
625: */
626: #define PetscIncompleteLLDestroy(lnk,bt) \
627:   (PetscFree(lnk) || PetscBTDestroy(bt) || (0))

630: #endif /* __PETSCSYS_H */