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