1 /* 2 Provides access to system related and general utility routines. 3 */ 4 #if !defined(__PETSCSYS_H) 5 #define __PETSCSYS_H 6 #include "petsc.h" 7 PETSC_EXTERN_CXX_BEGIN 8 9 EXTERN PetscErrorCode PetscGetArchType(char[],size_t); 10 EXTERN PetscErrorCode PetscGetHostName(char[],size_t); 11 EXTERN PetscErrorCode PetscGetUserName(char[],size_t); 12 EXTERN PetscErrorCode PetscGetProgramName(char[],size_t); 13 EXTERN PetscErrorCode PetscSetProgramName(const char[]); 14 EXTERN PetscErrorCode PetscGetDate(char[],size_t); 15 16 EXTERN PetscErrorCode PetscSortInt(PetscInt,PetscInt[]); 17 EXTERN PetscErrorCode PetscSortIntWithPermutation(PetscInt,const PetscInt[],PetscInt[]); 18 EXTERN PetscErrorCode PetscSortStrWithPermutation(PetscInt,const char*[],PetscInt[]); 19 EXTERN PetscErrorCode PetscSortIntWithArray(PetscInt,PetscInt[],PetscInt[]); 20 EXTERN PetscErrorCode PetscSortIntWithScalarArray(PetscInt,PetscInt[],PetscScalar[]); 21 EXTERN PetscErrorCode PetscSortReal(PetscInt,PetscReal[]); 22 EXTERN PetscErrorCode PetscSortRealWithPermutation(PetscInt,const PetscReal[],PetscInt[]); 23 24 EXTERN PetscErrorCode PetscSetDisplay(void); 25 EXTERN PetscErrorCode PetscGetDisplay(char[],size_t); 26 27 extern PetscCookie PETSC_RANDOM_COOKIE; 28 29 typedef enum { RANDOM_DEFAULT,RANDOM_DEFAULT_REAL, 30 RANDOM_DEFAULT_IMAGINARY } PetscRandomType; 31 32 /*S 33 PetscRandom - Abstract PETSc object that manages generating random numbers 34 35 Level: intermediate 36 37 Concepts: random numbers 38 39 .seealso: PetscRandomCreate(), PetscRandomGetValue() 40 S*/ 41 typedef struct _p_PetscRandom* PetscRandom; 42 43 EXTERN PetscErrorCode PetscRandomCreate(MPI_Comm,PetscRandomType,PetscRandom*); 44 EXTERN PetscErrorCode PetscRandomGetValue(PetscRandom,PetscScalar*); 45 EXTERN PetscErrorCode PetscRandomSetInterval(PetscRandom,PetscScalar,PetscScalar); 46 EXTERN PetscErrorCode PetscRandomDestroy(PetscRandom); 47 48 EXTERN PetscErrorCode PetscGetFullPath(const char[],char[],size_t); 49 EXTERN PetscErrorCode PetscGetRelativePath(const char[],char[],size_t); 50 EXTERN PetscErrorCode PetscGetWorkingDirectory(char[],size_t); 51 EXTERN PetscErrorCode PetscGetRealPath(char[],char[]); 52 EXTERN PetscErrorCode PetscGetHomeDirectory(char[],size_t); 53 EXTERN PetscErrorCode PetscTestFile(const char[],char,PetscTruth*); 54 EXTERN PetscErrorCode PetscTestDirectory(const char[],char,PetscTruth*); 55 EXTERN PetscErrorCode PetscBinaryRead(int,void*,PetscInt,PetscDataType); 56 EXTERN PetscErrorCode PetscSynchronizedBinaryRead(MPI_Comm,int,void*,PetscInt,PetscDataType); 57 EXTERN PetscErrorCode PetscBinaryWrite(int,void*,PetscInt,PetscDataType,PetscTruth); 58 EXTERN PetscErrorCode PetscBinaryOpen(const char[],int,int *); 59 EXTERN PetscErrorCode PetscBinaryClose(int); 60 EXTERN PetscErrorCode PetscSharedTmp(MPI_Comm,PetscTruth *); 61 EXTERN PetscErrorCode PetscSharedWorkingDirectory(MPI_Comm,PetscTruth *); 62 EXTERN PetscErrorCode PetscGetTmp(MPI_Comm,char *,size_t); 63 EXTERN PetscErrorCode PetscFileRetrieve(MPI_Comm,const char *,char *,size_t,PetscTruth*); 64 EXTERN PetscErrorCode PetscLs(MPI_Comm,const char[],char*,size_t,PetscTruth*); 65 EXTERN PetscErrorCode PetscDLLibraryCCAAppend(MPI_Comm,PetscDLLibraryList*,const char[]); 66 67 /* 68 In binary files variables are stored using the following lengths, 69 regardless of how they are stored in memory on any one particular 70 machine. Use these rather then sizeof() in computing sizes for 71 PetscBinarySeek(). 72 */ 73 #define PETSC_BINARY_INT_SIZE (32/8) 74 #define PETSC_BINARY_FLOAT_SIZE (32/8) 75 #define PETSC_BINARY_CHAR_SIZE (8/8) 76 #define PETSC_BINARY_SHORT_SIZE (16/8) 77 #define PETSC_BINARY_DOUBLE_SIZE (64/8) 78 #define PETSC_BINARY_SCALAR_SIZE sizeof(PetscScalar) 79 80 /*E 81 PetscBinarySeekType - argument to PetscBinarySeek() 82 83 Level: advanced 84 85 .seealso: PetscBinarySeek(), PetscSynchronizedBinarySeek() 86 E*/ 87 typedef enum {PETSC_BINARY_SEEK_SET = 0,PETSC_BINARY_SEEK_CUR = 1,PETSC_BINARY_SEEK_END = 2} PetscBinarySeekType; 88 EXTERN PetscErrorCode PetscBinarySeek(int,off_t,PetscBinarySeekType,off_t*); 89 EXTERN PetscErrorCode PetscSynchronizedBinarySeek(MPI_Comm,int,off_t,PetscBinarySeekType,off_t*); 90 91 EXTERN PetscErrorCode PetscSetDebugger(const char[],PetscTruth); 92 EXTERN PetscErrorCode PetscSetDefaultDebugger(void); 93 EXTERN PetscErrorCode PetscSetDebuggerFromString(char*); 94 EXTERN PetscErrorCode PetscAttachDebugger(void); 95 EXTERN PetscErrorCode PetscStopForDebugger(void); 96 97 EXTERN PetscErrorCode PetscGatherNumberOfMessages(MPI_Comm,PetscMPIInt*,PetscMPIInt*,PetscMPIInt*); 98 EXTERN PetscErrorCode PetscGatherMessageLengths(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt**,PetscMPIInt**); 99 EXTERN PetscErrorCode PetscGatherMessageLengths2(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt*,PetscMPIInt**,PetscMPIInt**,PetscMPIInt**); 100 EXTERN PetscErrorCode PetscPostIrecvInt(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt*,PetscInt***,MPI_Request**); 101 EXTERN PetscErrorCode PetscPostIrecvScalar(MPI_Comm,PetscMPIInt,PetscMPIInt,PetscMPIInt*,PetscMPIInt*,PetscScalar***,MPI_Request**); 102 103 EXTERN PetscErrorCode PetscSSEIsEnabled(MPI_Comm,PetscTruth *,PetscTruth *); 104 105 /* Parallel communication routines */ 106 /*E 107 InsertMode - Whether entries are inserted or added into vectors or matrices 108 109 Level: beginner 110 111 .seealso: VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(), 112 VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), 113 MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd() 114 E*/ 115 typedef enum {NOT_SET_VALUES, INSERT_VALUES, ADD_VALUES, MAX_VALUES} InsertMode; 116 117 /*M 118 INSERT_VALUES - Put a value into a vector or matrix, overwrites any previous value 119 120 Level: beginner 121 122 .seealso: InsertMode, VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(), 123 VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), ADD_VALUES, INSERT_VALUES, 124 MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd() 125 126 M*/ 127 128 /*M 129 ADD_VALUES - Adds a value into a vector or matrix, if there previously was no value, just puts the 130 value into that location 131 132 Level: beginner 133 134 .seealso: InsertMode, VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(), 135 VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), ADD_VALUES, INSERT_VALUES, 136 MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd() 137 138 M*/ 139 140 /*M 141 MAX_VALUES - Puts the maximum of the scattered/gathered value and the current value into each location 142 143 Level: beginner 144 145 .seealso: InsertMode, VecScatterBegin(), VecScatterEnd(), ADD_VALUES, INSERT_VALUES 146 147 M*/ 148 149 /*E 150 ScatterMode - Determines the direction of a scatter 151 152 Level: beginner 153 154 .seealso: VecScatter, VecScatterBegin(), VecScatterEnd() 155 E*/ 156 typedef enum {SCATTER_FORWARD=0, SCATTER_REVERSE=1, SCATTER_FORWARD_LOCAL=2, SCATTER_REVERSE_LOCAL=3, SCATTER_LOCAL=2} ScatterMode; 157 158 /*M 159 SCATTER_FORWARD - Scatters the values as dictated by the VecScatterCreate() call 160 161 Level: beginner 162 163 .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_REVERSE, SCATTER_FORWARD_LOCAL, 164 SCATTER_REVERSE_LOCAL 165 166 M*/ 167 168 /*M 169 SCATTER_REVERSE - Moves the values in the opposite direction then the directions indicated in 170 in the VecScatterCreate() 171 172 Level: beginner 173 174 .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_FORWARD, SCATTER_FORWARD_LOCAL, 175 SCATTER_REVERSE_LOCAL 176 177 M*/ 178 179 /*M 180 SCATTER_FORWARD_LOCAL - Scatters the values as dictated by the VecScatterCreate() call except NO parallel communication 181 is done. Any variables that have be moved between processes are ignored 182 183 Level: developer 184 185 .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_REVERSE, SCATTER_FORWARD, 186 SCATTER_REVERSE_LOCAL 187 188 M*/ 189 190 /*M 191 SCATTER_REVERSE_LOCAL - Moves the values in the opposite direction then the directions indicated in 192 in the VecScatterCreate() except NO parallel communication 193 is done. Any variables that have be moved between processes are ignored 194 195 Level: developer 196 197 .seealso: VecScatter, ScatterMode, VecScatterCreate(), VecScatterBegin(), VecScatterEnd(), SCATTER_FORWARD, SCATTER_FORWARD_LOCAL, 198 SCATTER_REVERSE 199 200 M*/ 201 202 /* 203 Create and initialize a linked list 204 Input Parameters: 205 idx_start - starting index of the list 206 lnk_max - max value of lnk indicating the end of the list 207 nlnk - max length of the list 208 Output Parameters: 209 lnk - list initialized 210 bt - PetscBT (bitarray) with all bits set to false 211 */ 212 #define PetscLLCreate(idx_start,lnk_max,nlnk,lnk,bt) \ 213 (PetscMalloc(nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,0)) 214 215 /* 216 Add a index set into a sorted linked list 217 Input Parameters: 218 nidx - number of input indices 219 indices - interger array 220 idx_start - starting index of the list 221 lnk - linked list(an integer array) that is created 222 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 223 output Parameters: 224 nlnk - number of newly added indices 225 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 226 bt - updated PetscBT (bitarray) 227 */ 228 #define PetscLLAdd(nidx,indices,idx_start,nlnk,lnk,bt) 0;\ 229 {\ 230 int _k,_entry,_location,_lnkdata;\ 231 nlnk = 0;\ 232 _lnkdata = idx_start;\ 233 for (_k=0; _k<nidx; _k++){\ 234 _entry = indices[_k];\ 235 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 236 /* search for insertion location */\ 237 /* start from the beginning if _entry < previous _entry */\ 238 if (_k && _entry < _lnkdata) _lnkdata = idx_start;\ 239 do {\ 240 _location = _lnkdata;\ 241 _lnkdata = lnk[_location];\ 242 } while (_entry > _lnkdata);\ 243 /* insertion location is found, add entry into lnk */\ 244 lnk[_location] = _entry;\ 245 lnk[_entry] = _lnkdata;\ 246 nlnk++;\ 247 _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\ 248 }\ 249 }\ 250 } 251 252 /* 253 Add a SORTED index set into a sorted linked list 254 Input Parameters: 255 nidx - number of input indices 256 indices - sorted interger array 257 idx_start - starting index of the list 258 lnk - linked list(an integer array) that is created 259 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 260 output Parameters: 261 nlnk - number of newly added indices 262 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 263 bt - updated PetscBT (bitarray) 264 */ 265 #define PetscLLAddSorted(nidx,indices,idx_start,nlnk,lnk,bt) 0;\ 266 {\ 267 int _k,_entry,_location,_lnkdata;\ 268 nlnk = 0;\ 269 _lnkdata = idx_start;\ 270 for (_k=0; _k<nidx; _k++){\ 271 _entry = indices[_k];\ 272 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 273 /* search for insertion location */\ 274 do {\ 275 _location = _lnkdata;\ 276 _lnkdata = lnk[_location];\ 277 } while (_entry > _lnkdata);\ 278 /* insertion location is found, add entry into lnk */\ 279 lnk[_location] = _entry;\ 280 lnk[_entry] = _lnkdata;\ 281 nlnk++;\ 282 _lnkdata = _entry; /* next search starts from here */\ 283 }\ 284 }\ 285 } 286 287 /* 288 Add a SORTED index set into a sorted linked list used for LUFactorSymbolic() 289 Same as PetscLLAddSorted() with an additional operation: 290 count the number of input indices that are no larger than 'diag' 291 Input Parameters: 292 nidx - number of input indices 293 indices - sorted interger array 294 idx_start - starting index of the list 295 lnk - linked list(an integer array) that is created 296 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 297 diag - index of the active row in LUFactorSymbolic 298 output Parameters: 299 nlnk - number of newly added indices 300 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 301 bt - updated PetscBT (bitarray) 302 nzbd - number of input indices that are no larger than 'diag' 303 */ 304 #define PetscLLAddSortedLU(nidx,indices,idx_start,nlnk,lnk,bt,diag,nzbd) 0;\ 305 {\ 306 int _k,_entry,_location,_lnkdata;\ 307 nlnk = 0;\ 308 _lnkdata = idx_start;\ 309 nzbd = 0;\ 310 for (_k=0; _k<nidx; _k++){\ 311 _entry = indices[_k];\ 312 if (_entry <= diag) nzbd++;\ 313 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 314 /* search for insertion location */\ 315 do {\ 316 _location = _lnkdata;\ 317 _lnkdata = lnk[_location];\ 318 } while (_entry > _lnkdata);\ 319 /* insertion location is found, add entry into lnk */\ 320 lnk[_location] = _entry;\ 321 lnk[_entry] = _lnkdata;\ 322 nlnk++;\ 323 _lnkdata = _entry; /* next search starts from here */\ 324 }\ 325 }\ 326 } 327 328 /* 329 Copy data on the list into an array, then initialize the list 330 Input Parameters: 331 idx_start - starting index of the list 332 lnk_max - max value of lnk indicating the end of the list 333 nlnk - number of data on the list to be copied 334 lnk - linked list 335 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 336 output Parameters: 337 indices - array that contains the copied data 338 lnk - linked list that is cleaned and initialize 339 bt - PetscBT (bitarray) with all bits set to false 340 */ 341 #define PetscLLClean(idx_start,lnk_max,nlnk,lnk,indices,bt) 0;\ 342 {\ 343 int _j,_idx=idx_start;\ 344 for (_j=0; _j<nlnk; _j++){\ 345 _idx = lnk[_idx];\ 346 *(indices+_j) = _idx;\ 347 PetscBTClear(bt,_idx);\ 348 }\ 349 lnk[idx_start] = lnk_max;\ 350 } 351 /* 352 Free memories used by the list 353 */ 354 #define PetscLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(bt)) 355 356 /* Routines below are used for incomplete matrix factorization */ 357 /* 358 Create and initialize a linked list and its levels 359 Input Parameters: 360 idx_start - starting index of the list 361 lnk_max - max value of lnk indicating the end of the list 362 nlnk - max length of the list 363 Output Parameters: 364 lnk - list initialized 365 lnk_lvl - array of size nlnk for storing levels of lnk 366 bt - PetscBT (bitarray) with all bits set to false 367 */ 368 #define PetscIncompleteLLCreate(idx_start,lnk_max,nlnk,lnk,lnk_lvl,bt)\ 369 (PetscMalloc(2*nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,lnk_lvl = lnk + nlnk,0)) 370 371 /* 372 Initialize a sorted linked list used for ILU and ICC 373 Input Parameters: 374 nidx - number of input idx 375 idx - interger array used for storing column indices 376 idx_start - starting index of the list 377 perm - indices of an IS 378 lnk - linked list(an integer array) that is created 379 lnklvl - levels of lnk 380 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 381 output Parameters: 382 nlnk - number of newly added idx 383 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx 384 lnklvl - levels of lnk 385 bt - updated PetscBT (bitarray) 386 */ 387 #define PetscIncompleteLLInit(nidx,idx,idx_start,perm,nlnk,lnk,lnklvl,bt) 0;\ 388 {\ 389 int _k,_entry,_location,_lnkdata;\ 390 nlnk = 0;\ 391 _lnkdata = idx_start;\ 392 for (_k=0; _k<nidx; _k++){\ 393 _entry = perm[idx[_k]];\ 394 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 395 /* search for insertion location */\ 396 if (_k && _entry < _lnkdata) _lnkdata = idx_start;\ 397 do {\ 398 _location = _lnkdata;\ 399 _lnkdata = lnk[_location];\ 400 } while (_entry > _lnkdata);\ 401 /* insertion location is found, add entry into lnk */\ 402 lnk[_location] = _entry;\ 403 lnk[_entry] = _lnkdata;\ 404 lnklvl[_entry] = 0;\ 405 nlnk++;\ 406 _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\ 407 }\ 408 }\ 409 } 410 411 /* 412 Add a SORTED index set into a sorted linked list for ILU 413 Input Parameters: 414 nidx - number of input indices 415 idx - sorted interger array used for storing column indices 416 level - level of fill, e.g., ICC(level) 417 idxlvl - level of idx 418 idx_start - starting index of the list 419 lnk - linked list(an integer array) that is created 420 lnklvl - levels of lnk 421 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 422 prow - the row number of idx 423 output Parameters: 424 nlnk - number of newly added idx 425 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx 426 lnklvl - levels of lnk 427 bt - updated PetscBT (bitarray) 428 429 Note: the level of factor(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(i,prow)+lvl(prow,j)+1) 430 where idx = non-zero columns of U(prow,prow+1:n-1), prow<i 431 */ 432 #define PetscILULLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,lnklvl_prow) 0;\ 433 {\ 434 int _k,_entry,_location,_lnkdata,_incrlev,_lnklvl_prow=lnklvl[prow];\ 435 nlnk = 0;\ 436 _lnkdata = idx_start;\ 437 for (_k=0; _k<nidx; _k++){\ 438 _incrlev = idxlvl[_k] + _lnklvl_prow + 1;\ 439 if (_incrlev > level) continue;\ 440 _entry = idx[_k];\ 441 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 442 /* search for insertion location */\ 443 do {\ 444 _location = _lnkdata;\ 445 _lnkdata = lnk[_location];\ 446 } while (_entry > _lnkdata);\ 447 /* insertion location is found, add entry into lnk */\ 448 lnk[_location] = _entry;\ 449 lnk[_entry] = _lnkdata;\ 450 lnklvl[_entry] = _incrlev;\ 451 nlnk++;\ 452 _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\ 453 } else { /* existing entry: update lnklvl */\ 454 if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\ 455 }\ 456 }\ 457 } 458 459 /* 460 Add a index set into a sorted linked list 461 Input Parameters: 462 nidx - number of input idx 463 idx - interger array used for storing column indices 464 level - level of fill, e.g., ICC(level) 465 idxlvl - level of idx 466 idx_start - starting index of the list 467 lnk - linked list(an integer array) that is created 468 lnklvl - levels of lnk 469 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 470 output Parameters: 471 nlnk - number of newly added idx 472 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx 473 lnklvl - levels of lnk 474 bt - updated PetscBT (bitarray) 475 */ 476 #define PetscIncompleteLLAdd(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\ 477 {\ 478 int _k,_entry,_location,_lnkdata,_incrlev;\ 479 nlnk = 0;\ 480 _lnkdata = idx_start;\ 481 for (_k=0; _k<nidx; _k++){\ 482 _incrlev = idxlvl[_k] + 1;\ 483 if (_incrlev > level) continue;\ 484 _entry = idx[_k];\ 485 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 486 /* search for insertion location */\ 487 if (_k && _entry < _lnkdata) _lnkdata = idx_start;\ 488 do {\ 489 _location = _lnkdata;\ 490 _lnkdata = lnk[_location];\ 491 } while (_entry > _lnkdata);\ 492 /* insertion location is found, add entry into lnk */\ 493 lnk[_location] = _entry;\ 494 lnk[_entry] = _lnkdata;\ 495 lnklvl[_entry] = _incrlev;\ 496 nlnk++;\ 497 _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\ 498 } else { /* existing entry: update lnklvl */\ 499 if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\ 500 }\ 501 }\ 502 } 503 504 /* 505 Add a SORTED index set into a sorted linked list 506 Input Parameters: 507 nidx - number of input indices 508 idx - sorted interger array used for storing column indices 509 level - level of fill, e.g., ICC(level) 510 idxlvl - level of idx 511 idx_start - starting index of the list 512 lnk - linked list(an integer array) that is created 513 lnklvl - levels of lnk 514 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 515 output Parameters: 516 nlnk - number of newly added idx 517 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx 518 lnklvl - levels of lnk 519 bt - updated PetscBT (bitarray) 520 */ 521 #define PetscIncompleteLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt) 0;\ 522 {\ 523 int _k,_entry,_location,_lnkdata,_incrlev;\ 524 nlnk = 0;\ 525 _lnkdata = idx_start;\ 526 for (_k=0; _k<nidx; _k++){\ 527 _incrlev = idxlvl[_k] + 1;\ 528 if (_incrlev > level) continue;\ 529 _entry = idx[_k];\ 530 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 531 /* search for insertion location */\ 532 do {\ 533 _location = _lnkdata;\ 534 _lnkdata = lnk[_location];\ 535 } while (_entry > _lnkdata);\ 536 /* insertion location is found, add entry into lnk */\ 537 lnk[_location] = _entry;\ 538 lnk[_entry] = _lnkdata;\ 539 lnklvl[_entry] = _incrlev;\ 540 nlnk++;\ 541 _lnkdata = _entry; /* next search starts from here */\ 542 } else { /* existing entry: update lnklvl */\ 543 if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\ 544 }\ 545 }\ 546 } 547 548 /* 549 Add a SORTED index set into a sorted linked list for ICC 550 Input Parameters: 551 nidx - number of input indices 552 idx - sorted interger array used for storing column indices 553 level - level of fill, e.g., ICC(level) 554 idxlvl - level of idx 555 idx_start - starting index of the list 556 lnk - linked list(an integer array) that is created 557 lnklvl - levels of lnk 558 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 559 idxlvl_prow - idxlvl[prow], where prow is the row number of the idx 560 output Parameters: 561 nlnk - number of newly added indices 562 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from idx 563 lnklvl - levels of lnk 564 bt - updated PetscBT (bitarray) 565 Note: the level of U(i,j) is set as lvl(i,j) = min{ lvl(i,j), lvl(prow,i)+lvl(prow,j)+1) 566 where idx = non-zero columns of U(prow,prow+1:n-1), prow<i 567 */ 568 #define PetscICCLLAddSorted(nidx,idx,level,idxlvl,idx_start,nlnk,lnk,lnklvl,bt,idxlvl_prow) 0;\ 569 {\ 570 int _k,_entry,_location,_lnkdata,_incrlev;\ 571 nlnk = 0;\ 572 _lnkdata = idx_start;\ 573 for (_k=0; _k<nidx; _k++){\ 574 _incrlev = idxlvl[_k] + idxlvl_prow + 1;\ 575 if (_incrlev > level) continue;\ 576 _entry = idx[_k];\ 577 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 578 /* search for insertion location */\ 579 do {\ 580 _location = _lnkdata;\ 581 _lnkdata = lnk[_location];\ 582 } while (_entry > _lnkdata);\ 583 /* insertion location is found, add entry into lnk */\ 584 lnk[_location] = _entry;\ 585 lnk[_entry] = _lnkdata;\ 586 lnklvl[_entry] = _incrlev;\ 587 nlnk++;\ 588 _lnkdata = _entry; /* next search starts from here */\ 589 } else { /* existing entry: update lnklvl */\ 590 if (lnklvl[_entry] > _incrlev) lnklvl[_entry] = _incrlev;\ 591 }\ 592 }\ 593 } 594 595 /* 596 Copy data on the list into an array, then initialize the list 597 Input Parameters: 598 idx_start - starting index of the list 599 lnk_max - max value of lnk indicating the end of the list 600 nlnk - number of data on the list to be copied 601 lnk - linked list 602 lnklvl - level of lnk 603 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 604 output Parameters: 605 indices - array that contains the copied data 606 lnk - linked list that is cleaned and initialize 607 lnklvl - level of lnk that is reinitialized 608 bt - PetscBT (bitarray) with all bits set to false 609 */ 610 #define PetscIncompleteLLClean(idx_start,lnk_max,nlnk,lnk,lnklvl,indices,indiceslvl,bt) 0;\ 611 {\ 612 int _j,_idx=idx_start;\ 613 for (_j=0; _j<nlnk; _j++){\ 614 _idx = lnk[_idx];\ 615 *(indices+_j) = _idx;\ 616 *(indiceslvl+_j) = lnklvl[_idx];\ 617 lnklvl[_idx] = -1;\ 618 PetscBTClear(bt,_idx);\ 619 }\ 620 lnk[idx_start] = lnk_max;\ 621 } 622 /* 623 Free memories used by the list 624 */ 625 #define PetscIncompleteLLDestroy(lnk,bt) \ 626 (PetscFree(lnk) || PetscBTDestroy(bt) || (0)) 627 628 PETSC_EXTERN_CXX_END 629 #endif /* __PETSCSYS_H */ 630