1 /* Define Feature test macros to make sure atoll is available (SVr4, POSIX.1-2001, 4.3BSD, C99), not in (C89 and POSIX.1-1996) */ 2 #define PETSC_DESIRE_FEATURE_TEST_MACROS /* for atoll() */ 3 4 /* 5 These routines simplify the use of command line, file options, etc., and are used to manipulate the options database. 6 This provides the low-level interface, the high level interface is in aoptions.c 7 8 Some routines use regular malloc and free because it cannot know what malloc is requested with the 9 options database until it has already processed the input. 10 */ 11 12 #include <petsc/private/petscimpl.h> /*I "petscsys.h" I*/ 13 #include <petscviewer.h> 14 #include <ctype.h> 15 #if defined(PETSC_HAVE_MALLOC_H) 16 #include <malloc.h> 17 #endif 18 #if defined(PETSC_HAVE_STRINGS_H) 19 #include <strings.h> /* strcasecmp */ 20 #endif 21 22 #if defined(PETSC_HAVE_STRCASECMP) 23 #define PetscOptNameCmp(a, b) strcasecmp(a, b) 24 #elif defined(PETSC_HAVE_STRICMP) 25 #define PetscOptNameCmp(a, b) stricmp(a, b) 26 #else 27 #define PetscOptNameCmp(a, b) Error_strcasecmp_not_found 28 #endif 29 30 #include <petsc/private/hashtable.h> 31 32 /* This assumes ASCII encoding and ignores locale settings */ 33 /* Using tolower() is about 2X slower in microbenchmarks */ 34 static inline int PetscToLower(int c) { 35 return ((c >= 'A') & (c <= 'Z')) ? c + 'a' - 'A' : c; 36 } 37 38 /* Bob Jenkins's one at a time hash function (case-insensitive) */ 39 static inline unsigned int PetscOptHash(const char key[]) { 40 unsigned int hash = 0; 41 while (*key) { 42 hash += PetscToLower(*key++); 43 hash += hash << 10; 44 hash ^= hash >> 6; 45 } 46 hash += hash << 3; 47 hash ^= hash >> 11; 48 hash += hash << 15; 49 return hash; 50 } 51 52 static inline int PetscOptEqual(const char a[], const char b[]) { 53 return !PetscOptNameCmp(a, b); 54 } 55 56 KHASH_INIT(HO, kh_cstr_t, int, 1, PetscOptHash, PetscOptEqual) 57 58 /* 59 This table holds all the options set by the user. For simplicity, we use a static size database 60 */ 61 #define MAXOPTNAME PETSC_MAX_OPTION_NAME 62 #define MAXOPTIONS 512 63 #define MAXALIASES 25 64 #define MAXPREFIXES 25 65 #define MAXOPTIONSMONITORS 5 66 67 struct _n_PetscOptions { 68 PetscOptions previous; 69 int N; /* number of options */ 70 char *names[MAXOPTIONS]; /* option names */ 71 char *values[MAXOPTIONS]; /* option values */ 72 PetscBool used[MAXOPTIONS]; /* flag option use */ 73 PetscBool precedentProcessed; 74 75 /* Hash table */ 76 khash_t(HO) * ht; 77 78 /* Prefixes */ 79 int prefixind; 80 int prefixstack[MAXPREFIXES]; 81 char prefix[MAXOPTNAME]; 82 83 /* Aliases */ 84 int Naliases; /* number or aliases */ 85 char *aliases1[MAXALIASES]; /* aliased */ 86 char *aliases2[MAXALIASES]; /* aliasee */ 87 88 /* Help */ 89 PetscBool help; /* flag whether "-help" is in the database */ 90 PetscBool help_intro; /* flag whether "-help intro" is in the database */ 91 92 /* Monitors */ 93 PetscBool monitorFromOptions, monitorCancel; 94 PetscErrorCode (*monitor[MAXOPTIONSMONITORS])(const char[], const char[], void *); /* returns control to user after */ 95 PetscErrorCode (*monitordestroy[MAXOPTIONSMONITORS])(void **); /* */ 96 void *monitorcontext[MAXOPTIONSMONITORS]; /* to pass arbitrary user data into monitor */ 97 PetscInt numbermonitors; /* to, for instance, detect options being set */ 98 }; 99 100 static PetscOptions defaultoptions = NULL; /* the options database routines query this object for options */ 101 102 /* list of options which preceed others, i.e., are processed in PetscOptionsProcessPrecedentFlags() */ 103 /* these options can only take boolean values, the code will crash if given a non-boolean value */ 104 static const char *precedentOptions[] = {"-petsc_ci", "-options_monitor", "-options_monitor_cancel", "-help", "-skip_petscrc"}; 105 enum PetscPrecedentOption { 106 PO_CI_ENABLE, 107 PO_OPTIONS_MONITOR, 108 PO_OPTIONS_MONITOR_CANCEL, 109 PO_HELP, 110 PO_SKIP_PETSCRC, 111 PO_NUM 112 }; 113 114 static PetscErrorCode PetscOptionsSetValue_Private(PetscOptions, const char[], const char[], int *); 115 116 /* 117 Options events monitor 118 */ 119 static PetscErrorCode PetscOptionsMonitor(PetscOptions options, const char name[], const char value[]) { 120 PetscFunctionBegin; 121 if (!value) value = ""; 122 if (options->monitorFromOptions) PetscCall(PetscOptionsMonitorDefault(name, value, NULL)); 123 for (PetscInt i = 0; i < options->numbermonitors; i++) PetscCall((*options->monitor[i])(name, value, options->monitorcontext[i])); 124 PetscFunctionReturn(0); 125 } 126 127 /*@ 128 PetscOptionsCreate - Creates an empty options database. 129 130 Logically collective 131 132 Output Parameter: 133 . options - Options database object 134 135 Level: advanced 136 137 Note: 138 Though PETSc has a concept of multiple options database the current code uses a single default `PetscOptions` object 139 140 Developer Notes: 141 We may want eventually to pass a `MPI_Comm` to determine the ownership of the object 142 143 This object never got developed after being introduced, it is not clear that supporting multiple `PetscOptions` objects is useful 144 145 .seealso: `PetscOptionsDestroy()`, `PetscOptionsPush()`, `PetscOptionsPop()`, `PetscOptionsInsert()`, `PetscOptionsSetValue()` 146 @*/ 147 PetscErrorCode PetscOptionsCreate(PetscOptions *options) { 148 PetscFunctionBegin; 149 PetscValidPointer(options, 1); 150 *options = (PetscOptions)calloc(1, sizeof(**options)); 151 PetscCheck(*options, PETSC_COMM_SELF, PETSC_ERR_MEM, "Failed to allocate the options database"); 152 PetscFunctionReturn(0); 153 } 154 155 /*@ 156 PetscOptionsDestroy - Destroys an option database. 157 158 Logically collective on whatever communicator was associated with the call to `PetscOptionsCreate()` 159 160 Input Parameter: 161 . options - the `PetscOptions` object 162 163 Level: advanced 164 165 .seealso: `PetscOptionsInsert()`, `PetscOptionsPush()`, `PetscOptionsPop()`, `PetscOptionsInsert()`, `PetscOptionsSetValue()` 166 @*/ 167 PetscErrorCode PetscOptionsDestroy(PetscOptions *options) { 168 PetscFunctionBegin; 169 if (!*options) PetscFunctionReturn(0); 170 PetscCheck(!(*options)->previous, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "You are destroying an option that has been used with PetscOptionsPush() but does not have a corresponding PetscOptionsPop()"); 171 PetscCall(PetscOptionsClear(*options)); 172 /* XXX what about monitors ? */ 173 free(*options); 174 *options = NULL; 175 PetscFunctionReturn(0); 176 } 177 178 /* 179 PetscOptionsCreateDefault - Creates the default global options database 180 */ 181 PetscErrorCode PetscOptionsCreateDefault(void) { 182 PetscFunctionBegin; 183 if (PetscUnlikely(!defaultoptions)) PetscCall(PetscOptionsCreate(&defaultoptions)); 184 PetscFunctionReturn(0); 185 } 186 187 /*@ 188 PetscOptionsPush - Push a new `PetscOptions` object as the default provider of options 189 Allows using different parts of a code to use different options databases 190 191 Logically Collective 192 193 Input Parameter: 194 . opt - the options obtained with `PetscOptionsCreate()` 195 196 Notes: 197 Use `PetscOptionsPop()` to return to the previous default options database 198 199 The collectivity of this routine is complex; only the MPI ranks that call this routine will 200 have the affect of these options. If some processes that create objects call this routine and others do 201 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 202 on different ranks. 203 204 Developer Note: 205 Though this functionality has been provided it has never been used in PETSc and might be removed. 206 207 Level: advanced 208 209 .seealso: `PetscOptionsPop()`, `PetscOptionsCreate()`, `PetscOptionsInsert()`, `PetscOptionsSetValue()`, `PetscOptionsLeft()` 210 @*/ 211 PetscErrorCode PetscOptionsPush(PetscOptions opt) { 212 PetscFunctionBegin; 213 PetscCall(PetscOptionsCreateDefault()); 214 opt->previous = defaultoptions; 215 defaultoptions = opt; 216 PetscFunctionReturn(0); 217 } 218 219 /*@ 220 PetscOptionsPop - Pop the most recent `PetscOptionsPush() `to return to the previous default options 221 222 Logically collective on whatever communicator was associated with the call to `PetscOptionsCreate()` 223 224 Level: advanced 225 226 .seealso: `PetscOptionsPop()`, `PetscOptionsCreate()`, `PetscOptionsInsert()`, `PetscOptionsSetValue()`, `PetscOptionsLeft()` 227 @*/ 228 PetscErrorCode PetscOptionsPop(void) { 229 PetscOptions current = defaultoptions; 230 231 PetscFunctionBegin; 232 PetscCheck(defaultoptions, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Missing default options"); 233 PetscCheck(defaultoptions->previous, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscOptionsPop() called too many times"); 234 defaultoptions = defaultoptions->previous; 235 current->previous = NULL; 236 PetscFunctionReturn(0); 237 } 238 239 /* 240 PetscOptionsDestroyDefault - Destroys the default global options database 241 */ 242 PetscErrorCode PetscOptionsDestroyDefault(void) { 243 PetscFunctionBegin; 244 if (!defaultoptions) PetscFunctionReturn(0); 245 /* Destroy any options that the user forgot to pop */ 246 while (defaultoptions->previous) { 247 PetscOptions tmp = defaultoptions; 248 249 PetscCall(PetscOptionsPop()); 250 PetscCall(PetscOptionsDestroy(&tmp)); 251 } 252 PetscCall(PetscOptionsDestroy(&defaultoptions)); 253 PetscFunctionReturn(0); 254 } 255 256 /*@C 257 PetscOptionsValidKey - PETSc Options database keys must begin with one or two dashes (-) followed by a letter. 258 259 Not collective 260 261 Input Parameter: 262 . key - string to check if valid 263 264 Output Parameter: 265 . valid - `PETSC_TRUE` if a valid key 266 267 Level: intermediate 268 @*/ 269 PetscErrorCode PetscOptionsValidKey(const char key[], PetscBool *valid) { 270 char *ptr; 271 272 PetscFunctionBegin; 273 if (key) PetscValidCharPointer(key, 1); 274 PetscValidBoolPointer(valid, 2); 275 *valid = PETSC_FALSE; 276 if (!key) PetscFunctionReturn(0); 277 if (key[0] != '-') PetscFunctionReturn(0); 278 if (key[1] == '-') key++; 279 if (!isalpha((int)key[1])) PetscFunctionReturn(0); 280 (void)strtod(key, &ptr); 281 if (ptr != key && !(*ptr == '_' || isalnum((int)*ptr))) PetscFunctionReturn(0); 282 *valid = PETSC_TRUE; 283 PetscFunctionReturn(0); 284 } 285 286 /*@C 287 PetscOptionsInsertString - Inserts options into the database from a string 288 289 Logically Collective 290 291 Input Parameters: 292 + options - options object 293 - in_str - string that contains options separated by blanks 294 295 Level: intermediate 296 297 The collectivity of this routine is complex; only the MPI ranks that call this routine will 298 have the affect of these options. If some processes that create objects call this routine and others do 299 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 300 on different ranks. 301 302 Contributed by Boyana Norris 303 304 .seealso: `PetscOptionsSetValue()`, `PetscOptionsView()`, `PetscOptionsHasName()`, `PetscOptionsGetInt()`, 305 `PetscOptionsGetReal()`, `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsBool()`, 306 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 307 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 308 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 309 `PetscOptionsFList()`, `PetscOptionsEList()`, `PetscOptionsInsertFile()` 310 @*/ 311 PetscErrorCode PetscOptionsInsertString(PetscOptions options, const char in_str[]) { 312 MPI_Comm comm = PETSC_COMM_SELF; 313 char *first, *second; 314 PetscToken token; 315 316 PetscFunctionBegin; 317 PetscCall(PetscTokenCreate(in_str, ' ', &token)); 318 PetscCall(PetscTokenFind(token, &first)); 319 while (first) { 320 PetscBool isfile, isfileyaml, isstringyaml, ispush, ispop, key; 321 PetscCall(PetscStrcasecmp(first, "-options_file", &isfile)); 322 PetscCall(PetscStrcasecmp(first, "-options_file_yaml", &isfileyaml)); 323 PetscCall(PetscStrcasecmp(first, "-options_string_yaml", &isstringyaml)); 324 PetscCall(PetscStrcasecmp(first, "-prefix_push", &ispush)); 325 PetscCall(PetscStrcasecmp(first, "-prefix_pop", &ispop)); 326 PetscCall(PetscOptionsValidKey(first, &key)); 327 if (!key) { 328 PetscCall(PetscTokenFind(token, &first)); 329 } else if (isfile) { 330 PetscCall(PetscTokenFind(token, &second)); 331 PetscCall(PetscOptionsInsertFile(comm, options, second, PETSC_TRUE)); 332 PetscCall(PetscTokenFind(token, &first)); 333 } else if (isfileyaml) { 334 PetscCall(PetscTokenFind(token, &second)); 335 PetscCall(PetscOptionsInsertFileYAML(comm, options, second, PETSC_TRUE)); 336 PetscCall(PetscTokenFind(token, &first)); 337 } else if (isstringyaml) { 338 PetscCall(PetscTokenFind(token, &second)); 339 PetscCall(PetscOptionsInsertStringYAML(options, second)); 340 PetscCall(PetscTokenFind(token, &first)); 341 } else if (ispush) { 342 PetscCall(PetscTokenFind(token, &second)); 343 PetscCall(PetscOptionsPrefixPush(options, second)); 344 PetscCall(PetscTokenFind(token, &first)); 345 } else if (ispop) { 346 PetscCall(PetscOptionsPrefixPop(options)); 347 PetscCall(PetscTokenFind(token, &first)); 348 } else { 349 PetscCall(PetscTokenFind(token, &second)); 350 PetscCall(PetscOptionsValidKey(second, &key)); 351 if (!key) { 352 PetscCall(PetscOptionsSetValue(options, first, second)); 353 PetscCall(PetscTokenFind(token, &first)); 354 } else { 355 PetscCall(PetscOptionsSetValue(options, first, NULL)); 356 first = second; 357 } 358 } 359 } 360 PetscCall(PetscTokenDestroy(&token)); 361 PetscFunctionReturn(0); 362 } 363 364 /* 365 Returns a line (ended by a \n, \r or null character of any length. Result should be freed with free() 366 */ 367 static char *Petscgetline(FILE *f) { 368 size_t size = 0; 369 size_t len = 0; 370 size_t last = 0; 371 char *buf = NULL; 372 373 if (feof(f)) return NULL; 374 do { 375 size += 1024; /* BUFSIZ is defined as "the optimal read size for this platform" */ 376 buf = (char *)realloc((void *)buf, size); /* realloc(NULL,n) is the same as malloc(n) */ 377 /* Actually do the read. Note that fgets puts a terminal '\0' on the 378 end of the string, so we make sure we overwrite this */ 379 if (!fgets(buf + len, 1024, f)) buf[len] = 0; 380 PetscStrlen(buf, &len); 381 last = len - 1; 382 } while (!feof(f) && buf[last] != '\n' && buf[last] != '\r'); 383 if (len) return buf; 384 free(buf); 385 return NULL; 386 } 387 388 static PetscErrorCode PetscOptionsFilename(MPI_Comm comm, const char file[], char filename[PETSC_MAX_PATH_LEN], PetscBool *yaml) { 389 char fname[PETSC_MAX_PATH_LEN + 8], path[PETSC_MAX_PATH_LEN + 8], *tail; 390 391 PetscFunctionBegin; 392 *yaml = PETSC_FALSE; 393 PetscCall(PetscStrreplace(comm, file, fname, sizeof(fname))); 394 PetscCall(PetscFixFilename(fname, path)); 395 PetscCall(PetscStrendswith(path, ":yaml", yaml)); 396 if (*yaml) { 397 PetscCall(PetscStrrchr(path, ':', &tail)); 398 tail[-1] = 0; /* remove ":yaml" suffix from path */ 399 } 400 PetscCall(PetscStrncpy(filename, path, PETSC_MAX_PATH_LEN)); 401 /* check for standard YAML and JSON filename extensions */ 402 if (!*yaml) PetscCall(PetscStrendswith(filename, ".yaml", yaml)); 403 if (!*yaml) PetscCall(PetscStrendswith(filename, ".yml", yaml)); 404 if (!*yaml) PetscCall(PetscStrendswith(filename, ".json", yaml)); 405 if (!*yaml) { /* check file contents */ 406 PetscMPIInt rank; 407 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 408 if (rank == 0) { 409 FILE *fh = fopen(filename, "r"); 410 if (fh) { 411 char buf[6] = ""; 412 if (fread(buf, 1, 6, fh) > 0) { 413 PetscCall(PetscStrncmp(buf, "%YAML ", 6, yaml)); /* check for '%YAML' tag */ 414 if (!*yaml) PetscCall(PetscStrncmp(buf, "---", 3, yaml)); /* check for document start */ 415 } 416 (void)fclose(fh); 417 } 418 } 419 PetscCallMPI(MPI_Bcast(yaml, 1, MPIU_BOOL, 0, comm)); 420 } 421 PetscFunctionReturn(0); 422 } 423 424 static PetscErrorCode PetscOptionsInsertFilePetsc(MPI_Comm comm, PetscOptions options, const char file[], PetscBool require) { 425 char *string, *vstring = NULL, *astring = NULL, *packed = NULL; 426 char *tokens[4]; 427 size_t i, len, bytes; 428 FILE *fd; 429 PetscToken token = NULL; 430 int err; 431 char *cmatch; 432 const char cmt = '#'; 433 PetscInt line = 1; 434 PetscMPIInt rank, cnt = 0, acnt = 0, counts[2]; 435 PetscBool isdir, alias = PETSC_FALSE, valid; 436 437 PetscFunctionBegin; 438 PetscCall(PetscMemzero(tokens, sizeof(tokens))); 439 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 440 if (rank == 0) { 441 char fpath[PETSC_MAX_PATH_LEN]; 442 char fname[PETSC_MAX_PATH_LEN]; 443 444 PetscCall(PetscStrreplace(PETSC_COMM_SELF, file, fpath, sizeof(fpath))); 445 PetscCall(PetscFixFilename(fpath, fname)); 446 447 fd = fopen(fname, "r"); 448 PetscCall(PetscTestDirectory(fname, 'r', &isdir)); 449 PetscCheck(!isdir || !require, PETSC_COMM_SELF, PETSC_ERR_USER, "Specified options file %s is a directory", fname); 450 if (fd && !isdir) { 451 PetscSegBuffer vseg, aseg; 452 PetscCall(PetscSegBufferCreate(1, 4000, &vseg)); 453 PetscCall(PetscSegBufferCreate(1, 2000, &aseg)); 454 455 /* the following line will not work when opening initial files (like .petscrc) since info is not yet set */ 456 PetscCall(PetscInfo(NULL, "Opened options file %s\n", file)); 457 458 while ((string = Petscgetline(fd))) { 459 /* eliminate comments from each line */ 460 PetscCall(PetscStrchr(string, cmt, &cmatch)); 461 if (cmatch) *cmatch = 0; 462 PetscCall(PetscStrlen(string, &len)); 463 /* replace tabs, ^M, \n with " " */ 464 for (i = 0; i < len; i++) { 465 if (string[i] == '\t' || string[i] == '\r' || string[i] == '\n') string[i] = ' '; 466 } 467 PetscCall(PetscTokenCreate(string, ' ', &token)); 468 PetscCall(PetscTokenFind(token, &tokens[0])); 469 if (!tokens[0]) { 470 goto destroy; 471 } else if (!tokens[0][0]) { /* if token 0 is empty (string begins with spaces), redo */ 472 PetscCall(PetscTokenFind(token, &tokens[0])); 473 } 474 for (i = 1; i < 4; i++) PetscCall(PetscTokenFind(token, &tokens[i])); 475 if (!tokens[0]) { 476 goto destroy; 477 } else if (tokens[0][0] == '-') { 478 PetscCall(PetscOptionsValidKey(tokens[0], &valid)); 479 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Error in options file %s line %" PetscInt_FMT ": invalid option %s", fname, line, tokens[0]); 480 PetscCall(PetscStrlen(tokens[0], &len)); 481 PetscCall(PetscSegBufferGet(vseg, len + 1, &vstring)); 482 PetscCall(PetscArraycpy(vstring, tokens[0], len)); 483 vstring[len] = ' '; 484 if (tokens[1]) { 485 PetscCall(PetscOptionsValidKey(tokens[1], &valid)); 486 PetscCheck(!valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Error in options file %s line %" PetscInt_FMT ": cannot specify two options per line (%s %s)", fname, line, tokens[0], tokens[1]); 487 PetscCall(PetscStrlen(tokens[1], &len)); 488 PetscCall(PetscSegBufferGet(vseg, len + 3, &vstring)); 489 vstring[0] = '"'; 490 PetscCall(PetscArraycpy(vstring + 1, tokens[1], len)); 491 vstring[len + 1] = '"'; 492 vstring[len + 2] = ' '; 493 } 494 } else { 495 PetscCall(PetscStrcasecmp(tokens[0], "alias", &alias)); 496 if (alias) { 497 PetscCall(PetscOptionsValidKey(tokens[1], &valid)); 498 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Error in options file %s line %" PetscInt_FMT ": invalid aliased option %s", fname, line, tokens[1]); 499 PetscCheck(tokens[2], PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Error in options file %s line %" PetscInt_FMT ": alias missing for %s", fname, line, tokens[1]); 500 PetscCall(PetscOptionsValidKey(tokens[2], &valid)); 501 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Error in options file %s line %" PetscInt_FMT ": invalid aliasee option %s", fname, line, tokens[2]); 502 PetscCall(PetscStrlen(tokens[1], &len)); 503 PetscCall(PetscSegBufferGet(aseg, len + 1, &astring)); 504 PetscCall(PetscArraycpy(astring, tokens[1], len)); 505 astring[len] = ' '; 506 507 PetscCall(PetscStrlen(tokens[2], &len)); 508 PetscCall(PetscSegBufferGet(aseg, len + 1, &astring)); 509 PetscCall(PetscArraycpy(astring, tokens[2], len)); 510 astring[len] = ' '; 511 } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown first token in options file %s line %" PetscInt_FMT ": %s", fname, line, tokens[0]); 512 } 513 { 514 const char *extraToken = alias ? tokens[3] : tokens[2]; 515 PetscCheck(!extraToken, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Error in options file %s line %" PetscInt_FMT ": extra token %s", fname, line, extraToken); 516 } 517 destroy: 518 free(string); 519 PetscCall(PetscTokenDestroy(&token)); 520 alias = PETSC_FALSE; 521 line++; 522 } 523 err = fclose(fd); 524 PetscCheck(!err, PETSC_COMM_SELF, PETSC_ERR_SYS, "fclose() failed on file %s", fname); 525 PetscCall(PetscSegBufferGetSize(aseg, &bytes)); /* size without null termination */ 526 PetscCall(PetscMPIIntCast(bytes, &acnt)); 527 PetscCall(PetscSegBufferGet(aseg, 1, &astring)); 528 astring[0] = 0; 529 PetscCall(PetscSegBufferGetSize(vseg, &bytes)); /* size without null termination */ 530 PetscCall(PetscMPIIntCast(bytes, &cnt)); 531 PetscCall(PetscSegBufferGet(vseg, 1, &vstring)); 532 vstring[0] = 0; 533 PetscCall(PetscMalloc1(2 + acnt + cnt, &packed)); 534 PetscCall(PetscSegBufferExtractTo(aseg, packed)); 535 PetscCall(PetscSegBufferExtractTo(vseg, packed + acnt + 1)); 536 PetscCall(PetscSegBufferDestroy(&aseg)); 537 PetscCall(PetscSegBufferDestroy(&vseg)); 538 } else PetscCheck(!require, PETSC_COMM_SELF, PETSC_ERR_USER, "Unable to open options file %s", fname); 539 } 540 541 counts[0] = acnt; 542 counts[1] = cnt; 543 err = MPI_Bcast(counts, 2, MPI_INT, 0, comm); 544 PetscCheck(!err, PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in first MPI collective call, could be caused by using an incorrect mpiexec or a network problem, it can be caused by having VPN running: see https://petsc.org/release/faq/"); 545 acnt = counts[0]; 546 cnt = counts[1]; 547 if (rank) PetscCall(PetscMalloc1(2 + acnt + cnt, &packed)); 548 if (acnt || cnt) { 549 PetscCallMPI(MPI_Bcast(packed, 2 + acnt + cnt, MPI_CHAR, 0, comm)); 550 astring = packed; 551 vstring = packed + acnt + 1; 552 } 553 554 if (acnt) { 555 PetscCall(PetscTokenCreate(astring, ' ', &token)); 556 PetscCall(PetscTokenFind(token, &tokens[0])); 557 while (tokens[0]) { 558 PetscCall(PetscTokenFind(token, &tokens[1])); 559 PetscCall(PetscOptionsSetAlias(options, tokens[0], tokens[1])); 560 PetscCall(PetscTokenFind(token, &tokens[0])); 561 } 562 PetscCall(PetscTokenDestroy(&token)); 563 } 564 565 if (cnt) PetscCall(PetscOptionsInsertString(options, vstring)); 566 PetscCall(PetscFree(packed)); 567 PetscFunctionReturn(0); 568 } 569 570 /*@C 571 PetscOptionsInsertFile - Inserts options into the database from a file. 572 573 Collective 574 575 Input Parameters: 576 + comm - the processes that will share the options (usually `PETSC_COMM_WORLD`) 577 . options - options database, use NULL for default global database 578 . file - name of file, 579 ".yml" and ".yaml" filename extensions are inserted as YAML options, 580 append ":yaml" to filename to force YAML options. 581 - require - if `PETSC_TRUE` will generate an error if the file does not exist 582 583 Notes: 584 Use # for lines that are comments and which should be ignored. 585 Usually, instead of using this command, one should list the file name in the call to `PetscInitialize()`, this insures that certain options 586 such as -log_view or -malloc_debug are processed properly. This routine only sets options into the options database that will be processed by later 587 calls to XXXSetFromOptions() it should not be used for options listed under PetscInitialize(). 588 The collectivity of this routine is complex; only the MPI ranks in comm will 589 have the affect of these options. If some processes that create objects call this routine and others do 590 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 591 on different ranks. 592 593 Level: developer 594 595 .seealso: `PetscOptionsSetValue()`, `PetscOptionsView()`, `PetscOptionsHasName()`, `PetscOptionsGetInt()`, 596 `PetscOptionsGetReal()`, `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsBool()`, 597 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 598 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 599 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 600 `PetscOptionsFList()`, `PetscOptionsEList()` 601 @*/ 602 PetscErrorCode PetscOptionsInsertFile(MPI_Comm comm, PetscOptions options, const char file[], PetscBool require) { 603 char filename[PETSC_MAX_PATH_LEN]; 604 PetscBool yaml; 605 606 PetscFunctionBegin; 607 PetscCall(PetscOptionsFilename(comm, file, filename, &yaml)); 608 if (yaml) { 609 PetscCall(PetscOptionsInsertFileYAML(comm, options, filename, require)); 610 } else { 611 PetscCall(PetscOptionsInsertFilePetsc(comm, options, filename, require)); 612 } 613 PetscFunctionReturn(0); 614 } 615 616 /*@C 617 PetscOptionsInsertArgs - Inserts options into the database from a array of strings 618 619 Logically Collective 620 621 Input Parameters: 622 + options - options object 623 . argc - the array length 624 - args - the string array 625 626 Level: intermediate 627 628 .seealso: `PetscOptions`, `PetscOptionsInsertString()`, `PetscOptionsInsertFile()` 629 @*/ 630 PetscErrorCode PetscOptionsInsertArgs(PetscOptions options, int argc, char *args[]) { 631 MPI_Comm comm = PETSC_COMM_WORLD; 632 int left = PetscMax(argc, 0); 633 char *const *eargs = args; 634 635 PetscFunctionBegin; 636 while (left) { 637 PetscBool isfile, isfileyaml, isstringyaml, ispush, ispop, key; 638 PetscCall(PetscStrcasecmp(eargs[0], "-options_file", &isfile)); 639 PetscCall(PetscStrcasecmp(eargs[0], "-options_file_yaml", &isfileyaml)); 640 PetscCall(PetscStrcasecmp(eargs[0], "-options_string_yaml", &isstringyaml)); 641 PetscCall(PetscStrcasecmp(eargs[0], "-prefix_push", &ispush)); 642 PetscCall(PetscStrcasecmp(eargs[0], "-prefix_pop", &ispop)); 643 PetscCall(PetscOptionsValidKey(eargs[0], &key)); 644 if (!key) { 645 eargs++; 646 left--; 647 } else if (isfile) { 648 PetscCheck(left > 1 && eargs[1][0] != '-', PETSC_COMM_SELF, PETSC_ERR_USER, "Missing filename for -options_file filename option"); 649 PetscCall(PetscOptionsInsertFile(comm, options, eargs[1], PETSC_TRUE)); 650 eargs += 2; 651 left -= 2; 652 } else if (isfileyaml) { 653 PetscCheck(left > 1 && eargs[1][0] != '-', PETSC_COMM_SELF, PETSC_ERR_USER, "Missing filename for -options_file_yaml filename option"); 654 PetscCall(PetscOptionsInsertFileYAML(comm, options, eargs[1], PETSC_TRUE)); 655 eargs += 2; 656 left -= 2; 657 } else if (isstringyaml) { 658 PetscCheck(left > 1 && eargs[1][0] != '-', PETSC_COMM_SELF, PETSC_ERR_USER, "Missing string for -options_string_yaml string option"); 659 PetscCall(PetscOptionsInsertStringYAML(options, eargs[1])); 660 eargs += 2; 661 left -= 2; 662 } else if (ispush) { 663 PetscCheck(left > 1, PETSC_COMM_SELF, PETSC_ERR_USER, "Missing prefix for -prefix_push option"); 664 PetscCheck(eargs[1][0] != '-', PETSC_COMM_SELF, PETSC_ERR_USER, "Missing prefix for -prefix_push option (prefixes cannot start with '-')"); 665 PetscCall(PetscOptionsPrefixPush(options, eargs[1])); 666 eargs += 2; 667 left -= 2; 668 } else if (ispop) { 669 PetscCall(PetscOptionsPrefixPop(options)); 670 eargs++; 671 left--; 672 } else { 673 PetscBool nextiskey = PETSC_FALSE; 674 if (left >= 2) PetscCall(PetscOptionsValidKey(eargs[1], &nextiskey)); 675 if (left < 2 || nextiskey) { 676 PetscCall(PetscOptionsSetValue(options, eargs[0], NULL)); 677 eargs++; 678 left--; 679 } else { 680 PetscCall(PetscOptionsSetValue(options, eargs[0], eargs[1])); 681 eargs += 2; 682 left -= 2; 683 } 684 } 685 } 686 PetscFunctionReturn(0); 687 } 688 689 static inline PetscErrorCode PetscOptionsStringToBoolIfSet_Private(enum PetscPrecedentOption opt, const char *val[], PetscBool set[], PetscBool *flg) { 690 PetscFunctionBegin; 691 if (set[opt]) { 692 PetscCall(PetscOptionsStringToBool(val[opt], flg)); 693 } else *flg = PETSC_FALSE; 694 PetscFunctionReturn(0); 695 } 696 697 /* Process options with absolute precedence, these are only processed from the command line, not the environment or files */ 698 static PetscErrorCode PetscOptionsProcessPrecedentFlags(PetscOptions options, int argc, char *args[], PetscBool *skip_petscrc, PetscBool *skip_petscrc_set) { 699 const char *const *opt = precedentOptions; 700 const size_t n = PO_NUM; 701 size_t o; 702 int a; 703 const char **val; 704 char **cval; 705 PetscBool *set, unneeded; 706 707 PetscFunctionBegin; 708 PetscCall(PetscCalloc2(n, &cval, n, &set)); 709 val = (const char **)cval; 710 711 /* Look for options possibly set using PetscOptionsSetValue beforehand */ 712 for (o = 0; o < n; o++) PetscCall(PetscOptionsFindPair(options, NULL, opt[o], &val[o], &set[o])); 713 714 /* Loop through all args to collect last occurring value of each option */ 715 for (a = 1; a < argc; a++) { 716 PetscBool valid, eq; 717 718 PetscCall(PetscOptionsValidKey(args[a], &valid)); 719 if (!valid) continue; 720 for (o = 0; o < n; o++) { 721 PetscCall(PetscStrcasecmp(args[a], opt[o], &eq)); 722 if (eq) { 723 set[o] = PETSC_TRUE; 724 if (a == argc - 1 || !args[a + 1] || !args[a + 1][0] || args[a + 1][0] == '-') val[o] = NULL; 725 else val[o] = args[a + 1]; 726 break; 727 } 728 } 729 } 730 731 /* Process flags */ 732 PetscCall(PetscStrcasecmp(val[PO_HELP], "intro", &options->help_intro)); 733 if (options->help_intro) options->help = PETSC_TRUE; 734 else PetscCall(PetscOptionsStringToBoolIfSet_Private(PO_HELP, val, set, &options->help)); 735 PetscCall(PetscOptionsStringToBoolIfSet_Private(PO_CI_ENABLE, val, set, &unneeded)); 736 /* need to manage PO_CI_ENABLE option before the PetscOptionsMonitor is turned on, so its setting is not monitored */ 737 if (set[PO_CI_ENABLE]) PetscCall(PetscOptionsSetValue_Private(options, opt[PO_CI_ENABLE], val[PO_CI_ENABLE], &a)); 738 PetscCall(PetscOptionsStringToBoolIfSet_Private(PO_OPTIONS_MONITOR_CANCEL, val, set, &options->monitorCancel)); 739 PetscCall(PetscOptionsStringToBoolIfSet_Private(PO_OPTIONS_MONITOR, val, set, &options->monitorFromOptions)); 740 PetscCall(PetscOptionsStringToBoolIfSet_Private(PO_SKIP_PETSCRC, val, set, skip_petscrc)); 741 *skip_petscrc_set = set[PO_SKIP_PETSCRC]; 742 743 /* Store precedent options in database and mark them as used */ 744 for (o = 1; o < n; o++) { 745 if (set[o]) { 746 PetscCall(PetscOptionsSetValue_Private(options, opt[o], val[o], &a)); 747 options->used[a] = PETSC_TRUE; 748 } 749 } 750 PetscCall(PetscFree2(cval, set)); 751 options->precedentProcessed = PETSC_TRUE; 752 PetscFunctionReturn(0); 753 } 754 755 static inline PetscErrorCode PetscOptionsSkipPrecedent(PetscOptions options, const char name[], PetscBool *flg) { 756 PetscFunctionBegin; 757 PetscValidBoolPointer(flg, 3); 758 *flg = PETSC_FALSE; 759 if (options->precedentProcessed) { 760 for (int i = 0; i < PO_NUM; ++i) { 761 if (!PetscOptNameCmp(precedentOptions[i], name)) { 762 /* check if precedent option has been set already */ 763 PetscCall(PetscOptionsFindPair(options, NULL, name, NULL, flg)); 764 if (*flg) break; 765 } 766 } 767 } 768 PetscFunctionReturn(0); 769 } 770 771 /*@C 772 PetscOptionsInsert - Inserts into the options database from the command line, 773 the environmental variable and a file. 774 775 Collective on `PETSC_COMM_WORLD` 776 777 Input Parameters: 778 + options - options database or NULL for the default global database 779 . argc - count of number of command line arguments 780 . args - the command line arguments 781 - file - [optional] PETSc database file, append ":yaml" to filename to specify YAML options format. 782 Use NULL or empty string to not check for code specific file. 783 Also checks ~/.petscrc, .petscrc and petscrc. 784 Use -skip_petscrc in the code specific file (or command line) to skip ~/.petscrc, .petscrc and petscrc files. 785 786 Options Database Keys: 787 + -options_file <filename> - read options from a file 788 - -options_file_yaml <filename> - read options from a YAML file 789 790 Notes: 791 Since PetscOptionsInsert() is automatically called by `PetscInitialize()`, 792 the user does not typically need to call this routine. `PetscOptionsInsert()` 793 can be called several times, adding additional entries into the database. 794 795 See `PetscInitialize()` for options related to option database monitoring. 796 797 Level: advanced 798 799 .seealso: `PetscOptionsDestroy()`, `PetscOptionsView()`, `PetscOptionsInsertString()`, `PetscOptionsInsertFile()`, 800 `PetscInitialize()` 801 @*/ 802 PetscErrorCode PetscOptionsInsert(PetscOptions options, int *argc, char ***args, const char file[]) { 803 MPI_Comm comm = PETSC_COMM_WORLD; 804 PetscMPIInt rank; 805 PetscBool hasArgs = (argc && *argc) ? PETSC_TRUE : PETSC_FALSE; 806 PetscBool skipPetscrc = PETSC_FALSE, skipPetscrcSet = PETSC_FALSE; 807 808 PetscFunctionBegin; 809 PetscCheck(!hasArgs || (args && *args), comm, PETSC_ERR_ARG_NULL, "*argc > 1 but *args not given"); 810 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 811 812 if (!options) { 813 PetscCall(PetscOptionsCreateDefault()); 814 options = defaultoptions; 815 } 816 if (hasArgs) { 817 /* process options with absolute precedence */ 818 PetscCall(PetscOptionsProcessPrecedentFlags(options, *argc, *args, &skipPetscrc, &skipPetscrcSet)); 819 PetscCall(PetscOptionsGetBool(NULL, NULL, "-petsc_ci", &PetscCIEnabled, NULL)); 820 } 821 if (file && file[0]) { 822 PetscCall(PetscOptionsInsertFile(comm, options, file, PETSC_TRUE)); 823 /* if -skip_petscrc has not been set from command line, check whether it has been set in the file */ 824 if (!skipPetscrcSet) PetscCall(PetscOptionsGetBool(options, NULL, "-skip_petscrc", &skipPetscrc, NULL)); 825 } 826 if (!skipPetscrc) { 827 char filename[PETSC_MAX_PATH_LEN]; 828 PetscCall(PetscGetHomeDirectory(filename, sizeof(filename))); 829 PetscCallMPI(MPI_Bcast(filename, (int)sizeof(filename), MPI_CHAR, 0, comm)); 830 if (filename[0]) PetscCall(PetscStrcat(filename, "/.petscrc")); 831 PetscCall(PetscOptionsInsertFile(comm, options, filename, PETSC_FALSE)); 832 PetscCall(PetscOptionsInsertFile(comm, options, ".petscrc", PETSC_FALSE)); 833 PetscCall(PetscOptionsInsertFile(comm, options, "petscrc", PETSC_FALSE)); 834 } 835 836 /* insert environment options */ 837 { 838 char *eoptions = NULL; 839 size_t len = 0; 840 if (rank == 0) { 841 eoptions = (char *)getenv("PETSC_OPTIONS"); 842 PetscCall(PetscStrlen(eoptions, &len)); 843 } 844 PetscCallMPI(MPI_Bcast(&len, 1, MPIU_SIZE_T, 0, comm)); 845 if (len) { 846 if (rank) PetscCall(PetscMalloc1(len + 1, &eoptions)); 847 PetscCallMPI(MPI_Bcast(eoptions, len, MPI_CHAR, 0, comm)); 848 if (rank) eoptions[len] = 0; 849 PetscCall(PetscOptionsInsertString(options, eoptions)); 850 if (rank) PetscCall(PetscFree(eoptions)); 851 } 852 } 853 854 /* insert YAML environment options */ 855 { 856 char *eoptions = NULL; 857 size_t len = 0; 858 if (rank == 0) { 859 eoptions = (char *)getenv("PETSC_OPTIONS_YAML"); 860 PetscCall(PetscStrlen(eoptions, &len)); 861 } 862 PetscCallMPI(MPI_Bcast(&len, 1, MPIU_SIZE_T, 0, comm)); 863 if (len) { 864 if (rank) PetscCall(PetscMalloc1(len + 1, &eoptions)); 865 PetscCallMPI(MPI_Bcast(eoptions, len, MPI_CHAR, 0, comm)); 866 if (rank) eoptions[len] = 0; 867 PetscCall(PetscOptionsInsertStringYAML(options, eoptions)); 868 if (rank) PetscCall(PetscFree(eoptions)); 869 } 870 } 871 872 /* insert command line options here because they take precedence over arguments in petscrc/environment */ 873 if (hasArgs) PetscCall(PetscOptionsInsertArgs(options, *argc - 1, *args + 1)); 874 PetscCall(PetscOptionsGetBool(NULL, NULL, "-petsc_ci_portable_error_output", &PetscCIEnabledPortableErrorOutput, NULL)); 875 PetscFunctionReturn(0); 876 } 877 878 /* These options are not printed with PetscOptionsView() or PetscOptionsMonitor() when PetscCIEnabled is on */ 879 /* TODO: get the list from the test harness, do not have it hardwired here. Maybe from gmakegentest.py */ 880 static const char *PetscCIOptions[] = { 881 "malloc_debug", 882 "malloc_dump", 883 "malloc_test", 884 "nox", 885 "nox_warning", 886 "display", 887 "saws_port_auto_select", 888 "saws_port_auto_select_silent", 889 "vecscatter_mpi1", 890 "check_pointer_intensity", 891 "cuda_initialize", 892 "error_output_stdout", 893 "use_gpu_aware_mpi", 894 "checkfunctionlist", 895 "petsc_ci", 896 "petsc_ci_portable_error_output", 897 }; 898 899 static PetscBool PetscCIOption(const char *name) { 900 PetscInt idx; 901 PetscBool found; 902 903 if (!PetscCIEnabled) return PETSC_FALSE; 904 PetscEListFind(PETSC_STATIC_ARRAY_LENGTH(PetscCIOptions), PetscCIOptions, name, &idx, &found); 905 return found; 906 } 907 908 /*@C 909 PetscOptionsView - Prints the options that have been loaded. This is 910 useful for debugging purposes. 911 912 Logically Collective on viewer 913 914 Input Parameters: 915 + options - options database, use NULL for default global database 916 - viewer - must be an `PETSCVIEWERASCII` viewer 917 918 Options Database Key: 919 . -options_view - Activates `PetscOptionsView()` within `PetscFinalize()` 920 921 Note: 922 Only the rank zero process of the `MPI_Comm` used to create view prints the option values. Other processes 923 may have different values but they are not printed. 924 925 Level: advanced 926 927 .seealso: `PetscOptionsAllUsed()` 928 @*/ 929 PetscErrorCode PetscOptionsView(PetscOptions options, PetscViewer viewer) { 930 PetscInt i, N = 0; 931 PetscBool isascii; 932 933 PetscFunctionBegin; 934 if (viewer) PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 2); 935 options = options ? options : defaultoptions; 936 if (!viewer) viewer = PETSC_VIEWER_STDOUT_WORLD; 937 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii)); 938 PetscCheck(isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Only supports ASCII viewer"); 939 940 for (i = 0; i < options->N; i++) { 941 if (PetscCIOption(options->names[i])) continue; 942 N++; 943 } 944 945 if (!N) { 946 PetscCall(PetscViewerASCIIPrintf(viewer, "#No PETSc Option Table entries\n")); 947 PetscFunctionReturn(0); 948 } 949 950 PetscCall(PetscViewerASCIIPrintf(viewer, "#PETSc Option Table entries:\n")); 951 for (i = 0; i < options->N; i++) { 952 if (PetscCIOption(options->names[i])) continue; 953 if (options->values[i]) { 954 PetscCall(PetscViewerASCIIPrintf(viewer, "-%s %s\n", options->names[i], options->values[i])); 955 } else { 956 PetscCall(PetscViewerASCIIPrintf(viewer, "-%s\n", options->names[i])); 957 } 958 } 959 PetscCall(PetscViewerASCIIPrintf(viewer, "#End of PETSc Option Table entries\n")); 960 PetscFunctionReturn(0); 961 } 962 963 /* 964 Called by error handlers to print options used in run 965 */ 966 PetscErrorCode PetscOptionsLeftError(void) { 967 PetscInt i, nopt = 0; 968 969 for (i = 0; i < defaultoptions->N; i++) { 970 if (!defaultoptions->used[i]) { 971 if (PetscCIOption(defaultoptions->names[i])) continue; 972 nopt++; 973 } 974 } 975 if (nopt) { 976 (*PetscErrorPrintf)("WARNING! There are option(s) set that were not used! Could be the program crashed before they were used or a spelling mistake, etc!\n"); 977 for (i = 0; i < defaultoptions->N; i++) { 978 if (!defaultoptions->used[i]) { 979 if (PetscCIOption(defaultoptions->names[i])) continue; 980 if (defaultoptions->values[i]) (*PetscErrorPrintf)("Option left: name:-%s value: %s\n", defaultoptions->names[i], defaultoptions->values[i]); 981 else (*PetscErrorPrintf)("Option left: name:-%s (no value)\n", defaultoptions->names[i]); 982 } 983 } 984 } 985 return 0; 986 } 987 988 PETSC_EXTERN PetscErrorCode PetscOptionsViewError(void) { 989 PetscInt i, N = 0; 990 PetscOptions options = defaultoptions; 991 992 for (i = 0; i < options->N; i++) { 993 if (PetscCIOption(options->names[i])) continue; 994 N++; 995 } 996 997 if (N) { 998 (*PetscErrorPrintf)("PETSc Option Table entries:\n"); 999 } else { 1000 (*PetscErrorPrintf)("No PETSc Option Table entries\n"); 1001 } 1002 for (i = 0; i < options->N; i++) { 1003 if (PetscCIOption(options->names[i])) continue; 1004 if (options->values[i]) { 1005 (*PetscErrorPrintf)("-%s %s\n", options->names[i], options->values[i]); 1006 } else { 1007 (*PetscErrorPrintf)("-%s\n", options->names[i]); 1008 } 1009 } 1010 return 0; 1011 } 1012 1013 /*@C 1014 PetscOptionsPrefixPush - Designate a prefix to be used by all options insertions to follow. 1015 1016 Logically Collective 1017 1018 Input Parameters: 1019 + options - options database, or NULL for the default global database 1020 - prefix - The string to append to the existing prefix 1021 1022 Options Database Keys: 1023 + -prefix_push <some_prefix_> - push the given prefix 1024 - -prefix_pop - pop the last prefix 1025 1026 Notes: 1027 It is common to use this in conjunction with -options_file as in 1028 1029 $ -prefix_push system1_ -options_file system1rc -prefix_pop -prefix_push system2_ -options_file system2rc -prefix_pop 1030 1031 where the files no longer require all options to be prefixed with -system2_. 1032 1033 The collectivity of this routine is complex; only the MPI ranks that call this routine will 1034 have the affect of these options. If some processes that create objects call this routine and others do 1035 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1036 on different ranks. 1037 1038 Level: advanced 1039 1040 .seealso: `PetscOptionsPrefixPop()`, `PetscOptionsPush()`, `PetscOptionsPop()`, `PetscOptionsCreate()`, `PetscOptionsSetValue()` 1041 @*/ 1042 PetscErrorCode PetscOptionsPrefixPush(PetscOptions options, const char prefix[]) { 1043 size_t n; 1044 PetscInt start; 1045 char key[MAXOPTNAME + 1]; 1046 PetscBool valid; 1047 1048 PetscFunctionBegin; 1049 PetscValidCharPointer(prefix, 2); 1050 options = options ? options : defaultoptions; 1051 PetscCheck(options->prefixind < MAXPREFIXES, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Maximum depth of prefix stack %d exceeded, recompile \n src/sys/objects/options.c with larger value for MAXPREFIXES", MAXPREFIXES); 1052 key[0] = '-'; /* keys must start with '-' */ 1053 PetscCall(PetscStrncpy(key + 1, prefix, sizeof(key) - 1)); 1054 PetscCall(PetscOptionsValidKey(key, &valid)); 1055 if (!valid && options->prefixind > 0 && isdigit((int)prefix[0])) valid = PETSC_TRUE; /* If the prefix stack is not empty, make numbers a valid prefix */ 1056 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_USER, "Given prefix \"%s\" not valid (the first character must be a letter%s, do not include leading '-')", prefix, options->prefixind ? " or digit" : ""); 1057 start = options->prefixind ? options->prefixstack[options->prefixind - 1] : 0; 1058 PetscCall(PetscStrlen(prefix, &n)); 1059 PetscCheck(n + 1 <= sizeof(options->prefix) - start, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Maximum prefix length %zu exceeded", sizeof(options->prefix)); 1060 PetscCall(PetscArraycpy(options->prefix + start, prefix, n + 1)); 1061 options->prefixstack[options->prefixind++] = start + n; 1062 PetscFunctionReturn(0); 1063 } 1064 1065 /*@C 1066 PetscOptionsPrefixPop - Remove the latest options prefix, see `PetscOptionsPrefixPush()` for details 1067 1068 Logically Collective on the `MPI_Comm` used when called `PetscOptionsPrefixPush()` 1069 1070 Input Parameter: 1071 . options - options database, or NULL for the default global database 1072 1073 Level: advanced 1074 1075 .seealso: `PetscOptionsPrefixPush()`, `PetscOptionsPush()`, `PetscOptionsPop()`, `PetscOptionsCreate()`, `PetscOptionsSetValue()` 1076 @*/ 1077 PetscErrorCode PetscOptionsPrefixPop(PetscOptions options) { 1078 PetscInt offset; 1079 1080 PetscFunctionBegin; 1081 options = options ? options : defaultoptions; 1082 PetscCheck(options->prefixind >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "More prefixes popped than pushed"); 1083 options->prefixind--; 1084 offset = options->prefixind ? options->prefixstack[options->prefixind - 1] : 0; 1085 options->prefix[offset] = 0; 1086 PetscFunctionReturn(0); 1087 } 1088 1089 /*@C 1090 PetscOptionsClear - Removes all options form the database leaving it empty. 1091 1092 Logically Collective 1093 1094 Input Parameter: 1095 . options - options database, use NULL for the default global database 1096 1097 The collectivity of this routine is complex; only the MPI ranks that call this routine will 1098 have the affect of these options. If some processes that create objects call this routine and others do 1099 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1100 on different ranks. 1101 1102 Level: developer 1103 1104 .seealso: `PetscOptionsInsert()` 1105 @*/ 1106 PetscErrorCode PetscOptionsClear(PetscOptions options) { 1107 PetscInt i; 1108 1109 PetscFunctionBegin; 1110 options = options ? options : defaultoptions; 1111 if (!options) PetscFunctionReturn(0); 1112 1113 for (i = 0; i < options->N; i++) { 1114 if (options->names[i]) free(options->names[i]); 1115 if (options->values[i]) free(options->values[i]); 1116 } 1117 options->N = 0; 1118 1119 for (i = 0; i < options->Naliases; i++) { 1120 free(options->aliases1[i]); 1121 free(options->aliases2[i]); 1122 } 1123 options->Naliases = 0; 1124 1125 /* destroy hash table */ 1126 kh_destroy(HO, options->ht); 1127 options->ht = NULL; 1128 1129 options->prefixind = 0; 1130 options->prefix[0] = 0; 1131 options->help = PETSC_FALSE; 1132 PetscFunctionReturn(0); 1133 } 1134 1135 /*@C 1136 PetscOptionsSetAlias - Makes a key and alias for another key 1137 1138 Logically Collective 1139 1140 Input Parameters: 1141 + options - options database, or NULL for default global database 1142 . newname - the alias 1143 - oldname - the name that alias will refer to 1144 1145 Level: advanced 1146 1147 The collectivity of this routine is complex; only the MPI ranks that call this routine will 1148 have the affect of these options. If some processes that create objects call this routine and others do 1149 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1150 on different ranks. 1151 1152 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, `OptionsHasName()`, 1153 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 1154 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 1155 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 1156 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 1157 `PetscOptionsFList()`, `PetscOptionsEList()` 1158 @*/ 1159 PetscErrorCode PetscOptionsSetAlias(PetscOptions options, const char newname[], const char oldname[]) { 1160 PetscInt n; 1161 size_t len; 1162 PetscBool valid; 1163 1164 PetscFunctionBegin; 1165 PetscValidCharPointer(newname, 2); 1166 PetscValidCharPointer(oldname, 3); 1167 options = options ? options : defaultoptions; 1168 PetscCall(PetscOptionsValidKey(newname, &valid)); 1169 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid aliased option %s", newname); 1170 PetscCall(PetscOptionsValidKey(oldname, &valid)); 1171 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid aliasee option %s", oldname); 1172 1173 n = options->Naliases; 1174 PetscCheck(n < MAXALIASES, PETSC_COMM_SELF, PETSC_ERR_MEM, "You have defined to many PETSc options aliases, limit %d recompile \n src/sys/objects/options.c with larger value for MAXALIASES", MAXALIASES); 1175 1176 newname++; 1177 oldname++; 1178 PetscCall(PetscStrlen(newname, &len)); 1179 options->aliases1[n] = (char *)malloc((len + 1) * sizeof(char)); 1180 PetscCall(PetscStrcpy(options->aliases1[n], newname)); 1181 PetscCall(PetscStrlen(oldname, &len)); 1182 options->aliases2[n] = (char *)malloc((len + 1) * sizeof(char)); 1183 PetscCall(PetscStrcpy(options->aliases2[n], oldname)); 1184 options->Naliases++; 1185 PetscFunctionReturn(0); 1186 } 1187 1188 /*@C 1189 PetscOptionsSetValue - Sets an option name-value pair in the options 1190 database, overriding whatever is already present. 1191 1192 Logically Collective 1193 1194 Input Parameters: 1195 + options - options database, use NULL for the default global database 1196 . name - name of option, this SHOULD have the - prepended 1197 - value - the option value (not used for all options, so can be NULL) 1198 1199 Level: intermediate 1200 1201 Note: 1202 This function can be called BEFORE `PetscInitialize()` 1203 1204 The collectivity of this routine is complex; only the MPI ranks that call this routine will 1205 have the affect of these options. If some processes that create objects call this routine and others do 1206 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1207 on different ranks. 1208 1209 Developers Note: Uses malloc() directly because PETSc may not be initialized yet. 1210 1211 .seealso: `PetscOptionsInsert()`, `PetscOptionsClearValue()` 1212 @*/ 1213 PetscErrorCode PetscOptionsSetValue(PetscOptions options, const char name[], const char value[]) { 1214 PetscFunctionBegin; 1215 PetscCall(PetscOptionsSetValue_Private(options, name, value, NULL)); 1216 PetscFunctionReturn(0); 1217 } 1218 1219 static PetscErrorCode PetscOptionsSetValue_Private(PetscOptions options, const char name[], const char value[], int *pos) { 1220 size_t len; 1221 int N, n, i; 1222 char **names; 1223 char fullname[MAXOPTNAME] = ""; 1224 PetscBool flg; 1225 1226 PetscFunctionBegin; 1227 if (!options) { 1228 PetscCall(PetscOptionsCreateDefault()); 1229 options = defaultoptions; 1230 } 1231 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "name %s must start with '-'", name); 1232 1233 PetscCall(PetscOptionsSkipPrecedent(options, name, &flg)); 1234 if (flg) PetscFunctionReturn(0); 1235 1236 name++; /* skip starting dash */ 1237 1238 if (options->prefixind > 0) { 1239 strncpy(fullname, options->prefix, sizeof(fullname)); 1240 fullname[sizeof(fullname) - 1] = 0; 1241 strncat(fullname, name, sizeof(fullname) - strlen(fullname) - 1); 1242 fullname[sizeof(fullname) - 1] = 0; 1243 name = fullname; 1244 } 1245 1246 /* check against aliases */ 1247 N = options->Naliases; 1248 for (i = 0; i < N; i++) { 1249 int result = PetscOptNameCmp(options->aliases1[i], name); 1250 if (!result) { 1251 name = options->aliases2[i]; 1252 break; 1253 } 1254 } 1255 1256 /* slow search */ 1257 N = n = options->N; 1258 names = options->names; 1259 for (i = 0; i < N; i++) { 1260 int result = PetscOptNameCmp(names[i], name); 1261 if (!result) { 1262 n = i; 1263 goto setvalue; 1264 } else if (result > 0) { 1265 n = i; 1266 break; 1267 } 1268 } 1269 PetscCheck(N < MAXOPTIONS, PETSC_COMM_SELF, PETSC_ERR_MEM, "Number of options %d < max number of options %d, can not allocate enough space", N, MAXOPTIONS); 1270 1271 /* shift remaining values up 1 */ 1272 for (i = N; i > n; i--) { 1273 options->names[i] = options->names[i - 1]; 1274 options->values[i] = options->values[i - 1]; 1275 options->used[i] = options->used[i - 1]; 1276 } 1277 options->names[n] = NULL; 1278 options->values[n] = NULL; 1279 options->used[n] = PETSC_FALSE; 1280 options->N++; 1281 1282 /* destroy hash table */ 1283 kh_destroy(HO, options->ht); 1284 options->ht = NULL; 1285 1286 /* set new name */ 1287 len = strlen(name); 1288 options->names[n] = (char *)malloc((len + 1) * sizeof(char)); 1289 PetscCheck(options->names[n], PETSC_COMM_SELF, PETSC_ERR_MEM, "Failed to allocate option name"); 1290 strcpy(options->names[n], name); 1291 1292 setvalue: 1293 /* set new value */ 1294 if (options->values[n]) free(options->values[n]); 1295 len = value ? strlen(value) : 0; 1296 if (len) { 1297 options->values[n] = (char *)malloc((len + 1) * sizeof(char)); 1298 if (!options->values[n]) return PETSC_ERR_MEM; 1299 strcpy(options->values[n], value); 1300 } else { 1301 options->values[n] = NULL; 1302 } 1303 1304 /* handle -help so that it can be set from anywhere */ 1305 if (!PetscOptNameCmp(name, "help")) { 1306 options->help = PETSC_TRUE; 1307 options->help_intro = (value && !PetscOptNameCmp(value, "intro")) ? PETSC_TRUE : PETSC_FALSE; 1308 options->used[n] = PETSC_TRUE; 1309 } 1310 1311 PetscCall(PetscOptionsMonitor(options, name, value)); 1312 if (pos) *pos = n; 1313 PetscFunctionReturn(0); 1314 } 1315 1316 /*@C 1317 PetscOptionsClearValue - Clears an option name-value pair in the options 1318 database, overriding whatever is already present. 1319 1320 Logically Collective 1321 1322 Input Parameters: 1323 + options - options database, use NULL for the default global database 1324 - name - name of option, this SHOULD have the - prepended 1325 1326 Level: intermediate 1327 1328 Note: 1329 The collectivity of this routine is complex; only the MPI ranks that call this routine will 1330 have the affect of these options. If some processes that create objects call this routine and others do 1331 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1332 on different ranks. 1333 1334 .seealso: `PetscOptionsInsert()` 1335 @*/ 1336 PetscErrorCode PetscOptionsClearValue(PetscOptions options, const char name[]) { 1337 int N, n, i; 1338 char **names; 1339 1340 PetscFunctionBegin; 1341 options = options ? options : defaultoptions; 1342 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Name must begin with '-': Instead %s", name); 1343 if (!PetscOptNameCmp(name, "-help")) options->help = options->help_intro = PETSC_FALSE; 1344 1345 name++; /* skip starting dash */ 1346 1347 /* slow search */ 1348 N = n = options->N; 1349 names = options->names; 1350 for (i = 0; i < N; i++) { 1351 int result = PetscOptNameCmp(names[i], name); 1352 if (!result) { 1353 n = i; 1354 break; 1355 } else if (result > 0) { 1356 n = N; 1357 break; 1358 } 1359 } 1360 if (n == N) PetscFunctionReturn(0); /* it was not present */ 1361 1362 /* remove name and value */ 1363 if (options->names[n]) free(options->names[n]); 1364 if (options->values[n]) free(options->values[n]); 1365 /* shift remaining values down 1 */ 1366 for (i = n; i < N - 1; i++) { 1367 options->names[i] = options->names[i + 1]; 1368 options->values[i] = options->values[i + 1]; 1369 options->used[i] = options->used[i + 1]; 1370 } 1371 options->N--; 1372 1373 /* destroy hash table */ 1374 kh_destroy(HO, options->ht); 1375 options->ht = NULL; 1376 1377 PetscCall(PetscOptionsMonitor(options, name, NULL)); 1378 PetscFunctionReturn(0); 1379 } 1380 1381 /*@C 1382 PetscOptionsFindPair - Gets an option name-value pair from the options database. 1383 1384 Not Collective 1385 1386 Input Parameters: 1387 + options - options database, use NULL for the default global database 1388 . pre - the string to prepend to the name or NULL, this SHOULD NOT have the "-" prepended 1389 - name - name of option, this SHOULD have the "-" prepended 1390 1391 Output Parameters: 1392 + value - the option value (optional, not used for all options) 1393 - set - whether the option is set (optional) 1394 1395 Note: 1396 Each process may find different values or no value depending on how options were inserted into the database 1397 1398 Level: developer 1399 1400 .seealso: `PetscOptionsSetValue()`, `PetscOptionsClearValue()` 1401 @*/ 1402 PetscErrorCode PetscOptionsFindPair(PetscOptions options, const char pre[], const char name[], const char *value[], PetscBool *set) { 1403 char buf[MAXOPTNAME]; 1404 PetscBool usehashtable = PETSC_TRUE; 1405 PetscBool matchnumbers = PETSC_TRUE; 1406 1407 PetscFunctionBegin; 1408 options = options ? options : defaultoptions; 1409 PetscCheck(!pre || !PetscUnlikely(pre[0] == '-'), PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Prefix cannot begin with '-': Instead %s", pre); 1410 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Name must begin with '-': Instead %s", name); 1411 1412 name++; /* skip starting dash */ 1413 1414 /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */ 1415 if (pre && pre[0]) { 1416 char *ptr = buf; 1417 if (name[0] == '-') { 1418 *ptr++ = '-'; 1419 name++; 1420 } 1421 PetscCall(PetscStrncpy(ptr, pre, buf + sizeof(buf) - ptr)); 1422 PetscCall(PetscStrlcat(buf, name, sizeof(buf))); 1423 name = buf; 1424 } 1425 1426 if (PetscDefined(USE_DEBUG)) { 1427 PetscBool valid; 1428 char key[MAXOPTNAME + 1] = "-"; 1429 PetscCall(PetscStrncpy(key + 1, name, sizeof(key) - 1)); 1430 PetscCall(PetscOptionsValidKey(key, &valid)); 1431 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid option '%s' obtained from pre='%s' and name='%s'", key, pre ? pre : "", name); 1432 } 1433 1434 if (!options->ht && usehashtable) { 1435 int i, ret; 1436 khiter_t it; 1437 khash_t(HO) * ht; 1438 ht = kh_init(HO); 1439 PetscCheck(ht, PETSC_COMM_SELF, PETSC_ERR_MEM, "Hash table allocation failed"); 1440 ret = kh_resize(HO, ht, options->N * 2); /* twice the required size to reduce risk of collisions */ 1441 PetscCheck(!ret, PETSC_COMM_SELF, PETSC_ERR_MEM, "Hash table allocation failed"); 1442 for (i = 0; i < options->N; i++) { 1443 it = kh_put(HO, ht, options->names[i], &ret); 1444 PetscCheck(ret == 1, PETSC_COMM_SELF, PETSC_ERR_MEM, "Hash table allocation failed"); 1445 kh_val(ht, it) = i; 1446 } 1447 options->ht = ht; 1448 } 1449 1450 if (usehashtable) { /* fast search */ 1451 khash_t(HO) *ht = options->ht; 1452 khiter_t it = kh_get(HO, ht, name); 1453 if (it != kh_end(ht)) { 1454 int i = kh_val(ht, it); 1455 options->used[i] = PETSC_TRUE; 1456 if (value) *value = options->values[i]; 1457 if (set) *set = PETSC_TRUE; 1458 PetscFunctionReturn(0); 1459 } 1460 } else { /* slow search */ 1461 int i, N = options->N; 1462 for (i = 0; i < N; i++) { 1463 int result = PetscOptNameCmp(options->names[i], name); 1464 if (!result) { 1465 options->used[i] = PETSC_TRUE; 1466 if (value) *value = options->values[i]; 1467 if (set) *set = PETSC_TRUE; 1468 PetscFunctionReturn(0); 1469 } else if (result > 0) { 1470 break; 1471 } 1472 } 1473 } 1474 1475 /* 1476 The following block slows down all lookups in the most frequent path (most lookups are unsuccessful). 1477 Maybe this special lookup mode should be enabled on request with a push/pop API. 1478 The feature of matching _%d_ used sparingly in the codebase. 1479 */ 1480 if (matchnumbers) { 1481 int i, j, cnt = 0, locs[16], loce[16]; 1482 /* determine the location and number of all _%d_ in the key */ 1483 for (i = 0; name[i]; i++) { 1484 if (name[i] == '_') { 1485 for (j = i + 1; name[j]; j++) { 1486 if (name[j] >= '0' && name[j] <= '9') continue; 1487 if (name[j] == '_' && j > i + 1) { /* found a number */ 1488 locs[cnt] = i + 1; 1489 loce[cnt++] = j + 1; 1490 } 1491 i = j - 1; 1492 break; 1493 } 1494 } 1495 } 1496 for (i = 0; i < cnt; i++) { 1497 PetscBool found; 1498 char opt[MAXOPTNAME + 1] = "-", tmp[MAXOPTNAME]; 1499 PetscCall(PetscStrncpy(tmp, name, PetscMin((size_t)(locs[i] + 1), sizeof(tmp)))); 1500 PetscCall(PetscStrlcat(opt, tmp, sizeof(opt))); 1501 PetscCall(PetscStrlcat(opt, name + loce[i], sizeof(opt))); 1502 PetscCall(PetscOptionsFindPair(options, NULL, opt, value, &found)); 1503 if (found) { 1504 if (set) *set = PETSC_TRUE; 1505 PetscFunctionReturn(0); 1506 } 1507 } 1508 } 1509 1510 if (set) *set = PETSC_FALSE; 1511 PetscFunctionReturn(0); 1512 } 1513 1514 /* Check whether any option begins with pre+name */ 1515 PETSC_EXTERN PetscErrorCode PetscOptionsFindPairPrefix_Private(PetscOptions options, const char pre[], const char name[], const char *value[], PetscBool *set) { 1516 char buf[MAXOPTNAME]; 1517 int numCnt = 0, locs[16], loce[16]; 1518 1519 PetscFunctionBegin; 1520 options = options ? options : defaultoptions; 1521 PetscCheck(!pre || pre[0] != '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Prefix cannot begin with '-': Instead %s", pre); 1522 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Name must begin with '-': Instead %s", name); 1523 1524 name++; /* skip starting dash */ 1525 1526 /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */ 1527 if (pre && pre[0]) { 1528 char *ptr = buf; 1529 if (name[0] == '-') { 1530 *ptr++ = '-'; 1531 name++; 1532 } 1533 PetscCall(PetscStrncpy(ptr, pre, sizeof(buf) + (size_t)(ptr - buf))); 1534 PetscCall(PetscStrlcat(buf, name, sizeof(buf))); 1535 name = buf; 1536 } 1537 1538 if (PetscDefined(USE_DEBUG)) { 1539 PetscBool valid; 1540 char key[MAXOPTNAME + 1] = "-"; 1541 PetscCall(PetscStrncpy(key + 1, name, sizeof(key) - 1)); 1542 PetscCall(PetscOptionsValidKey(key, &valid)); 1543 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid option '%s' obtained from pre='%s' and name='%s'", key, pre ? pre : "", name); 1544 } 1545 1546 /* determine the location and number of all _%d_ in the key */ 1547 { 1548 int i, j; 1549 for (i = 0; name[i]; i++) { 1550 if (name[i] == '_') { 1551 for (j = i + 1; name[j]; j++) { 1552 if (name[j] >= '0' && name[j] <= '9') continue; 1553 if (name[j] == '_' && j > i + 1) { /* found a number */ 1554 locs[numCnt] = i + 1; 1555 loce[numCnt++] = j + 1; 1556 } 1557 i = j - 1; 1558 break; 1559 } 1560 } 1561 } 1562 } 1563 1564 { /* slow search */ 1565 int c, i; 1566 size_t len; 1567 PetscBool match; 1568 1569 for (c = -1; c < numCnt; ++c) { 1570 char opt[MAXOPTNAME + 1] = "", tmp[MAXOPTNAME]; 1571 1572 if (c < 0) { 1573 PetscCall(PetscStrcpy(opt, name)); 1574 } else { 1575 PetscCall(PetscStrncpy(tmp, name, PetscMin((size_t)(locs[c] + 1), sizeof(tmp)))); 1576 PetscCall(PetscStrlcat(opt, tmp, sizeof(opt))); 1577 PetscCall(PetscStrlcat(opt, name + loce[c], sizeof(opt))); 1578 } 1579 PetscCall(PetscStrlen(opt, &len)); 1580 for (i = 0; i < options->N; i++) { 1581 PetscCall(PetscStrncmp(options->names[i], opt, len, &match)); 1582 if (match) { 1583 options->used[i] = PETSC_TRUE; 1584 if (value) *value = options->values[i]; 1585 if (set) *set = PETSC_TRUE; 1586 PetscFunctionReturn(0); 1587 } 1588 } 1589 } 1590 } 1591 1592 if (set) *set = PETSC_FALSE; 1593 PetscFunctionReturn(0); 1594 } 1595 1596 /*@C 1597 PetscOptionsReject - Generates an error if a certain option is given. 1598 1599 Not Collective 1600 1601 Input Parameters: 1602 + options - options database, use NULL for default global database 1603 . pre - the option prefix (may be NULL) 1604 . name - the option name one is seeking 1605 - mess - error message (may be NULL) 1606 1607 Level: advanced 1608 1609 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, `OptionsHasName()`, 1610 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 1611 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 1612 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 1613 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 1614 `PetscOptionsFList()`, `PetscOptionsEList()` 1615 @*/ 1616 PetscErrorCode PetscOptionsReject(PetscOptions options, const char pre[], const char name[], const char mess[]) { 1617 PetscBool flag = PETSC_FALSE; 1618 1619 PetscFunctionBegin; 1620 PetscCall(PetscOptionsHasName(options, pre, name, &flag)); 1621 if (flag) { 1622 PetscCheck(!mess || !mess[0], PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Program has disabled option: -%s%s with %s", pre ? pre : "", name + 1, mess); 1623 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Program has disabled option: -%s%s", pre ? pre : "", name + 1); 1624 } 1625 PetscFunctionReturn(0); 1626 } 1627 1628 /*@C 1629 PetscOptionsHasHelp - Determines whether the "-help" option is in the database. 1630 1631 Not Collective 1632 1633 Input Parameter: 1634 . options - options database, use NULL for default global database 1635 1636 Output Parameter: 1637 . set - `PETSC_TRUE` if found else `PETSC_FALSE`. 1638 1639 Level: advanced 1640 1641 .seealso: `PetscOptionsHasName()` 1642 @*/ 1643 PetscErrorCode PetscOptionsHasHelp(PetscOptions options, PetscBool *set) { 1644 PetscFunctionBegin; 1645 PetscValidBoolPointer(set, 2); 1646 options = options ? options : defaultoptions; 1647 *set = options->help; 1648 PetscFunctionReturn(0); 1649 } 1650 1651 PetscErrorCode PetscOptionsHasHelpIntro_Internal(PetscOptions options, PetscBool *set) { 1652 PetscFunctionBegin; 1653 PetscValidBoolPointer(set, 2); 1654 options = options ? options : defaultoptions; 1655 *set = options->help_intro; 1656 PetscFunctionReturn(0); 1657 } 1658 1659 /*@C 1660 PetscOptionsHasName - Determines whether a certain option is given in the database. This returns true whether the option is a number, string or Boolean, even 1661 if its value is set to false. 1662 1663 Not Collective 1664 1665 Input Parameters: 1666 + options - options database, use NULL for default global database 1667 . pre - string to prepend to the name or NULL 1668 - name - the option one is seeking 1669 1670 Output Parameter: 1671 . set - `PETSC_TRUE` if found else `PETSC_FALSE`. 1672 1673 Level: beginner 1674 1675 Note: 1676 In many cases you probably want to use `PetscOptionsGetBool()` instead of calling this, to allowing toggling values. 1677 1678 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 1679 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 1680 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 1681 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 1682 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 1683 `PetscOptionsFList()`, `PetscOptionsEList()` 1684 @*/ 1685 PetscErrorCode PetscOptionsHasName(PetscOptions options, const char pre[], const char name[], PetscBool *set) { 1686 const char *value; 1687 PetscBool flag; 1688 1689 PetscFunctionBegin; 1690 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 1691 if (set) *set = flag; 1692 PetscFunctionReturn(0); 1693 } 1694 1695 /*@C 1696 PetscOptionsGetAll - Lists all the options the program was run with in a single string. 1697 1698 Not Collective 1699 1700 Input Parameter: 1701 . options - the options database, use NULL for the default global database 1702 1703 Output Parameter: 1704 . copts - pointer where string pointer is stored 1705 1706 Notes: 1707 The array and each entry in the array should be freed with `PetscFree()` 1708 1709 Each process may have different values depending on how the options were inserted into the database 1710 1711 Level: advanced 1712 1713 .seealso: `PetscOptionsAllUsed()`, `PetscOptionsView()`, `PetscOptionsPush()`, `PetscOptionsPop()` 1714 @*/ 1715 PetscErrorCode PetscOptionsGetAll(PetscOptions options, char *copts[]) { 1716 PetscInt i; 1717 size_t len = 1, lent = 0; 1718 char *coptions = NULL; 1719 1720 PetscFunctionBegin; 1721 PetscValidPointer(copts, 2); 1722 options = options ? options : defaultoptions; 1723 /* count the length of the required string */ 1724 for (i = 0; i < options->N; i++) { 1725 PetscCall(PetscStrlen(options->names[i], &lent)); 1726 len += 2 + lent; 1727 if (options->values[i]) { 1728 PetscCall(PetscStrlen(options->values[i], &lent)); 1729 len += 1 + lent; 1730 } 1731 } 1732 PetscCall(PetscMalloc1(len, &coptions)); 1733 coptions[0] = 0; 1734 for (i = 0; i < options->N; i++) { 1735 PetscCall(PetscStrcat(coptions, "-")); 1736 PetscCall(PetscStrcat(coptions, options->names[i])); 1737 PetscCall(PetscStrcat(coptions, " ")); 1738 if (options->values[i]) { 1739 PetscCall(PetscStrcat(coptions, options->values[i])); 1740 PetscCall(PetscStrcat(coptions, " ")); 1741 } 1742 } 1743 *copts = coptions; 1744 PetscFunctionReturn(0); 1745 } 1746 1747 /*@C 1748 PetscOptionsUsed - Indicates if PETSc has used a particular option set in the database 1749 1750 Not Collective 1751 1752 Input Parameters: 1753 + options - options database, use NULL for default global database 1754 - name - string name of option 1755 1756 Output Parameter: 1757 . used - `PETSC_TRUE` if the option was used, otherwise false, including if option was not found in options database 1758 1759 Level: advanced 1760 1761 Note: 1762 The value returned may be different on each process and depends on which options have been processed 1763 on the given process 1764 1765 .seealso: `PetscOptionsView()`, `PetscOptionsLeft()`, `PetscOptionsAllUsed()` 1766 @*/ 1767 PetscErrorCode PetscOptionsUsed(PetscOptions options, const char *name, PetscBool *used) { 1768 PetscInt i; 1769 1770 PetscFunctionBegin; 1771 PetscValidCharPointer(name, 2); 1772 PetscValidBoolPointer(used, 3); 1773 options = options ? options : defaultoptions; 1774 *used = PETSC_FALSE; 1775 for (i = 0; i < options->N; i++) { 1776 PetscCall(PetscStrcasecmp(options->names[i], name, used)); 1777 if (*used) { 1778 *used = options->used[i]; 1779 break; 1780 } 1781 } 1782 PetscFunctionReturn(0); 1783 } 1784 1785 /*@ 1786 PetscOptionsAllUsed - Returns a count of the number of options in the 1787 database that have never been selected. 1788 1789 Not Collective 1790 1791 Input Parameter: 1792 . options - options database, use NULL for default global database 1793 1794 Output Parameter: 1795 . N - count of options not used 1796 1797 Level: advanced 1798 1799 Note: 1800 The value returned may be different on each process and depends on which options have been processed 1801 on the given process 1802 1803 .seealso: `PetscOptionsView()` 1804 @*/ 1805 PetscErrorCode PetscOptionsAllUsed(PetscOptions options, PetscInt *N) { 1806 PetscInt i, n = 0; 1807 1808 PetscFunctionBegin; 1809 PetscValidIntPointer(N, 2); 1810 options = options ? options : defaultoptions; 1811 for (i = 0; i < options->N; i++) { 1812 if (!options->used[i]) n++; 1813 } 1814 *N = n; 1815 PetscFunctionReturn(0); 1816 } 1817 1818 /*@ 1819 PetscOptionsLeft - Prints to screen any options that were set and never used. 1820 1821 Not Collective 1822 1823 Input Parameter: 1824 . options - options database; use NULL for default global database 1825 1826 Options Database Key: 1827 . -options_left - activates `PetscOptionsAllUsed()` within `PetscFinalize()` 1828 1829 Notes: 1830 This is rarely used directly, it is called by `PetscFinalize()` in debug more or if -options_left 1831 is passed otherwise to help users determine possible mistakes in their usage of options. This 1832 only prints values on process zero of `PETSC_COMM_WORLD`. 1833 1834 Other processes depending the objects 1835 used may have different options that are left unused. 1836 1837 Level: advanced 1838 1839 .seealso: `PetscOptionsAllUsed()` 1840 @*/ 1841 PetscErrorCode PetscOptionsLeft(PetscOptions options) { 1842 PetscInt i; 1843 PetscInt cnt = 0; 1844 PetscOptions toptions; 1845 1846 PetscFunctionBegin; 1847 toptions = options ? options : defaultoptions; 1848 for (i = 0; i < toptions->N; i++) { 1849 if (!toptions->used[i]) { 1850 if (PetscCIOption(toptions->names[i])) continue; 1851 if (toptions->values[i]) { 1852 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Option left: name:-%s value: %s\n", toptions->names[i], toptions->values[i])); 1853 } else { 1854 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Option left: name:-%s (no value)\n", toptions->names[i])); 1855 } 1856 } 1857 } 1858 if (!options) { 1859 toptions = defaultoptions; 1860 while (toptions->previous) { 1861 cnt++; 1862 toptions = toptions->previous; 1863 } 1864 if (cnt) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Option left: You may have forgotten some calls to PetscOptionsPop(),\n PetscOptionsPop() has been called %" PetscInt_FMT " less times than PetscOptionsPush()\n", cnt)); 1865 } 1866 PetscFunctionReturn(0); 1867 } 1868 1869 /*@C 1870 PetscOptionsLeftGet - Returns all options that were set and never used. 1871 1872 Not Collective 1873 1874 Input Parameter: 1875 . options - options database, use NULL for default global database 1876 1877 Output Parameters: 1878 + N - count of options not used 1879 . names - names of options not used 1880 - values - values of options not used 1881 1882 Level: advanced 1883 1884 Notes: 1885 Users should call `PetscOptionsLeftRestore()` to free the memory allocated in this routine 1886 1887 The value returned may be different on each process and depends on which options have been processed 1888 on the given process 1889 1890 .seealso: `PetscOptionsAllUsed()`, `PetscOptionsLeft()` 1891 @*/ 1892 PetscErrorCode PetscOptionsLeftGet(PetscOptions options, PetscInt *N, char **names[], char **values[]) { 1893 PetscInt i, n; 1894 1895 PetscFunctionBegin; 1896 if (N) PetscValidIntPointer(N, 2); 1897 if (names) PetscValidPointer(names, 3); 1898 if (values) PetscValidPointer(values, 4); 1899 options = options ? options : defaultoptions; 1900 1901 /* The number of unused PETSc options */ 1902 n = 0; 1903 for (i = 0; i < options->N; i++) { 1904 if (PetscCIOption(options->names[i])) continue; 1905 if (!options->used[i]) n++; 1906 } 1907 if (N) *N = n; 1908 if (names) PetscCall(PetscMalloc1(n, names)); 1909 if (values) PetscCall(PetscMalloc1(n, values)); 1910 1911 n = 0; 1912 if (names || values) { 1913 for (i = 0; i < options->N; i++) { 1914 if (!options->used[i]) { 1915 if (PetscCIOption(options->names[i])) continue; 1916 if (names) (*names)[n] = options->names[i]; 1917 if (values) (*values)[n] = options->values[i]; 1918 n++; 1919 } 1920 } 1921 } 1922 PetscFunctionReturn(0); 1923 } 1924 1925 /*@C 1926 PetscOptionsLeftRestore - Free memory for the unused PETSc options obtained using `PetscOptionsLeftGet()`. 1927 1928 Not Collective 1929 1930 Input Parameters: 1931 + options - options database, use NULL for default global database 1932 . names - names of options not used 1933 - values - values of options not used 1934 1935 Level: advanced 1936 1937 .seealso: `PetscOptionsAllUsed()`, `PetscOptionsLeft()`, `PetscOptionsLeftGet()` 1938 @*/ 1939 PetscErrorCode PetscOptionsLeftRestore(PetscOptions options, PetscInt *N, char **names[], char **values[]) { 1940 PetscFunctionBegin; 1941 if (N) PetscValidIntPointer(N, 2); 1942 if (names) PetscValidPointer(names, 3); 1943 if (values) PetscValidPointer(values, 4); 1944 if (N) *N = 0; 1945 if (names) PetscCall(PetscFree(*names)); 1946 if (values) PetscCall(PetscFree(*values)); 1947 PetscFunctionReturn(0); 1948 } 1949 1950 /*@C 1951 PetscOptionsMonitorDefault - Print all options set value events using the supplied `PetscViewer`. 1952 1953 Logically Collective on ctx 1954 1955 Input Parameters: 1956 + name - option name string 1957 . value - option value string 1958 - ctx - an ASCII viewer or NULL 1959 1960 Level: intermediate 1961 1962 Notes: 1963 If ctx is NULL, `PetscPrintf()` is used. 1964 The first MPI rank in the `PetscViewer` viewer actually prints the values, other 1965 processes may have different values set 1966 1967 If `PetscCIEnabled` then do not print the test harness options 1968 1969 .seealso: `PetscOptionsMonitorSet()` 1970 @*/ 1971 PetscErrorCode PetscOptionsMonitorDefault(const char name[], const char value[], void *ctx) { 1972 PetscFunctionBegin; 1973 if (PetscCIOption(name)) PetscFunctionReturn(0); 1974 1975 if (ctx) { 1976 PetscViewer viewer = (PetscViewer)ctx; 1977 if (!value) { 1978 PetscCall(PetscViewerASCIIPrintf(viewer, "Removing option: %s\n", name)); 1979 } else if (!value[0]) { 1980 PetscCall(PetscViewerASCIIPrintf(viewer, "Setting option: %s (no value)\n", name)); 1981 } else { 1982 PetscCall(PetscViewerASCIIPrintf(viewer, "Setting option: %s = %s\n", name, value)); 1983 } 1984 } else { 1985 MPI_Comm comm = PETSC_COMM_WORLD; 1986 if (!value) { 1987 PetscCall(PetscPrintf(comm, "Removing option: %s\n", name)); 1988 } else if (!value[0]) { 1989 PetscCall(PetscPrintf(comm, "Setting option: %s (no value)\n", name)); 1990 } else { 1991 PetscCall(PetscPrintf(comm, "Setting option: %s = %s\n", name, value)); 1992 } 1993 } 1994 PetscFunctionReturn(0); 1995 } 1996 1997 /*@C 1998 PetscOptionsMonitorSet - Sets an ADDITIONAL function to be called at every method that 1999 modified the PETSc options database. 2000 2001 Not Collective 2002 2003 Input Parameters: 2004 + monitor - pointer to function (if this is NULL, it turns off monitoring 2005 . mctx - [optional] context for private data for the 2006 monitor routine (use NULL if no context is desired) 2007 - monitordestroy - [optional] routine that frees monitor context 2008 (may be NULL) 2009 2010 Calling Sequence of monitor: 2011 $ monitor (const char name[], const char value[], void *mctx) 2012 2013 + name - option name string 2014 . value - option value string 2015 - mctx - optional monitoring context, as set by `PetscOptionsMonitorSet()` 2016 2017 Options Database Keys: 2018 See `PetscInitialize()` for options related to option database monitoring. 2019 2020 Notes: 2021 The default is to do nothing. To print the name and value of options 2022 being inserted into the database, use `PetscOptionsMonitorDefault()` as the monitoring routine, 2023 with a null monitoring context. 2024 2025 Several different monitoring routines may be set by calling 2026 `PetscOptionsMonitorSet()` multiple times; all will be called in the 2027 order in which they were set. 2028 2029 Level: intermediate 2030 2031 .seealso: `PetscOptionsMonitorDefault()`, `PetscInitialize()` 2032 @*/ 2033 PetscErrorCode PetscOptionsMonitorSet(PetscErrorCode (*monitor)(const char name[], const char value[], void *), void *mctx, PetscErrorCode (*monitordestroy)(void **)) { 2034 PetscOptions options = defaultoptions; 2035 2036 PetscFunctionBegin; 2037 if (options->monitorCancel) PetscFunctionReturn(0); 2038 PetscCheck(options->numbermonitors < MAXOPTIONSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many PetscOptions monitors set"); 2039 options->monitor[options->numbermonitors] = monitor; 2040 options->monitordestroy[options->numbermonitors] = monitordestroy; 2041 options->monitorcontext[options->numbermonitors++] = (void *)mctx; 2042 PetscFunctionReturn(0); 2043 } 2044 2045 /* 2046 PetscOptionsStringToBool - Converts string to PetscBool, handles cases like "yes", "no", "true", "false", "0", "1", "off", "on". 2047 Empty string is considered as true. 2048 */ 2049 PetscErrorCode PetscOptionsStringToBool(const char value[], PetscBool *a) { 2050 PetscBool istrue, isfalse; 2051 size_t len; 2052 2053 PetscFunctionBegin; 2054 /* PetscStrlen() returns 0 for NULL or "" */ 2055 PetscCall(PetscStrlen(value, &len)); 2056 if (!len) { 2057 *a = PETSC_TRUE; 2058 PetscFunctionReturn(0); 2059 } 2060 PetscCall(PetscStrcasecmp(value, "TRUE", &istrue)); 2061 if (istrue) { 2062 *a = PETSC_TRUE; 2063 PetscFunctionReturn(0); 2064 } 2065 PetscCall(PetscStrcasecmp(value, "YES", &istrue)); 2066 if (istrue) { 2067 *a = PETSC_TRUE; 2068 PetscFunctionReturn(0); 2069 } 2070 PetscCall(PetscStrcasecmp(value, "1", &istrue)); 2071 if (istrue) { 2072 *a = PETSC_TRUE; 2073 PetscFunctionReturn(0); 2074 } 2075 PetscCall(PetscStrcasecmp(value, "on", &istrue)); 2076 if (istrue) { 2077 *a = PETSC_TRUE; 2078 PetscFunctionReturn(0); 2079 } 2080 PetscCall(PetscStrcasecmp(value, "FALSE", &isfalse)); 2081 if (isfalse) { 2082 *a = PETSC_FALSE; 2083 PetscFunctionReturn(0); 2084 } 2085 PetscCall(PetscStrcasecmp(value, "NO", &isfalse)); 2086 if (isfalse) { 2087 *a = PETSC_FALSE; 2088 PetscFunctionReturn(0); 2089 } 2090 PetscCall(PetscStrcasecmp(value, "0", &isfalse)); 2091 if (isfalse) { 2092 *a = PETSC_FALSE; 2093 PetscFunctionReturn(0); 2094 } 2095 PetscCall(PetscStrcasecmp(value, "off", &isfalse)); 2096 if (isfalse) { 2097 *a = PETSC_FALSE; 2098 PetscFunctionReturn(0); 2099 } 2100 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown logical value: %s", value); 2101 } 2102 2103 /* 2104 PetscOptionsStringToInt - Converts a string to an integer value. Handles special cases such as "default" and "decide" 2105 */ 2106 PetscErrorCode PetscOptionsStringToInt(const char name[], PetscInt *a) { 2107 size_t len; 2108 PetscBool decide, tdefault, mouse; 2109 2110 PetscFunctionBegin; 2111 PetscCall(PetscStrlen(name, &len)); 2112 PetscCheck(len, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "character string of length zero has no numerical value"); 2113 2114 PetscCall(PetscStrcasecmp(name, "PETSC_DEFAULT", &tdefault)); 2115 if (!tdefault) PetscCall(PetscStrcasecmp(name, "DEFAULT", &tdefault)); 2116 PetscCall(PetscStrcasecmp(name, "PETSC_DECIDE", &decide)); 2117 if (!decide) PetscCall(PetscStrcasecmp(name, "DECIDE", &decide)); 2118 PetscCall(PetscStrcasecmp(name, "mouse", &mouse)); 2119 2120 if (tdefault) *a = PETSC_DEFAULT; 2121 else if (decide) *a = PETSC_DECIDE; 2122 else if (mouse) *a = -1; 2123 else { 2124 char *endptr; 2125 long strtolval; 2126 2127 strtolval = strtol(name, &endptr, 10); 2128 PetscCheck((size_t)(endptr - name) == len, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s has no integer value (do not include . in it)", name); 2129 2130 #if defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE_ATOLL) 2131 (void)strtolval; 2132 *a = atoll(name); 2133 #elif defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE___INT64) 2134 (void)strtolval; 2135 *a = _atoi64(name); 2136 #else 2137 *a = (PetscInt)strtolval; 2138 #endif 2139 } 2140 PetscFunctionReturn(0); 2141 } 2142 2143 #if defined(PETSC_USE_REAL___FLOAT128) 2144 #include <quadmath.h> 2145 #endif 2146 2147 static PetscErrorCode PetscStrtod(const char name[], PetscReal *a, char **endptr) { 2148 PetscFunctionBegin; 2149 #if defined(PETSC_USE_REAL___FLOAT128) 2150 *a = strtoflt128(name, endptr); 2151 #else 2152 *a = (PetscReal)strtod(name, endptr); 2153 #endif 2154 PetscFunctionReturn(0); 2155 } 2156 2157 static PetscErrorCode PetscStrtoz(const char name[], PetscScalar *a, char **endptr, PetscBool *isImaginary) { 2158 PetscBool hasi = PETSC_FALSE; 2159 char *ptr; 2160 PetscReal strtoval; 2161 2162 PetscFunctionBegin; 2163 PetscCall(PetscStrtod(name, &strtoval, &ptr)); 2164 if (ptr == name) { 2165 strtoval = 1.; 2166 hasi = PETSC_TRUE; 2167 if (name[0] == 'i') { 2168 ptr++; 2169 } else if (name[0] == '+' && name[1] == 'i') { 2170 ptr += 2; 2171 } else if (name[0] == '-' && name[1] == 'i') { 2172 strtoval = -1.; 2173 ptr += 2; 2174 } 2175 } else if (*ptr == 'i') { 2176 hasi = PETSC_TRUE; 2177 ptr++; 2178 } 2179 *endptr = ptr; 2180 *isImaginary = hasi; 2181 if (hasi) { 2182 #if !defined(PETSC_USE_COMPLEX) 2183 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s contains imaginary but complex not supported ", name); 2184 #else 2185 *a = PetscCMPLX(0., strtoval); 2186 #endif 2187 } else { 2188 *a = strtoval; 2189 } 2190 PetscFunctionReturn(0); 2191 } 2192 2193 /* 2194 Converts a string to PetscReal value. Handles special cases like "default" and "decide" 2195 */ 2196 PetscErrorCode PetscOptionsStringToReal(const char name[], PetscReal *a) { 2197 size_t len; 2198 PetscBool match; 2199 char *endptr; 2200 2201 PetscFunctionBegin; 2202 PetscCall(PetscStrlen(name, &len)); 2203 PetscCheck(len, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "String of length zero has no numerical value"); 2204 2205 PetscCall(PetscStrcasecmp(name, "PETSC_DEFAULT", &match)); 2206 if (!match) PetscCall(PetscStrcasecmp(name, "DEFAULT", &match)); 2207 if (match) { 2208 *a = PETSC_DEFAULT; 2209 PetscFunctionReturn(0); 2210 } 2211 2212 PetscCall(PetscStrcasecmp(name, "PETSC_DECIDE", &match)); 2213 if (!match) PetscCall(PetscStrcasecmp(name, "DECIDE", &match)); 2214 if (match) { 2215 *a = PETSC_DECIDE; 2216 PetscFunctionReturn(0); 2217 } 2218 2219 PetscCall(PetscStrtod(name, a, &endptr)); 2220 PetscCheck((size_t)(endptr - name) == len, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s has no numeric value", name); 2221 PetscFunctionReturn(0); 2222 } 2223 2224 PetscErrorCode PetscOptionsStringToScalar(const char name[], PetscScalar *a) { 2225 PetscBool imag1; 2226 size_t len; 2227 PetscScalar val = 0.; 2228 char *ptr = NULL; 2229 2230 PetscFunctionBegin; 2231 PetscCall(PetscStrlen(name, &len)); 2232 PetscCheck(len, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "character string of length zero has no numerical value"); 2233 PetscCall(PetscStrtoz(name, &val, &ptr, &imag1)); 2234 #if defined(PETSC_USE_COMPLEX) 2235 if ((size_t)(ptr - name) < len) { 2236 PetscBool imag2; 2237 PetscScalar val2; 2238 2239 PetscCall(PetscStrtoz(ptr, &val2, &ptr, &imag2)); 2240 if (imag1) PetscCheck(imag2, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s: must specify imaginary component second", name); 2241 val = PetscCMPLX(PetscRealPart(val), PetscImaginaryPart(val2)); 2242 } 2243 #endif 2244 PetscCheck((size_t)(ptr - name) == len, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s has no numeric value ", name); 2245 *a = val; 2246 PetscFunctionReturn(0); 2247 } 2248 2249 /*@C 2250 PetscOptionsGetBool - Gets the Logical (true or false) value for a particular 2251 option in the database. 2252 2253 Not Collective 2254 2255 Input Parameters: 2256 + options - options database, use NULL for default global database 2257 . pre - the string to prepend to the name or NULL 2258 - name - the option one is seeking 2259 2260 Output Parameters: 2261 + ivalue - the logical value to return 2262 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2263 2264 Level: beginner 2265 2266 Notes: 2267 TRUE, true, YES, yes, nostring, and 1 all translate to `PETSC_TRUE` 2268 FALSE, false, NO, no, and 0 all translate to `PETSC_FALSE` 2269 2270 If the option is given, but no value is provided, then ivalue and set are both given the value `PETSC_TRUE`. That is -requested_bool 2271 is equivalent to -requested_bool true 2272 2273 If the user does not supply the option at all ivalue is NOT changed. Thus 2274 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2275 2276 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, 2277 `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsGetInt()`, `PetscOptionsBool()`, 2278 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2279 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2280 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2281 `PetscOptionsFList()`, `PetscOptionsEList()` 2282 @*/ 2283 PetscErrorCode PetscOptionsGetBool(PetscOptions options, const char pre[], const char name[], PetscBool *ivalue, PetscBool *set) { 2284 const char *value; 2285 PetscBool flag; 2286 2287 PetscFunctionBegin; 2288 PetscValidCharPointer(name, 3); 2289 if (ivalue) PetscValidBoolPointer(ivalue, 4); 2290 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2291 if (flag) { 2292 if (set) *set = PETSC_TRUE; 2293 PetscCall(PetscOptionsStringToBool(value, &flag)); 2294 if (ivalue) *ivalue = flag; 2295 } else { 2296 if (set) *set = PETSC_FALSE; 2297 } 2298 PetscFunctionReturn(0); 2299 } 2300 2301 /*@C 2302 PetscOptionsGetEList - Puts a list of option values that a single one may be selected from 2303 2304 Not Collective 2305 2306 Input Parameters: 2307 + options - options database, use NULL for default global database 2308 . pre - the string to prepend to the name or NULL 2309 . opt - option name 2310 . list - the possible choices (one of these must be selected, anything else is invalid) 2311 - ntext - number of choices 2312 2313 Output Parameters: 2314 + value - the index of the value to return (defaults to zero if the option name is given but no choice is listed) 2315 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2316 2317 Level: intermediate 2318 2319 Notes: 2320 If the user does not supply the option value is NOT changed. Thus 2321 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2322 2323 See `PetscOptionsFList()` for when the choices are given in a `PetscFunctionList` 2324 2325 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 2326 `PetscOptionsHasName()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2327 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2328 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2329 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2330 `PetscOptionsFList()`, `PetscOptionsEList()` 2331 @*/ 2332 PetscErrorCode PetscOptionsGetEList(PetscOptions options, const char pre[], const char opt[], const char *const *list, PetscInt ntext, PetscInt *value, PetscBool *set) { 2333 size_t alen, len = 0, tlen = 0; 2334 char *svalue; 2335 PetscBool aset, flg = PETSC_FALSE; 2336 PetscInt i; 2337 2338 PetscFunctionBegin; 2339 PetscValidCharPointer(opt, 3); 2340 for (i = 0; i < ntext; i++) { 2341 PetscCall(PetscStrlen(list[i], &alen)); 2342 if (alen > len) len = alen; 2343 tlen += len + 1; 2344 } 2345 len += 5; /* a little extra space for user mistypes */ 2346 PetscCall(PetscMalloc1(len, &svalue)); 2347 PetscCall(PetscOptionsGetString(options, pre, opt, svalue, len, &aset)); 2348 if (aset) { 2349 PetscCall(PetscEListFind(ntext, list, svalue, value, &flg)); 2350 if (!flg) { 2351 char *avail, *pavl; 2352 2353 PetscCall(PetscMalloc1(tlen, &avail)); 2354 pavl = avail; 2355 for (i = 0; i < ntext; i++) { 2356 PetscCall(PetscStrlen(list[i], &alen)); 2357 PetscCall(PetscStrcpy(pavl, list[i])); 2358 pavl += alen; 2359 PetscCall(PetscStrcpy(pavl, " ")); 2360 pavl += 1; 2361 } 2362 PetscCall(PetscStrtolower(avail)); 2363 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_USER, "Unknown option %s for -%s%s. Available options: %s", svalue, pre ? pre : "", opt + 1, avail); 2364 } 2365 if (set) *set = PETSC_TRUE; 2366 } else if (set) *set = PETSC_FALSE; 2367 PetscCall(PetscFree(svalue)); 2368 PetscFunctionReturn(0); 2369 } 2370 2371 /*@C 2372 PetscOptionsGetEnum - Gets the enum value for a particular option in the database. 2373 2374 Not Collective 2375 2376 Input Parameters: 2377 + options - options database, use NULL for default global database 2378 . pre - option prefix or NULL 2379 . opt - option name 2380 - list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null 2381 2382 Output Parameters: 2383 + value - the value to return 2384 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2385 2386 Level: beginner 2387 2388 Notes: 2389 If the user does not supply the option value is NOT changed. Thus 2390 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2391 2392 List is usually something like `PCASMTypes` or some other predefined list of enum names 2393 2394 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, `PetscOptionsGetInt()`, 2395 `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()` 2396 `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`, 2397 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2398 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2399 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2400 `PetscOptionsFList()`, `PetscOptionsEList()`, `PetscOptionsGetEList()`, `PetscOptionsEnum()` 2401 @*/ 2402 PetscErrorCode PetscOptionsGetEnum(PetscOptions options, const char pre[], const char opt[], const char *const *list, PetscEnum *value, PetscBool *set) { 2403 PetscInt ntext = 0, tval; 2404 PetscBool fset; 2405 2406 PetscFunctionBegin; 2407 PetscValidCharPointer(opt, 3); 2408 while (list[ntext++]) PetscCheck(ntext <= 50, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "List argument appears to be wrong or have more than 50 entries"); 2409 PetscCheck(ntext >= 3, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "List argument must have at least two entries: typename and type prefix"); 2410 ntext -= 3; 2411 PetscCall(PetscOptionsGetEList(options, pre, opt, list, ntext, &tval, &fset)); 2412 /* with PETSC_USE_64BIT_INDICES sizeof(PetscInt) != sizeof(PetscEnum) */ 2413 if (fset) *value = (PetscEnum)tval; 2414 if (set) *set = fset; 2415 PetscFunctionReturn(0); 2416 } 2417 2418 /*@C 2419 PetscOptionsGetInt - Gets the integer value for a particular option in the database. 2420 2421 Not Collective 2422 2423 Input Parameters: 2424 + options - options database, use NULL for default global database 2425 . pre - the string to prepend to the name or NULL 2426 - name - the option one is seeking 2427 2428 Output Parameters: 2429 + ivalue - the integer value to return 2430 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2431 2432 Level: beginner 2433 2434 Notes: 2435 If the user does not supply the option ivalue is NOT changed. Thus 2436 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2437 2438 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, 2439 `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()` 2440 `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`, 2441 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2442 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2443 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2444 `PetscOptionsFList()`, `PetscOptionsEList()` 2445 @*/ 2446 PetscErrorCode PetscOptionsGetInt(PetscOptions options, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set) { 2447 const char *value; 2448 PetscBool flag; 2449 2450 PetscFunctionBegin; 2451 PetscValidCharPointer(name, 3); 2452 PetscValidIntPointer(ivalue, 4); 2453 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2454 if (flag) { 2455 if (!value) { 2456 if (set) *set = PETSC_FALSE; 2457 } else { 2458 if (set) *set = PETSC_TRUE; 2459 PetscCall(PetscOptionsStringToInt(value, ivalue)); 2460 } 2461 } else { 2462 if (set) *set = PETSC_FALSE; 2463 } 2464 PetscFunctionReturn(0); 2465 } 2466 2467 /*@C 2468 PetscOptionsGetReal - Gets the double precision value for a particular 2469 option in the database. 2470 2471 Not Collective 2472 2473 Input Parameters: 2474 + options - options database, use NULL for default global database 2475 . pre - string to prepend to each name or NULL 2476 - name - the option one is seeking 2477 2478 Output Parameters: 2479 + dvalue - the double value to return 2480 - set - `PETSC_TRUE` if found, `PETSC_FALSE` if not found 2481 2482 Note: 2483 If the user does not supply the option dvalue is NOT changed. Thus 2484 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2485 2486 Level: beginner 2487 2488 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2489 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2490 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2491 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2492 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2493 `PetscOptionsFList()`, `PetscOptionsEList()` 2494 @*/ 2495 PetscErrorCode PetscOptionsGetReal(PetscOptions options, const char pre[], const char name[], PetscReal *dvalue, PetscBool *set) { 2496 const char *value; 2497 PetscBool flag; 2498 2499 PetscFunctionBegin; 2500 PetscValidCharPointer(name, 3); 2501 PetscValidRealPointer(dvalue, 4); 2502 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2503 if (flag) { 2504 if (!value) { 2505 if (set) *set = PETSC_FALSE; 2506 } else { 2507 if (set) *set = PETSC_TRUE; 2508 PetscCall(PetscOptionsStringToReal(value, dvalue)); 2509 } 2510 } else { 2511 if (set) *set = PETSC_FALSE; 2512 } 2513 PetscFunctionReturn(0); 2514 } 2515 2516 /*@C 2517 PetscOptionsGetScalar - Gets the scalar value for a particular 2518 option in the database. 2519 2520 Not Collective 2521 2522 Input Parameters: 2523 + options - options database, use NULL for default global database 2524 . pre - string to prepend to each name or NULL 2525 - name - the option one is seeking 2526 2527 Output Parameters: 2528 + dvalue - the double value to return 2529 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2530 2531 Level: beginner 2532 2533 Usage: 2534 A complex number 2+3i must be specified with NO spaces 2535 2536 Note: 2537 If the user does not supply the option dvalue is NOT changed. Thus 2538 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2539 2540 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2541 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2542 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2543 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2544 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2545 `PetscOptionsFList()`, `PetscOptionsEList()` 2546 @*/ 2547 PetscErrorCode PetscOptionsGetScalar(PetscOptions options, const char pre[], const char name[], PetscScalar *dvalue, PetscBool *set) { 2548 const char *value; 2549 PetscBool flag; 2550 2551 PetscFunctionBegin; 2552 PetscValidCharPointer(name, 3); 2553 PetscValidScalarPointer(dvalue, 4); 2554 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2555 if (flag) { 2556 if (!value) { 2557 if (set) *set = PETSC_FALSE; 2558 } else { 2559 #if !defined(PETSC_USE_COMPLEX) 2560 PetscCall(PetscOptionsStringToReal(value, dvalue)); 2561 #else 2562 PetscCall(PetscOptionsStringToScalar(value, dvalue)); 2563 #endif 2564 if (set) *set = PETSC_TRUE; 2565 } 2566 } else { /* flag */ 2567 if (set) *set = PETSC_FALSE; 2568 } 2569 PetscFunctionReturn(0); 2570 } 2571 2572 /*@C 2573 PetscOptionsGetString - Gets the string value for a particular option in 2574 the database. 2575 2576 Not Collective 2577 2578 Input Parameters: 2579 + options - options database, use NULL for default global database 2580 . pre - string to prepend to name or NULL 2581 . name - the option one is seeking 2582 - len - maximum length of the string including null termination 2583 2584 Output Parameters: 2585 + string - location to copy string 2586 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2587 2588 Level: beginner 2589 2590 Fortran Note: 2591 The Fortran interface is slightly different from the C/C++ 2592 interface (len is not used). Sample usage in Fortran follows 2593 .vb 2594 character *20 string 2595 PetscErrorCode ierr 2596 PetscBool set 2597 call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s',string,set,ierr) 2598 .ve 2599 2600 Note: 2601 if the option is given but no string is provided then an empty string is returned and set is given the value of `PETSC_TRUE` 2602 2603 If the user does not use the option then the string is not changed. Thus 2604 you should ALWAYS initialize the string if you access it without first checking if the set flag is true. 2605 2606 Even if the user provided no string (for example -optionname -someotheroption) the flag is set to PETSC_TRUE (and the string is fulled with nulls). 2607 2608 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 2609 `PetscOptionsHasName()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2610 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2611 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2612 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2613 `PetscOptionsFList()`, `PetscOptionsEList()` 2614 @*/ 2615 PetscErrorCode PetscOptionsGetString(PetscOptions options, const char pre[], const char name[], char string[], size_t len, PetscBool *set) { 2616 const char *value; 2617 PetscBool flag; 2618 2619 PetscFunctionBegin; 2620 PetscValidCharPointer(name, 3); 2621 PetscValidCharPointer(string, 4); 2622 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2623 if (!flag) { 2624 if (set) *set = PETSC_FALSE; 2625 } else { 2626 if (set) *set = PETSC_TRUE; 2627 if (value) PetscCall(PetscStrncpy(string, value, len)); 2628 else PetscCall(PetscArrayzero(string, len)); 2629 } 2630 PetscFunctionReturn(0); 2631 } 2632 2633 char *PetscOptionsGetStringMatlab(PetscOptions options, const char pre[], const char name[]) { 2634 const char *value; 2635 PetscBool flag; 2636 2637 PetscFunctionBegin; 2638 if (PetscOptionsFindPair(options, pre, name, &value, &flag)) PetscFunctionReturn(NULL); 2639 if (flag) PetscFunctionReturn((char *)value); 2640 PetscFunctionReturn(NULL); 2641 } 2642 2643 /*@C 2644 PetscOptionsGetBoolArray - Gets an array of Logical (true or false) values for a particular 2645 option in the database. The values must be separated with commas with no intervening spaces. 2646 2647 Not Collective 2648 2649 Input Parameters: 2650 + options - options database, use NULL for default global database 2651 . pre - string to prepend to each name or NULL 2652 - name - the option one is seeking 2653 2654 Output Parameters: 2655 + dvalue - the integer values to return 2656 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2657 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2658 2659 Level: beginner 2660 2661 Note: 2662 TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE. FALSE, false, NO, no, and 0 all translate to PETSC_FALSE 2663 2664 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2665 `PetscOptionsGetString()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2666 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2667 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2668 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2669 `PetscOptionsFList()`, `PetscOptionsEList()` 2670 @*/ 2671 PetscErrorCode PetscOptionsGetBoolArray(PetscOptions options, const char pre[], const char name[], PetscBool dvalue[], PetscInt *nmax, PetscBool *set) { 2672 const char *svalue; 2673 char *value; 2674 PetscInt n = 0; 2675 PetscBool flag; 2676 PetscToken token; 2677 2678 PetscFunctionBegin; 2679 PetscValidCharPointer(name, 3); 2680 PetscValidBoolPointer(dvalue, 4); 2681 PetscValidIntPointer(nmax, 5); 2682 2683 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2684 if (!flag || !svalue) { 2685 if (set) *set = PETSC_FALSE; 2686 *nmax = 0; 2687 PetscFunctionReturn(0); 2688 } 2689 if (set) *set = PETSC_TRUE; 2690 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2691 PetscCall(PetscTokenFind(token, &value)); 2692 while (value && n < *nmax) { 2693 PetscCall(PetscOptionsStringToBool(value, dvalue)); 2694 PetscCall(PetscTokenFind(token, &value)); 2695 dvalue++; 2696 n++; 2697 } 2698 PetscCall(PetscTokenDestroy(&token)); 2699 *nmax = n; 2700 PetscFunctionReturn(0); 2701 } 2702 2703 /*@C 2704 PetscOptionsGetEnumArray - Gets an array of enum values for a particular option in the database. 2705 2706 Not Collective 2707 2708 Input Parameters: 2709 + options - options database, use NULL for default global database 2710 . pre - option prefix or NULL 2711 . name - option name 2712 - list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null 2713 2714 Output Parameters: 2715 + ivalue - the enum values to return 2716 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2717 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2718 2719 Level: beginner 2720 2721 Notes: 2722 The array must be passed as a comma separated list. 2723 2724 There must be no intervening spaces between the values. 2725 2726 list is usually something like `PCASMTypes` or some other predefined list of enum names. 2727 2728 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, `PetscOptionsGetInt()`, 2729 `PetscOptionsGetEnum()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()` 2730 `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`, `PetscOptionsName()`, 2731 `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, 2732 `PetscOptionsScalar()`, `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2733 `PetscOptionsFList()`, `PetscOptionsEList()`, `PetscOptionsGetEList()`, `PetscOptionsEnum()` 2734 @*/ 2735 PetscErrorCode PetscOptionsGetEnumArray(PetscOptions options, const char pre[], const char name[], const char *const *list, PetscEnum ivalue[], PetscInt *nmax, PetscBool *set) { 2736 const char *svalue; 2737 char *value; 2738 PetscInt n = 0; 2739 PetscEnum evalue; 2740 PetscBool flag; 2741 PetscToken token; 2742 2743 PetscFunctionBegin; 2744 PetscValidCharPointer(name, 3); 2745 PetscValidPointer(list, 4); 2746 PetscValidPointer(ivalue, 5); 2747 PetscValidIntPointer(nmax, 6); 2748 2749 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2750 if (!flag || !svalue) { 2751 if (set) *set = PETSC_FALSE; 2752 *nmax = 0; 2753 PetscFunctionReturn(0); 2754 } 2755 if (set) *set = PETSC_TRUE; 2756 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2757 PetscCall(PetscTokenFind(token, &value)); 2758 while (value && n < *nmax) { 2759 PetscCall(PetscEnumFind(list, value, &evalue, &flag)); 2760 PetscCheck(flag, PETSC_COMM_SELF, PETSC_ERR_USER, "Unknown enum value '%s' for -%s%s", svalue, pre ? pre : "", name + 1); 2761 ivalue[n++] = evalue; 2762 PetscCall(PetscTokenFind(token, &value)); 2763 } 2764 PetscCall(PetscTokenDestroy(&token)); 2765 *nmax = n; 2766 PetscFunctionReturn(0); 2767 } 2768 2769 /*@C 2770 PetscOptionsGetIntArray - Gets an array of integer values for a particular option in the database. 2771 2772 Not Collective 2773 2774 Input Parameters: 2775 + options - options database, use NULL for default global database 2776 . pre - string to prepend to each name or NULL 2777 - name - the option one is seeking 2778 2779 Output Parameters: 2780 + ivalue - the integer values to return 2781 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2782 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2783 2784 Level: beginner 2785 2786 Notes: 2787 The array can be passed as 2788 + a comma separated list - 0,1,2,3,4,5,6,7 2789 . a range (start\-end+1) - 0-8 2790 . a range with given increment (start\-end+1:inc) - 0-7:2 2791 - a combination of values and ranges separated by commas - 0,1-8,8-15:2 2792 2793 There must be no intervening spaces between the values. 2794 2795 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2796 `PetscOptionsGetString()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2797 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2798 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2799 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2800 `PetscOptionsFList()`, `PetscOptionsEList()` 2801 @*/ 2802 PetscErrorCode PetscOptionsGetIntArray(PetscOptions options, const char pre[], const char name[], PetscInt ivalue[], PetscInt *nmax, PetscBool *set) { 2803 const char *svalue; 2804 char *value; 2805 PetscInt n = 0, i, j, start, end, inc, nvalues; 2806 size_t len; 2807 PetscBool flag, foundrange; 2808 PetscToken token; 2809 2810 PetscFunctionBegin; 2811 PetscValidCharPointer(name, 3); 2812 PetscValidIntPointer(ivalue, 4); 2813 PetscValidIntPointer(nmax, 5); 2814 2815 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2816 if (!flag || !svalue) { 2817 if (set) *set = PETSC_FALSE; 2818 *nmax = 0; 2819 PetscFunctionReturn(0); 2820 } 2821 if (set) *set = PETSC_TRUE; 2822 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2823 PetscCall(PetscTokenFind(token, &value)); 2824 while (value && n < *nmax) { 2825 /* look for form d-D where d and D are integers */ 2826 foundrange = PETSC_FALSE; 2827 PetscCall(PetscStrlen(value, &len)); 2828 if (value[0] == '-') i = 2; 2829 else i = 1; 2830 for (; i < (int)len; i++) { 2831 if (value[i] == '-') { 2832 PetscCheck(i != (int)len - 1, PETSC_COMM_SELF, PETSC_ERR_USER, "Error in %" PetscInt_FMT "-th array entry %s", n, value); 2833 value[i] = 0; 2834 2835 PetscCall(PetscOptionsStringToInt(value, &start)); 2836 inc = 1; 2837 j = i + 1; 2838 for (; j < (int)len; j++) { 2839 if (value[j] == ':') { 2840 value[j] = 0; 2841 2842 PetscCall(PetscOptionsStringToInt(value + j + 1, &inc)); 2843 PetscCheck(inc > 0, PETSC_COMM_SELF, PETSC_ERR_USER, "Error in %" PetscInt_FMT "-th array entry,%s cannot have negative increment", n, value + j + 1); 2844 break; 2845 } 2846 } 2847 PetscCall(PetscOptionsStringToInt(value + i + 1, &end)); 2848 PetscCheck(end > start, PETSC_COMM_SELF, PETSC_ERR_USER, "Error in %" PetscInt_FMT "-th array entry, %s-%s cannot have decreasing list", n, value, value + i + 1); 2849 nvalues = (end - start) / inc + (end - start) % inc; 2850 PetscCheck(n + nvalues <= *nmax, PETSC_COMM_SELF, PETSC_ERR_USER, "Error in %" PetscInt_FMT "-th array entry, not enough space left in array (%" PetscInt_FMT ") to contain entire range from %" PetscInt_FMT " to %" PetscInt_FMT, n, *nmax - n, start, end); 2851 for (; start < end; start += inc) { 2852 *ivalue = start; 2853 ivalue++; 2854 n++; 2855 } 2856 foundrange = PETSC_TRUE; 2857 break; 2858 } 2859 } 2860 if (!foundrange) { 2861 PetscCall(PetscOptionsStringToInt(value, ivalue)); 2862 ivalue++; 2863 n++; 2864 } 2865 PetscCall(PetscTokenFind(token, &value)); 2866 } 2867 PetscCall(PetscTokenDestroy(&token)); 2868 *nmax = n; 2869 PetscFunctionReturn(0); 2870 } 2871 2872 /*@C 2873 PetscOptionsGetRealArray - Gets an array of double precision values for a 2874 particular option in the database. The values must be separated with commas with no intervening spaces. 2875 2876 Not Collective 2877 2878 Input Parameters: 2879 + options - options database, use NULL for default global database 2880 . pre - string to prepend to each name or NULL 2881 - name - the option one is seeking 2882 2883 Output Parameters: 2884 + dvalue - the double values to return 2885 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2886 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2887 2888 Level: beginner 2889 2890 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2891 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsBool()`, 2892 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2893 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2894 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2895 `PetscOptionsFList()`, `PetscOptionsEList()` 2896 @*/ 2897 PetscErrorCode PetscOptionsGetRealArray(PetscOptions options, const char pre[], const char name[], PetscReal dvalue[], PetscInt *nmax, PetscBool *set) { 2898 const char *svalue; 2899 char *value; 2900 PetscInt n = 0; 2901 PetscBool flag; 2902 PetscToken token; 2903 2904 PetscFunctionBegin; 2905 PetscValidCharPointer(name, 3); 2906 PetscValidRealPointer(dvalue, 4); 2907 PetscValidIntPointer(nmax, 5); 2908 2909 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2910 if (!flag || !svalue) { 2911 if (set) *set = PETSC_FALSE; 2912 *nmax = 0; 2913 PetscFunctionReturn(0); 2914 } 2915 if (set) *set = PETSC_TRUE; 2916 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2917 PetscCall(PetscTokenFind(token, &value)); 2918 while (value && n < *nmax) { 2919 PetscCall(PetscOptionsStringToReal(value, dvalue++)); 2920 PetscCall(PetscTokenFind(token, &value)); 2921 n++; 2922 } 2923 PetscCall(PetscTokenDestroy(&token)); 2924 *nmax = n; 2925 PetscFunctionReturn(0); 2926 } 2927 2928 /*@C 2929 PetscOptionsGetScalarArray - Gets an array of scalars for a 2930 particular option in the database. The values must be separated with commas with no intervening spaces. 2931 2932 Not Collective 2933 2934 Input Parameters: 2935 + options - options database, use NULL for default global database 2936 . pre - string to prepend to each name or NULL 2937 - name - the option one is seeking 2938 2939 Output Parameters: 2940 + dvalue - the scalar values to return 2941 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2942 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2943 2944 Level: beginner 2945 2946 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2947 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsBool()`, 2948 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2949 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2950 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2951 `PetscOptionsFList()`, `PetscOptionsEList()` 2952 @*/ 2953 PetscErrorCode PetscOptionsGetScalarArray(PetscOptions options, const char pre[], const char name[], PetscScalar dvalue[], PetscInt *nmax, PetscBool *set) { 2954 const char *svalue; 2955 char *value; 2956 PetscInt n = 0; 2957 PetscBool flag; 2958 PetscToken token; 2959 2960 PetscFunctionBegin; 2961 PetscValidCharPointer(name, 3); 2962 PetscValidScalarPointer(dvalue, 4); 2963 PetscValidIntPointer(nmax, 5); 2964 2965 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2966 if (!flag || !svalue) { 2967 if (set) *set = PETSC_FALSE; 2968 *nmax = 0; 2969 PetscFunctionReturn(0); 2970 } 2971 if (set) *set = PETSC_TRUE; 2972 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2973 PetscCall(PetscTokenFind(token, &value)); 2974 while (value && n < *nmax) { 2975 PetscCall(PetscOptionsStringToScalar(value, dvalue++)); 2976 PetscCall(PetscTokenFind(token, &value)); 2977 n++; 2978 } 2979 PetscCall(PetscTokenDestroy(&token)); 2980 *nmax = n; 2981 PetscFunctionReturn(0); 2982 } 2983 2984 /*@C 2985 PetscOptionsGetStringArray - Gets an array of string values for a particular 2986 option in the database. The values must be separated with commas with no intervening spaces. 2987 2988 Not Collective 2989 2990 Input Parameters: 2991 + options - options database, use NULL for default global database 2992 . pre - string to prepend to name or NULL 2993 - name - the option one is seeking 2994 2995 Output Parameters: 2996 + strings - location to copy strings 2997 . nmax - On input maximum number of strings, on output the actual number of strings found 2998 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2999 3000 Level: beginner 3001 3002 Notes: 3003 The nmax parameter is used for both input and output. 3004 3005 The user should pass in an array of pointers to char, to hold all the 3006 strings returned by this function. 3007 3008 The user is responsible for deallocating the strings that are 3009 returned. The Fortran interface for this routine is not supported. 3010 3011 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 3012 `PetscOptionsHasName()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 3013 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 3014 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 3015 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 3016 `PetscOptionsFList()`, `PetscOptionsEList()` 3017 @*/ 3018 PetscErrorCode PetscOptionsGetStringArray(PetscOptions options, const char pre[], const char name[], char *strings[], PetscInt *nmax, PetscBool *set) { 3019 const char *svalue; 3020 char *value; 3021 PetscInt n = 0; 3022 PetscBool flag; 3023 PetscToken token; 3024 3025 PetscFunctionBegin; 3026 PetscValidCharPointer(name, 3); 3027 PetscValidPointer(strings, 4); 3028 PetscValidIntPointer(nmax, 5); 3029 3030 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 3031 if (!flag || !svalue) { 3032 if (set) *set = PETSC_FALSE; 3033 *nmax = 0; 3034 PetscFunctionReturn(0); 3035 } 3036 if (set) *set = PETSC_TRUE; 3037 PetscCall(PetscTokenCreate(svalue, ',', &token)); 3038 PetscCall(PetscTokenFind(token, &value)); 3039 while (value && n < *nmax) { 3040 PetscCall(PetscStrallocpy(value, &strings[n])); 3041 PetscCall(PetscTokenFind(token, &value)); 3042 n++; 3043 } 3044 PetscCall(PetscTokenDestroy(&token)); 3045 *nmax = n; 3046 PetscFunctionReturn(0); 3047 } 3048 3049 /*@C 3050 PetscOptionsDeprecated - mark an option as deprecated, optionally replacing it with a new one 3051 3052 Prints a deprecation warning, unless an option is supplied to suppress. 3053 3054 Logically Collective 3055 3056 Input Parameters: 3057 + pre - string to prepend to name or NULL 3058 . oldname - the old, deprecated option 3059 . newname - the new option, or NULL if option is purely removed 3060 . version - a string describing the version of first deprecation, e.g. "3.9" 3061 - info - additional information string, or NULL. 3062 3063 Options Database Key: 3064 . -options_suppress_deprecated_warnings - do not print deprecation warnings 3065 3066 Notes: 3067 Must be called between `PetscOptionsBegin()` (or `PetscObjectOptionsBegin()`) and `PetscOptionsEnd()`. 3068 Only the proces of rank zero that owns the `PetscOptionsItems` are argument (managed by `PetscOptionsBegin()` or 3069 `PetscObjectOptionsBegin()` prints the information 3070 If newname is provided, the old option is replaced. Otherwise, it remains 3071 in the options database. 3072 If an option is not replaced, the info argument should be used to advise the user 3073 on how to proceed. 3074 There is a limit on the length of the warning printed, so very long strings 3075 provided as info may be truncated. 3076 3077 Level: developer 3078 3079 .seealso: `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsScalar()`, `PetscOptionsBool()`, `PetscOptionsString()`, `PetscOptionsSetValue()` 3080 @*/ 3081 PetscErrorCode PetscOptionsDeprecated_Private(PetscOptionItems *PetscOptionsObject, const char oldname[], const char newname[], const char version[], const char info[]) { 3082 PetscBool found, quiet; 3083 const char *value; 3084 const char *const quietopt = "-options_suppress_deprecated_warnings"; 3085 char msg[4096]; 3086 char *prefix = NULL; 3087 PetscOptions options = NULL; 3088 MPI_Comm comm = PETSC_COMM_SELF; 3089 3090 PetscFunctionBegin; 3091 PetscValidCharPointer(oldname, 2); 3092 PetscValidCharPointer(version, 4); 3093 if (PetscOptionsObject) { 3094 prefix = PetscOptionsObject->prefix; 3095 options = PetscOptionsObject->options; 3096 comm = PetscOptionsObject->comm; 3097 } 3098 PetscCall(PetscOptionsFindPair(options, prefix, oldname, &value, &found)); 3099 if (found) { 3100 if (newname) { 3101 if (prefix) PetscCall(PetscOptionsPrefixPush(options, prefix)); 3102 PetscCall(PetscOptionsSetValue(options, newname, value)); 3103 if (prefix) PetscCall(PetscOptionsPrefixPop(options)); 3104 PetscCall(PetscOptionsClearValue(options, oldname)); 3105 } 3106 quiet = PETSC_FALSE; 3107 PetscCall(PetscOptionsGetBool(options, NULL, quietopt, &quiet, NULL)); 3108 if (!quiet) { 3109 PetscCall(PetscStrcpy(msg, "** PETSc DEPRECATION WARNING ** : the option ")); 3110 PetscCall(PetscStrcat(msg, oldname)); 3111 PetscCall(PetscStrcat(msg, " is deprecated as of version ")); 3112 PetscCall(PetscStrcat(msg, version)); 3113 PetscCall(PetscStrcat(msg, " and will be removed in a future release.")); 3114 if (newname) { 3115 PetscCall(PetscStrcat(msg, " Please use the option ")); 3116 PetscCall(PetscStrcat(msg, newname)); 3117 PetscCall(PetscStrcat(msg, " instead.")); 3118 } 3119 if (info) { 3120 PetscCall(PetscStrcat(msg, " ")); 3121 PetscCall(PetscStrcat(msg, info)); 3122 } 3123 PetscCall(PetscStrcat(msg, " (Silence this warning with ")); 3124 PetscCall(PetscStrcat(msg, quietopt)); 3125 PetscCall(PetscStrcat(msg, ")\n")); 3126 PetscCall(PetscPrintf(comm, "%s", msg)); 3127 } 3128 } 3129 PetscFunctionReturn(0); 3130 } 3131