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 Add a index set into a sorted linked list 373 Input Parameters: 374 nidx - number of input indices 375 indices - interger array used for storing column indices 376 level - level of fill, e.g., ICC(level) 377 indices_lvl - level of indices 378 idx_start - starting index of the list 379 lnk - linked list(an integer array) that is created 380 lnk_lvl - levels of lnk 381 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 382 output Parameters: 383 nlnk - number of newly added indices 384 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 385 lnk_lvl - levels of lnk 386 bt - updated PetscBT (bitarray) 387 */ 388 #define PetscIncompleteLLAdd(nidx,indices,level,indices_lvl,idx_start,nlnk,lnk,lnk_lvl,bt) 0;\ 389 {\ 390 int _k,_entry,_location,_lnkdata,_incrlev;\ 391 nlnk = 0;\ 392 _lnkdata = idx_start;\ 393 for (_k=0; _k<nidx; _k++){\ 394 _incrlev = indices_lvl[_k] + 1;\ 395 if (_incrlev > level) {_k++; continue;} \ 396 _entry = indices[_k];\ 397 if (!PetscBTLookupSet(bt,_entry)){ /* new entry */\ 398 /* search for insertion location */\ 399 if (_k && _entry < _lnkdata) _lnkdata = idx_start;\ 400 do {\ 401 _location = _lnkdata;\ 402 _lnkdata = lnk[_location];\ 403 } while (_entry > _lnkdata);\ 404 /* insertion location is found, add entry into lnk */\ 405 lnk[_location] = _entry;\ 406 lnk[_entry] = _lnkdata;\ 407 lnk_lvl[_entry] = _incrlev;\ 408 nlnk++;\ 409 _lnkdata = _entry; /* next search starts from here if next_entry > _entry */\ 410 } else { /* existing entry: update lnk_lvl */\ 411 if (lnk_lvl[_entry] > _incrlev) lnk_lvl[_entry] = _incrlev;\ 412 }\ 413 }\ 414 } 415 416 /* 417 Add a SORTED index set into a sorted linked list 418 Input Parameters: 419 nidx - number of input indices 420 indices - sorted interger array used for storing column indices 421 level - level of fill, e.g., ICC(level) 422 indices_lvl - level of indices 423 idx_start - starting index of the list 424 lnk - linked list(an integer array) that is created 425 lnk_lvl - levels of lnk 426 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 427 output Parameters: 428 nlnk - number of newly added indices 429 lnk - the sorted(increasing order) linked list containing new and non-redundate entries from indices 430 lnk_lvl - levels of lnk 431 bt - updated PetscBT (bitarray) 432 */ 433 #define PetscIncompleteLLAddSorted(nidx,indices,level,indices_lvl,idx_start,nlnk,lnk,lnk_lvl,bt) 0;\ 434 {\ 435 int _k,_entry,_location,_lnkdata,_incrlev;\ 436 nlnk = 0;\ 437 _lnkdata = idx_start;\ 438 for (_k=0; _k<nidx; _k++){\ 439 _incrlev = indices_lvl[_k] + 1;\ 440 if (_incrlev > level) {_k++; continue;} \ 441 _entry = indices[_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 lnk_lvl[_entry] = _incrlev;\ 452 nlnk++;\ 453 _lnkdata = _entry; /* next search starts from here */\ 454 } else { /* existing entry: update lnk_lvl */\ 455 if (lnk_lvl[_entry] > _incrlev) lnk_lvl[_entry] = _incrlev;\ 456 }\ 457 }\ 458 } 459 460 /* 461 Copy data on the list into an array, then initialize the list 462 Input Parameters: 463 idx_start - starting index of the list 464 lnk_max - max value of lnk indicating the end of the list 465 nlnk - number of data on the list to be copied 466 lnk - linked list 467 lnk_lvl - level of lnk 468 bt - PetscBT (bitarray), bt[idx]=true marks idx is in lnk 469 output Parameters: 470 indices - array that contains the copied data 471 lnk -llinked list that is cleaned and initialize 472 lnk_lvl - level of lnk that is reinitialized 473 bt - PetscBT (bitarray) with all bits set to false 474 */ 475 #define PetscIncompleteLLClean(idx_start,lnk_max,nlnk,lnk,lnk_lvl,indices,indices_lvl,bt) 0;\ 476 {\ 477 int _j,_idx=idx_start;\ 478 for (_j=0; _j<nlnk; _j++){\ 479 _idx = lnk[_idx];\ 480 *(indices+_j) = _idx;\ 481 *(indices_lvl+_j) = lnk_lvl[_idx];\ 482 lnk_lvl[_idx] = -1;\ 483 PetscBTClear(bt,_idx);\ 484 }\ 485 lnk[idx_start] = lnk_max;\ 486 } 487 /* 488 Free memories used by the list 489 */ 490 #define PetscIncompleteLLDestroy(lnk,bt) \ 491 (PetscFree(lnk) || PetscBTDestroy(bt) || (0)) 492 493 PETSC_EXTERN_CXX_END 494 #endif /* __PETSCSYS_H */ 495