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 _k=nidx;\ 233 while (_k){/* assume indices are almost in increasing order, starting from its end saves computation */\ 234 _entry = indices[--_k];\ 235 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 236 /* search for insertion location */\ 237 _lnkdata = idx_start;\ 238 do {\ 239 _location = _lnkdata;\ 240 _lnkdata = lnk[_location];\ 241 } while (_entry > _lnkdata);\ 242 /* insertion location is found, add entry into lnk */\ 243 lnk[_location] = _entry;\ 244 lnk[_entry] = _lnkdata;\ 245 nlnk++;\ 246 }\ 247 }\ 248 } 249 250 /* 251 Add a SORTED index set into a sorted linked list 252 Input Parameters: 253 nidx - number of input indices 254 indices - sorted interger array 255 idx_start - starting index of the list 256 lnk - linked list(an integer array) that is created 257 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 258 output Parameters: 259 nlnk - number of newly added indices 260 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 261 bt - updated PetscBT (bitarray) 262 */ 263 #define PetscLLAddSorted(nidx,indices,idx_start,nlnk,lnk,bt) 0;\ 264 {\ 265 int _k,_entry,_location,_lnkdata;\ 266 nlnk = 0;\ 267 _lnkdata = idx_start;\ 268 for (_k=0; _k<nidx; _k++){\ 269 _entry = indices[_k];\ 270 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 271 /* search for insertion location */\ 272 do {\ 273 _location = _lnkdata;\ 274 _lnkdata = lnk[_location];\ 275 } while (_entry > _lnkdata);\ 276 /* insertion location is found, add entry into lnk */\ 277 lnk[_location] = _entry;\ 278 lnk[_entry] = _lnkdata;\ 279 nlnk++;\ 280 _lnkdata = _entry; /* next search starts from here */\ 281 }\ 282 }\ 283 } 284 285 /* 286 Add a SORTED index set into a sorted linked list used for LUFactorSymbolic() 287 Same as PetscLLAddSorted() with an additional operation: 288 count the number of input indices that are no larger than 'diag' 289 Input Parameters: 290 nidx - number of input indices 291 indices - sorted interger array 292 idx_start - starting index of the list 293 lnk - linked list(an integer array) that is created 294 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 295 diag - index of the active row in LUFactorSymbolic 296 output Parameters: 297 nlnk - number of newly added indices 298 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 299 bt - updated PetscBT (bitarray) 300 nzbd - number of input indices that are no larger than 'diag' 301 */ 302 #define PetscLLAddSortedLU(nidx,indices,idx_start,nlnk,lnk,bt,diag,nzbd) 0;\ 303 {\ 304 int _k,_entry,_location,_lnkdata;\ 305 nlnk = 0;\ 306 _lnkdata = idx_start;\ 307 nzbd = 0;\ 308 for (_k=0; _k<nidx; _k++){\ 309 _entry = indices[_k];\ 310 if (_entry <= diag) nzbd++;\ 311 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 312 /* search for insertion location */\ 313 do {\ 314 _location = _lnkdata;\ 315 _lnkdata = lnk[_location];\ 316 } while (_entry > _lnkdata);\ 317 /* insertion location is found, add entry into lnk */\ 318 lnk[_location] = _entry;\ 319 lnk[_entry] = _lnkdata;\ 320 nlnk++;\ 321 _lnkdata = _entry; /* next search starts from here */\ 322 }\ 323 }\ 324 } 325 326 /* 327 Copy data on the list into an array, then initialize the list 328 Input Parameters: 329 idx_start - starting index of the list 330 lnk_max - max value of lnk indicating the end of the list 331 nlnk - number of data on the list to be copied 332 lnk - linked list 333 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 334 output Parameters: 335 indices - array that contains the copied data 336 lnk - linked list that is cleaned and initialize 337 bt - PetscBT (bitarray) with all bits set to false 338 */ 339 #define PetscLLClean(idx_start,lnk_max,nlnk,lnk,indices,bt) 0;\ 340 {\ 341 int _j,_idx=idx_start;\ 342 for (_j=0; _j<nlnk; _j++){\ 343 _idx = lnk[_idx];\ 344 *(indices+_j) = _idx;\ 345 PetscBTClear(bt,_idx);\ 346 }\ 347 lnk[idx_start] = lnk_max;\ 348 } 349 /* 350 Free memories used by the list 351 */ 352 #define PetscLLDestroy(lnk,bt) (PetscFree(lnk) || PetscBTDestroy(bt)) 353 354 /* Routines below are used for incomplete matrix factorization */ 355 /* 356 Create and initialize a linked list and its levels 357 Input Parameters: 358 idx_start - starting index of the list 359 lnk_max - max value of lnk indicating the end of the list 360 nlnk - max length of the list 361 Output Parameters: 362 lnk - list initialized 363 lnk_lvl - array of size nlnk for storing levels of lnk 364 bt - PetscBT (bitarray) with all bits set to false 365 */ 366 #define PetscIncompleteLLCreate(idx_start,lnk_max,nlnk,lnk,lnk_lvl,bt)\ 367 (PetscMalloc(2*nlnk*sizeof(PetscInt),&lnk) || PetscBTCreate(nlnk,bt) || PetscBTMemzero(nlnk,bt) || (lnk[idx_start] = lnk_max,lnk_lvl = lnk + nlnk,0)) 368 369 /* 370 Add a index set into a sorted linked list 371 Input Parameters: 372 nidx - number of input indices 373 indices - interger array used for storing column indices 374 level - level of fill, e.g., ICC(level) 375 indices_lvl - level of indices 376 idx_start - starting index of the list 377 lnk - linked list(an integer array) that is created 378 lnk_lvl - levels of lnk 379 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 380 output Parameters: 381 nlnk - number of newly added indices 382 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 383 lnk_lvl - levels of lnk 384 bt - updated PetscBT (bitarray) 385 */ 386 #define PetscIncompleteLLAdd(nidx,indices,level,indices_lvl,idx_start,nlnk,lnk,lnk_lvl,bt) 0;\ 387 {\ 388 int _k,_entry,_location,_lnkdata,_incrlev;\ 389 nlnk = 0;\ 390 _k = nidx;\ 391 while (_k){/* assume indices are almost in increasing order, starting from its end saves computation */\ 392 _incrlev = indices_lvl[_k-1] + 1;\ 393 if (_incrlev > level) { --_k; continue;} \ 394 _entry = indices[--_k];\ 395 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 396 /* search for insertion location */\ 397 _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 nlnk++;\ 406 lnk_lvl[_entry] = _incrlev;\ 407 } else { /* existing entry: update lnk_lvl */\ 408 if (lnk_lvl[_entry] > _incrlev) lnk_lvl[_entry] = _incrlev;\ 409 }\ 410 }\ 411 } 412 413 /* 414 Add a SORTED index set into a sorted linked list 415 Input Parameters: 416 nidx - number of input indices 417 indices - sorted interger array used for storing column indices 418 level - level of fill, e.g., ICC(level) 419 indices_lvl - level of indices 420 idx_start - starting index of the list 421 lnk - linked list(an integer array) that is created 422 lnk_lvl - levels of lnk 423 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 424 output Parameters: 425 nlnk - number of newly added indices 426 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 427 lnk_lvl - levels of lnk 428 bt - updated PetscBT (bitarray) 429 */ 430 #define PetscIncompleteLLAddSorted(nidx,indices,level,indices_lvl,idx_start,nlnk,lnk,lnk_lvl,bt) 0;\ 431 {\ 432 int _k,_entry,_location,_lnkdata,_incrlev;\ 433 nlnk = 0;\ 434 _lnkdata = idx_start;\ 435 for (_k=0; _k<nidx; _k++){\ 436 _incrlev = indices_lvl[_k] + 1;\ 437 if (_incrlev > level) {_k++; continue;} \ 438 _entry = indices[_k];\ 439 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 440 /* search for insertion location */\ 441 do {\ 442 _location = _lnkdata;\ 443 _lnkdata = lnk[_location];\ 444 } while (_entry > _lnkdata);\ 445 /* insertion location is found, add entry into lnk */\ 446 lnk[_location] = _entry;\ 447 lnk[_entry] = _lnkdata;\ 448 lnk_lvl[_entry] = _incrlev;\ 449 nlnk++;\ 450 _lnkdata = _entry; /* next search starts from here */\ 451 } else { /* existing entry: update lnk_lvl */\ 452 if (lnk_lvl[_entry] > _incrlev) lnk_lvl[_entry] = _incrlev;\ 453 }\ 454 }\ 455 } 456 457 /* 458 Copy data on the list into an array, then initialize the list 459 Input Parameters: 460 idx_start - starting index of the list 461 lnk_max - max value of lnk indicating the end of the list 462 nlnk - number of data on the list to be copied 463 lnk - linked list 464 lnk_lvl - level of lnk 465 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 466 output Parameters: 467 indices - array that contains the copied data 468 lnk -llinked list that is cleaned and initialize 469 lnk_lvl - level of lnk that is reinitialized 470 bt - PetscBT (bitarray) with all bits set to false 471 */ 472 #define PetscIncompleteLLClean(idx_start,lnk_max,nlnk,lnk,lnk_lvl,indices,indices_lvl,bt) 0;\ 473 {\ 474 int _j,_idx=idx_start;\ 475 for (_j=0; _j<nlnk; _j++){\ 476 _idx = lnk[_idx];\ 477 *(indices+_j) = _idx;\ 478 *(indices_lvl+_j) = lnk_lvl[_idx];\ 479 lnk_lvl[_idx] = -1;\ 480 PetscBTClear(bt,_idx);\ 481 }\ 482 lnk[idx_start] = lnk_max;\ 483 } 484 /* 485 Free memories used by the list 486 */ 487 #define PetscIncompleteLLDestroy(lnk,bt) \ 488 (PetscFree(lnk) || PetscBTDestroy(bt) || (0)) 489 490 PETSC_EXTERN_CXX_END 491 #endif /* __PETSCSYS_H */ 492