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()`, `PetscOptionsInsert()`, `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 Note: 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[] = { 912 "malloc_debug", 913 "malloc_dump", 914 "malloc_test", 915 "malloc", 916 "nox", 917 "nox_warning", 918 "display", 919 "saws_port_auto_select", 920 "saws_port_auto_select_silent", 921 "vecscatter_mpi1", 922 "check_pointer_intensity", 923 "cuda_initialize", 924 "error_output_stdout", 925 "use_gpu_aware_mpi", 926 "checkfunctionlist", 927 "fp_trap", 928 "petsc_ci", 929 "petsc_ci_portable_error_output", 930 }; 931 932 static PetscBool PetscCIOption(const char *name) 933 { 934 PetscInt idx; 935 PetscBool found; 936 937 if (!PetscCIEnabled) return PETSC_FALSE; 938 PetscCallAbort(PETSC_COMM_SELF, PetscEListFind(PETSC_STATIC_ARRAY_LENGTH(PetscCIOptions), PetscCIOptions, name, &idx, &found)); 939 return found; 940 } 941 942 /*@C 943 PetscOptionsView - Prints the options that have been loaded. This is 944 useful for debugging purposes. 945 946 Logically Collective 947 948 Input Parameters: 949 + options - options database, use `NULL` for default global database 950 - viewer - must be an `PETSCVIEWERASCII` viewer 951 952 Options Database Key: 953 . -options_view - Activates `PetscOptionsView()` within `PetscFinalize()` 954 955 Level: advanced 956 957 Note: 958 Only the MPI rank 0 of the `MPI_Comm` used to create view prints the option values. Other processes 959 may have different values but they are not printed. 960 961 .seealso: `PetscOptionsAllUsed()` 962 @*/ 963 PetscErrorCode PetscOptionsView(PetscOptions options, PetscViewer viewer) 964 { 965 PetscInt i, N = 0; 966 PetscBool isascii; 967 968 PetscFunctionBegin; 969 if (viewer) PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 2); 970 options = options ? options : defaultoptions; 971 if (!viewer) viewer = PETSC_VIEWER_STDOUT_WORLD; 972 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii)); 973 PetscCheck(isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Only supports ASCII viewer"); 974 975 for (i = 0; i < options->N; i++) { 976 if (PetscCIOption(options->names[i])) continue; 977 N++; 978 } 979 980 if (!N) { 981 PetscCall(PetscViewerASCIIPrintf(viewer, "#No PETSc Option Table entries\n")); 982 PetscFunctionReturn(PETSC_SUCCESS); 983 } 984 985 PetscCall(PetscViewerASCIIPrintf(viewer, "#PETSc Option Table entries:\n")); 986 for (i = 0; i < options->N; i++) { 987 if (PetscCIOption(options->names[i])) continue; 988 if (options->values[i]) { 989 PetscCall(PetscViewerASCIIPrintf(viewer, "-%s %s", options->names[i], options->values[i])); 990 } else { 991 PetscCall(PetscViewerASCIIPrintf(viewer, "-%s", options->names[i])); 992 } 993 PetscCall(PetscViewerASCIIPrintf(viewer, " # (source: %s)\n", PetscOptionSources[options->source[i]])); 994 } 995 PetscCall(PetscViewerASCIIPrintf(viewer, "#End of PETSc Option Table entries\n")); 996 PetscFunctionReturn(PETSC_SUCCESS); 997 } 998 999 /* 1000 Called by error handlers to print options used in run 1001 */ 1002 PetscErrorCode PetscOptionsLeftError(void) 1003 { 1004 PetscInt i, nopt = 0; 1005 1006 for (i = 0; i < defaultoptions->N; i++) { 1007 if (!defaultoptions->used[i]) { 1008 if (PetscCIOption(defaultoptions->names[i])) continue; 1009 nopt++; 1010 } 1011 } 1012 if (nopt) { 1013 PetscCall((*PetscErrorPrintf)("WARNING! There are unused option(s) set! Could be the program crashed before usage or a spelling mistake, etc!\n")); 1014 for (i = 0; i < defaultoptions->N; i++) { 1015 if (!defaultoptions->used[i]) { 1016 if (PetscCIOption(defaultoptions->names[i])) continue; 1017 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]])); 1018 else PetscCall((*PetscErrorPrintf)(" Option left: name:-%s (no value) source: %s\n", defaultoptions->names[i], PetscOptionSources[defaultoptions->source[i]])); 1019 } 1020 } 1021 } 1022 return PETSC_SUCCESS; 1023 } 1024 1025 PETSC_EXTERN PetscErrorCode PetscOptionsViewError(void) 1026 { 1027 PetscInt i, N = 0; 1028 PetscOptions options = defaultoptions; 1029 1030 for (i = 0; i < options->N; i++) { 1031 if (PetscCIOption(options->names[i])) continue; 1032 N++; 1033 } 1034 1035 if (N) { 1036 PetscCall((*PetscErrorPrintf)("PETSc Option Table entries:\n")); 1037 } else { 1038 PetscCall((*PetscErrorPrintf)("No PETSc Option Table entries\n")); 1039 } 1040 for (i = 0; i < options->N; i++) { 1041 if (PetscCIOption(options->names[i])) continue; 1042 if (options->values[i]) { 1043 PetscCall((*PetscErrorPrintf)("-%s %s (source: %s)\n", options->names[i], options->values[i], PetscOptionSources[options->source[i]])); 1044 } else { 1045 PetscCall((*PetscErrorPrintf)("-%s (source: %s)\n", options->names[i], PetscOptionSources[options->source[i]])); 1046 } 1047 } 1048 return PETSC_SUCCESS; 1049 } 1050 1051 /*@C 1052 PetscOptionsPrefixPush - Designate a prefix to be used by all options insertions to follow. 1053 1054 Logically Collective 1055 1056 Input Parameters: 1057 + options - options database, or `NULL` for the default global database 1058 - prefix - The string to append to the existing prefix 1059 1060 Options Database Keys: 1061 + -prefix_push <some_prefix_> - push the given prefix 1062 - -prefix_pop - pop the last prefix 1063 1064 Level: advanced 1065 1066 Notes: 1067 It is common to use this in conjunction with `-options_file` as in 1068 1069 $ -prefix_push system1_ -options_file system1rc -prefix_pop -prefix_push system2_ -options_file system2rc -prefix_pop 1070 1071 where the files no longer require all options to be prefixed with `-system2_`. 1072 1073 The collectivity of this routine is complex; only the MPI processes that call this routine will 1074 have the affect of these options. If some processes that create objects call this routine and others do 1075 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1076 on different ranks. 1077 1078 .seealso: `PetscOptionsPrefixPop()`, `PetscOptionsPush()`, `PetscOptionsPop()`, `PetscOptionsCreate()`, `PetscOptionsSetValue()` 1079 @*/ 1080 PetscErrorCode PetscOptionsPrefixPush(PetscOptions options, const char prefix[]) 1081 { 1082 size_t n; 1083 PetscInt start; 1084 char key[PETSC_MAX_OPTION_NAME + 1]; 1085 PetscBool valid; 1086 1087 PetscFunctionBegin; 1088 PetscValidCharPointer(prefix, 2); 1089 options = options ? options : defaultoptions; 1090 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); 1091 key[0] = '-'; /* keys must start with '-' */ 1092 PetscCall(PetscStrncpy(key + 1, prefix, sizeof(key) - 1)); 1093 PetscCall(PetscOptionsValidKey(key, &valid)); 1094 if (!valid && options->prefixind > 0 && isdigit((int)prefix[0])) valid = PETSC_TRUE; /* If the prefix stack is not empty, make numbers a valid prefix */ 1095 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" : ""); 1096 start = options->prefixind ? options->prefixstack[options->prefixind - 1] : 0; 1097 PetscCall(PetscStrlen(prefix, &n)); 1098 PetscCheck(n + 1 <= sizeof(options->prefix) - start, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Maximum prefix length %zu exceeded", sizeof(options->prefix)); 1099 PetscCall(PetscArraycpy(options->prefix + start, prefix, n + 1)); 1100 options->prefixstack[options->prefixind++] = start + n; 1101 PetscFunctionReturn(PETSC_SUCCESS); 1102 } 1103 1104 /*@C 1105 PetscOptionsPrefixPop - Remove the latest options prefix, see `PetscOptionsPrefixPush()` for details 1106 1107 Logically Collective on the `MPI_Comm` used when called `PetscOptionsPrefixPush()` 1108 1109 Input Parameter: 1110 . options - options database, or `NULL` for the default global database 1111 1112 Level: advanced 1113 1114 .seealso: `PetscOptionsPrefixPush()`, `PetscOptionsPush()`, `PetscOptionsPop()`, `PetscOptionsCreate()`, `PetscOptionsSetValue()` 1115 @*/ 1116 PetscErrorCode PetscOptionsPrefixPop(PetscOptions options) 1117 { 1118 PetscInt offset; 1119 1120 PetscFunctionBegin; 1121 options = options ? options : defaultoptions; 1122 PetscCheck(options->prefixind >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "More prefixes popped than pushed"); 1123 options->prefixind--; 1124 offset = options->prefixind ? options->prefixstack[options->prefixind - 1] : 0; 1125 options->prefix[offset] = 0; 1126 PetscFunctionReturn(PETSC_SUCCESS); 1127 } 1128 1129 /*@C 1130 PetscOptionsClear - Removes all options form the database leaving it empty. 1131 1132 Logically Collective 1133 1134 Input Parameter: 1135 . options - options database, use `NULL` for the default global database 1136 1137 Level: developer 1138 1139 Note: 1140 The collectivity of this routine is complex; only the MPI processes that call this routine will 1141 have the affect of these options. If some processes that create objects call this routine and others do 1142 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1143 on different ranks. 1144 1145 .seealso: `PetscOptionsInsert()` 1146 @*/ 1147 PetscErrorCode PetscOptionsClear(PetscOptions options) 1148 { 1149 PetscInt i; 1150 1151 PetscFunctionBegin; 1152 options = options ? options : defaultoptions; 1153 if (!options) PetscFunctionReturn(PETSC_SUCCESS); 1154 1155 for (i = 0; i < options->N; i++) { 1156 if (options->names[i]) free(options->names[i]); 1157 if (options->values[i]) free(options->values[i]); 1158 } 1159 options->N = 0; 1160 free(options->names); 1161 free(options->values); 1162 free(options->used); 1163 free(options->source); 1164 options->names = NULL; 1165 options->values = NULL; 1166 options->used = NULL; 1167 options->source = NULL; 1168 options->Nalloc = 0; 1169 1170 for (i = 0; i < options->Na; i++) { 1171 free(options->aliases1[i]); 1172 free(options->aliases2[i]); 1173 } 1174 options->Na = 0; 1175 free(options->aliases1); 1176 free(options->aliases2); 1177 options->aliases1 = options->aliases2 = NULL; 1178 options->Naalloc = 0; 1179 1180 /* destroy hash table */ 1181 kh_destroy(HO, options->ht); 1182 options->ht = NULL; 1183 1184 options->prefixind = 0; 1185 options->prefix[0] = 0; 1186 options->help = PETSC_FALSE; 1187 options->help_intro = PETSC_FALSE; 1188 PetscFunctionReturn(PETSC_SUCCESS); 1189 } 1190 1191 /*@C 1192 PetscOptionsSetAlias - Makes a key and alias for another key 1193 1194 Logically Collective 1195 1196 Input Parameters: 1197 + options - options database, or `NULL` for default global database 1198 . newname - the alias 1199 - oldname - the name that alias will refer to 1200 1201 Level: advanced 1202 1203 Note: 1204 The collectivity of this routine is complex; only the MPI processes that call this routine will 1205 have the affect of these options. If some processes that create objects call this routine and others do 1206 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1207 on different ranks. 1208 1209 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, `OptionsHasName()`, 1210 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 1211 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 1212 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 1213 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 1214 `PetscOptionsFList()`, `PetscOptionsEList()` 1215 @*/ 1216 PetscErrorCode PetscOptionsSetAlias(PetscOptions options, const char newname[], const char oldname[]) 1217 { 1218 size_t len; 1219 PetscBool valid; 1220 1221 PetscFunctionBegin; 1222 PetscValidCharPointer(newname, 2); 1223 PetscValidCharPointer(oldname, 3); 1224 options = options ? options : defaultoptions; 1225 PetscCall(PetscOptionsValidKey(newname, &valid)); 1226 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid aliased option %s", newname); 1227 PetscCall(PetscOptionsValidKey(oldname, &valid)); 1228 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid aliasee option %s", oldname); 1229 1230 if (options->Na == options->Naalloc) { 1231 char **tmpA1, **tmpA2; 1232 1233 options->Naalloc = PetscMax(4, options->Naalloc * 2); 1234 tmpA1 = (char **)malloc(options->Naalloc * sizeof(char *)); 1235 tmpA2 = (char **)malloc(options->Naalloc * sizeof(char *)); 1236 for (int i = 0; i < options->Na; ++i) { 1237 tmpA1[i] = options->aliases1[i]; 1238 tmpA2[i] = options->aliases2[i]; 1239 } 1240 free(options->aliases1); 1241 free(options->aliases2); 1242 options->aliases1 = tmpA1; 1243 options->aliases2 = tmpA2; 1244 } 1245 newname++; 1246 oldname++; 1247 PetscCall(PetscStrlen(newname, &len)); 1248 options->aliases1[options->Na] = (char *)malloc((len + 1) * sizeof(char)); 1249 PetscCall(PetscStrncpy(options->aliases1[options->Na], newname, len + 1)); 1250 PetscCall(PetscStrlen(oldname, &len)); 1251 options->aliases2[options->Na] = (char *)malloc((len + 1) * sizeof(char)); 1252 PetscCall(PetscStrncpy(options->aliases2[options->Na], oldname, len + 1)); 1253 ++options->Na; 1254 PetscFunctionReturn(PETSC_SUCCESS); 1255 } 1256 1257 /*@C 1258 PetscOptionsSetValue - Sets an option name-value pair in the options 1259 database, overriding whatever is already present. 1260 1261 Logically Collective 1262 1263 Input Parameters: 1264 + options - options database, use `NULL` for the default global database 1265 . name - name of option, this SHOULD have the - prepended 1266 - value - the option value (not used for all options, so can be `NULL`) 1267 1268 Level: intermediate 1269 1270 Note: 1271 This function can be called BEFORE `PetscInitialize()` 1272 1273 The collectivity of this routine is complex; only the MPI processes that call this routine will 1274 have the affect of these options. If some processes that create objects call this routine and others do 1275 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1276 on different ranks. 1277 1278 Developers Note: 1279 Uses malloc() directly because PETSc may not be initialized yet. 1280 1281 .seealso: `PetscOptionsInsert()`, `PetscOptionsClearValue()` 1282 @*/ 1283 PetscErrorCode PetscOptionsSetValue(PetscOptions options, const char name[], const char value[]) 1284 { 1285 PetscFunctionBegin; 1286 PetscCall(PetscOptionsSetValue_Private(options, name, value, NULL, PETSC_OPT_CODE)); 1287 PetscFunctionReturn(PETSC_SUCCESS); 1288 } 1289 1290 PetscErrorCode PetscOptionsSetValue_Private(PetscOptions options, const char name[], const char value[], int *pos, PetscOptionSource source) 1291 { 1292 size_t len; 1293 int n, i; 1294 char **names; 1295 char fullname[PETSC_MAX_OPTION_NAME] = ""; 1296 PetscBool flg; 1297 1298 PetscFunctionBegin; 1299 if (!options) { 1300 PetscCall(PetscOptionsCreateDefault()); 1301 options = defaultoptions; 1302 } 1303 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "name %s must start with '-'", name); 1304 1305 PetscCall(PetscOptionsSkipPrecedent(options, name, &flg)); 1306 if (flg) PetscFunctionReturn(PETSC_SUCCESS); 1307 1308 name++; /* skip starting dash */ 1309 1310 if (options->prefixind > 0) { 1311 strncpy(fullname, options->prefix, sizeof(fullname)); 1312 fullname[sizeof(fullname) - 1] = 0; 1313 strncat(fullname, name, sizeof(fullname) - strlen(fullname) - 1); 1314 fullname[sizeof(fullname) - 1] = 0; 1315 name = fullname; 1316 } 1317 1318 /* check against aliases */ 1319 for (i = 0; i < options->Na; i++) { 1320 int result = PetscOptNameCmp(options->aliases1[i], name); 1321 if (!result) { 1322 name = options->aliases2[i]; 1323 break; 1324 } 1325 } 1326 1327 /* slow search */ 1328 n = options->N; 1329 names = options->names; 1330 for (i = 0; i < options->N; i++) { 1331 int result = PetscOptNameCmp(names[i], name); 1332 if (!result) { 1333 n = i; 1334 goto setvalue; 1335 } else if (result > 0) { 1336 n = i; 1337 break; 1338 } 1339 } 1340 if (options->N == options->Nalloc) { 1341 char **names, **values; 1342 PetscBool *used; 1343 PetscOptionSource *source; 1344 1345 options->Nalloc = PetscMax(10, options->Nalloc * 2); 1346 names = (char **)malloc(options->Nalloc * sizeof(char *)); 1347 values = (char **)malloc(options->Nalloc * sizeof(char *)); 1348 used = (PetscBool *)malloc(options->Nalloc * sizeof(PetscBool)); 1349 source = (PetscOptionSource *)malloc(options->Nalloc * sizeof(PetscOptionSource)); 1350 for (int i = 0; i < options->N; ++i) { 1351 names[i] = options->names[i]; 1352 values[i] = options->values[i]; 1353 used[i] = options->used[i]; 1354 source[i] = options->source[i]; 1355 } 1356 free(options->names); 1357 free(options->values); 1358 free(options->used); 1359 free(options->source); 1360 options->names = names; 1361 options->values = values; 1362 options->used = used; 1363 options->source = source; 1364 } 1365 1366 /* shift remaining values up 1 */ 1367 for (i = options->N; i > n; i--) { 1368 options->names[i] = options->names[i - 1]; 1369 options->values[i] = options->values[i - 1]; 1370 options->used[i] = options->used[i - 1]; 1371 options->source[i] = options->source[i - 1]; 1372 } 1373 options->names[n] = NULL; 1374 options->values[n] = NULL; 1375 options->used[n] = PETSC_FALSE; 1376 options->source[n] = PETSC_OPT_CODE; 1377 options->N++; 1378 1379 /* destroy hash table */ 1380 kh_destroy(HO, options->ht); 1381 options->ht = NULL; 1382 1383 /* set new name */ 1384 len = strlen(name); 1385 options->names[n] = (char *)malloc((len + 1) * sizeof(char)); 1386 PetscCheck(options->names[n], PETSC_COMM_SELF, PETSC_ERR_MEM, "Failed to allocate option name"); 1387 strcpy(options->names[n], name); 1388 1389 setvalue: 1390 /* set new value */ 1391 if (options->values[n]) free(options->values[n]); 1392 len = value ? strlen(value) : 0; 1393 if (len) { 1394 options->values[n] = (char *)malloc((len + 1) * sizeof(char)); 1395 if (!options->values[n]) return PETSC_ERR_MEM; 1396 strcpy(options->values[n], value); 1397 } else { 1398 options->values[n] = NULL; 1399 } 1400 options->source[n] = source; 1401 1402 /* handle -help so that it can be set from anywhere */ 1403 if (!PetscOptNameCmp(name, "help")) { 1404 options->help = PETSC_TRUE; 1405 options->help_intro = (value && !PetscOptNameCmp(value, "intro")) ? PETSC_TRUE : PETSC_FALSE; 1406 options->used[n] = PETSC_TRUE; 1407 } 1408 1409 PetscCall(PetscOptionsMonitor(options, name, value, source)); 1410 if (pos) *pos = n; 1411 PetscFunctionReturn(PETSC_SUCCESS); 1412 } 1413 1414 /*@C 1415 PetscOptionsClearValue - Clears an option name-value pair in the options 1416 database, overriding whatever is already present. 1417 1418 Logically Collective 1419 1420 Input Parameters: 1421 + options - options database, use `NULL` for the default global database 1422 - name - name of option, this SHOULD have the - prepended 1423 1424 Level: intermediate 1425 1426 Note: 1427 The collectivity of this routine is complex; only the MPI processes that call this routine will 1428 have the affect of these options. If some processes that create objects call this routine and others do 1429 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1430 on different ranks. 1431 1432 .seealso: `PetscOptionsInsert()` 1433 @*/ 1434 PetscErrorCode PetscOptionsClearValue(PetscOptions options, const char name[]) 1435 { 1436 int N, n, i; 1437 char **names; 1438 1439 PetscFunctionBegin; 1440 options = options ? options : defaultoptions; 1441 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Name must begin with '-': Instead %s", name); 1442 if (!PetscOptNameCmp(name, "-help")) options->help = options->help_intro = PETSC_FALSE; 1443 1444 name++; /* skip starting dash */ 1445 1446 /* slow search */ 1447 N = n = options->N; 1448 names = options->names; 1449 for (i = 0; i < N; i++) { 1450 int result = PetscOptNameCmp(names[i], name); 1451 if (!result) { 1452 n = i; 1453 break; 1454 } else if (result > 0) { 1455 n = N; 1456 break; 1457 } 1458 } 1459 if (n == N) PetscFunctionReturn(PETSC_SUCCESS); /* it was not present */ 1460 1461 /* remove name and value */ 1462 if (options->names[n]) free(options->names[n]); 1463 if (options->values[n]) free(options->values[n]); 1464 /* shift remaining values down 1 */ 1465 for (i = n; i < N - 1; i++) { 1466 options->names[i] = options->names[i + 1]; 1467 options->values[i] = options->values[i + 1]; 1468 options->used[i] = options->used[i + 1]; 1469 options->source[i] = options->source[i + 1]; 1470 } 1471 options->N--; 1472 1473 /* destroy hash table */ 1474 kh_destroy(HO, options->ht); 1475 options->ht = NULL; 1476 1477 PetscCall(PetscOptionsMonitor(options, name, NULL, PETSC_OPT_CODE)); 1478 PetscFunctionReturn(PETSC_SUCCESS); 1479 } 1480 1481 /*@C 1482 PetscOptionsFindPair - Gets an option name-value pair from the options database. 1483 1484 Not Collective 1485 1486 Input Parameters: 1487 + options - options database, use `NULL` for the default global database 1488 . pre - the string to prepend to the name or `NULL`, this SHOULD NOT have the "-" prepended 1489 - name - name of option, this SHOULD have the "-" prepended 1490 1491 Output Parameters: 1492 + value - the option value (optional, not used for all options) 1493 - set - whether the option is set (optional) 1494 1495 Level: developer 1496 1497 Note: 1498 Each process may find different values or no value depending on how options were inserted into the database 1499 1500 .seealso: `PetscOptionsSetValue()`, `PetscOptionsClearValue()` 1501 @*/ 1502 PetscErrorCode PetscOptionsFindPair(PetscOptions options, const char pre[], const char name[], const char *value[], PetscBool *set) 1503 { 1504 char buf[PETSC_MAX_OPTION_NAME]; 1505 PetscBool usehashtable = PETSC_TRUE; 1506 PetscBool matchnumbers = PETSC_TRUE; 1507 1508 PetscFunctionBegin; 1509 options = options ? options : defaultoptions; 1510 PetscCheck(!pre || !PetscUnlikely(pre[0] == '-'), PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Prefix cannot begin with '-': Instead %s", pre); 1511 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Name must begin with '-': Instead %s", name); 1512 1513 name++; /* skip starting dash */ 1514 1515 /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */ 1516 if (pre && pre[0]) { 1517 char *ptr = buf; 1518 if (name[0] == '-') { 1519 *ptr++ = '-'; 1520 name++; 1521 } 1522 PetscCall(PetscStrncpy(ptr, pre, buf + sizeof(buf) - ptr)); 1523 PetscCall(PetscStrlcat(buf, name, sizeof(buf))); 1524 name = buf; 1525 } 1526 1527 if (PetscDefined(USE_DEBUG)) { 1528 PetscBool valid; 1529 char key[PETSC_MAX_OPTION_NAME + 1] = "-"; 1530 PetscCall(PetscStrncpy(key + 1, name, sizeof(key) - 1)); 1531 PetscCall(PetscOptionsValidKey(key, &valid)); 1532 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid option '%s' obtained from pre='%s' and name='%s'", key, pre ? pre : "", name); 1533 } 1534 1535 if (!options->ht && usehashtable) { 1536 int i, ret; 1537 khiter_t it; 1538 khash_t(HO) *ht; 1539 ht = kh_init(HO); 1540 PetscCheck(ht, PETSC_COMM_SELF, PETSC_ERR_MEM, "Hash table allocation failed"); 1541 ret = kh_resize(HO, ht, options->N * 2); /* twice the required size to reduce risk of collisions */ 1542 PetscCheck(!ret, PETSC_COMM_SELF, PETSC_ERR_MEM, "Hash table allocation failed"); 1543 for (i = 0; i < options->N; i++) { 1544 it = kh_put(HO, ht, options->names[i], &ret); 1545 PetscCheck(ret == 1, PETSC_COMM_SELF, PETSC_ERR_MEM, "Hash table allocation failed"); 1546 kh_val(ht, it) = i; 1547 } 1548 options->ht = ht; 1549 } 1550 1551 if (usehashtable) { /* fast search */ 1552 khash_t(HO) *ht = options->ht; 1553 khiter_t it = kh_get(HO, ht, name); 1554 if (it != kh_end(ht)) { 1555 int i = kh_val(ht, it); 1556 options->used[i] = PETSC_TRUE; 1557 if (value) *value = options->values[i]; 1558 if (set) *set = PETSC_TRUE; 1559 PetscFunctionReturn(PETSC_SUCCESS); 1560 } 1561 } else { /* slow search */ 1562 int i, N = options->N; 1563 for (i = 0; i < N; i++) { 1564 int result = PetscOptNameCmp(options->names[i], name); 1565 if (!result) { 1566 options->used[i] = PETSC_TRUE; 1567 if (value) *value = options->values[i]; 1568 if (set) *set = PETSC_TRUE; 1569 PetscFunctionReturn(PETSC_SUCCESS); 1570 } else if (result > 0) { 1571 break; 1572 } 1573 } 1574 } 1575 1576 /* 1577 The following block slows down all lookups in the most frequent path (most lookups are unsuccessful). 1578 Maybe this special lookup mode should be enabled on request with a push/pop API. 1579 The feature of matching _%d_ used sparingly in the codebase. 1580 */ 1581 if (matchnumbers) { 1582 int i, j, cnt = 0, locs[16], loce[16]; 1583 /* determine the location and number of all _%d_ in the key */ 1584 for (i = 0; name[i]; i++) { 1585 if (name[i] == '_') { 1586 for (j = i + 1; name[j]; j++) { 1587 if (name[j] >= '0' && name[j] <= '9') continue; 1588 if (name[j] == '_' && j > i + 1) { /* found a number */ 1589 locs[cnt] = i + 1; 1590 loce[cnt++] = j + 1; 1591 } 1592 i = j - 1; 1593 break; 1594 } 1595 } 1596 } 1597 for (i = 0; i < cnt; i++) { 1598 PetscBool found; 1599 char opt[PETSC_MAX_OPTION_NAME + 1] = "-", tmp[PETSC_MAX_OPTION_NAME]; 1600 PetscCall(PetscStrncpy(tmp, name, PetscMin((size_t)(locs[i] + 1), sizeof(tmp)))); 1601 PetscCall(PetscStrlcat(opt, tmp, sizeof(opt))); 1602 PetscCall(PetscStrlcat(opt, name + loce[i], sizeof(opt))); 1603 PetscCall(PetscOptionsFindPair(options, NULL, opt, value, &found)); 1604 if (found) { 1605 if (set) *set = PETSC_TRUE; 1606 PetscFunctionReturn(PETSC_SUCCESS); 1607 } 1608 } 1609 } 1610 1611 if (set) *set = PETSC_FALSE; 1612 PetscFunctionReturn(PETSC_SUCCESS); 1613 } 1614 1615 /* Check whether any option begins with pre+name */ 1616 PETSC_EXTERN PetscErrorCode PetscOptionsFindPairPrefix_Private(PetscOptions options, const char pre[], const char name[], const char *value[], PetscBool *set) 1617 { 1618 char buf[PETSC_MAX_OPTION_NAME]; 1619 int numCnt = 0, locs[16], loce[16]; 1620 1621 PetscFunctionBegin; 1622 options = options ? options : defaultoptions; 1623 PetscCheck(!pre || pre[0] != '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Prefix cannot begin with '-': Instead %s", pre); 1624 PetscCheck(name[0] == '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Name must begin with '-': Instead %s", name); 1625 1626 name++; /* skip starting dash */ 1627 1628 /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */ 1629 if (pre && pre[0]) { 1630 char *ptr = buf; 1631 if (name[0] == '-') { 1632 *ptr++ = '-'; 1633 name++; 1634 } 1635 PetscCall(PetscStrncpy(ptr, pre, sizeof(buf) - ((ptr == buf) ? 0 : 1))); 1636 PetscCall(PetscStrlcat(buf, name, sizeof(buf))); 1637 name = buf; 1638 } 1639 1640 if (PetscDefined(USE_DEBUG)) { 1641 PetscBool valid; 1642 char key[PETSC_MAX_OPTION_NAME + 1] = "-"; 1643 PetscCall(PetscStrncpy(key + 1, name, sizeof(key) - 1)); 1644 PetscCall(PetscOptionsValidKey(key, &valid)); 1645 PetscCheck(valid, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid option '%s' obtained from pre='%s' and name='%s'", key, pre ? pre : "", name); 1646 } 1647 1648 /* determine the location and number of all _%d_ in the key */ 1649 { 1650 int i, j; 1651 for (i = 0; name[i]; i++) { 1652 if (name[i] == '_') { 1653 for (j = i + 1; name[j]; j++) { 1654 if (name[j] >= '0' && name[j] <= '9') continue; 1655 if (name[j] == '_' && j > i + 1) { /* found a number */ 1656 locs[numCnt] = i + 1; 1657 loce[numCnt++] = j + 1; 1658 } 1659 i = j - 1; 1660 break; 1661 } 1662 } 1663 } 1664 } 1665 1666 /* slow search */ 1667 for (int c = -1; c < numCnt; ++c) { 1668 char opt[PETSC_MAX_OPTION_NAME + 2] = ""; 1669 size_t len; 1670 1671 if (c < 0) { 1672 PetscCall(PetscStrncpy(opt, name, sizeof(opt))); 1673 } else { 1674 PetscCall(PetscStrncpy(opt, name, PetscMin((size_t)(locs[c] + 1), sizeof(opt)))); 1675 PetscCall(PetscStrlcat(opt, name + loce[c], sizeof(opt) - 1)); 1676 } 1677 PetscCall(PetscStrlen(opt, &len)); 1678 for (int i = 0; i < options->N; i++) { 1679 PetscBool match; 1680 1681 PetscCall(PetscStrncmp(options->names[i], opt, len, &match)); 1682 if (match) { 1683 options->used[i] = PETSC_TRUE; 1684 if (value) *value = options->values[i]; 1685 if (set) *set = PETSC_TRUE; 1686 PetscFunctionReturn(PETSC_SUCCESS); 1687 } 1688 } 1689 } 1690 1691 if (set) *set = PETSC_FALSE; 1692 PetscFunctionReturn(PETSC_SUCCESS); 1693 } 1694 1695 /*@C 1696 PetscOptionsReject - Generates an error if a certain option is given. 1697 1698 Not Collective 1699 1700 Input Parameters: 1701 + options - options database, use `NULL` for default global database 1702 . pre - the option prefix (may be `NULL`) 1703 . name - the option name one is seeking 1704 - mess - error message (may be `NULL`) 1705 1706 Level: advanced 1707 1708 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, `OptionsHasName()`, 1709 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 1710 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 1711 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 1712 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 1713 `PetscOptionsFList()`, `PetscOptionsEList()` 1714 @*/ 1715 PetscErrorCode PetscOptionsReject(PetscOptions options, const char pre[], const char name[], const char mess[]) 1716 { 1717 PetscBool flag = PETSC_FALSE; 1718 1719 PetscFunctionBegin; 1720 PetscCall(PetscOptionsHasName(options, pre, name, &flag)); 1721 if (flag) { 1722 PetscCheck(!mess || !mess[0], PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Program has disabled option: -%s%s with %s", pre ? pre : "", name + 1, mess); 1723 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Program has disabled option: -%s%s", pre ? pre : "", name + 1); 1724 } 1725 PetscFunctionReturn(PETSC_SUCCESS); 1726 } 1727 1728 /*@C 1729 PetscOptionsHasHelp - Determines whether the "-help" option is in the database. 1730 1731 Not Collective 1732 1733 Input Parameter: 1734 . options - options database, use `NULL` for default global database 1735 1736 Output Parameter: 1737 . set - `PETSC_TRUE` if found else `PETSC_FALSE`. 1738 1739 Level: advanced 1740 1741 .seealso: `PetscOptionsHasName()` 1742 @*/ 1743 PetscErrorCode PetscOptionsHasHelp(PetscOptions options, PetscBool *set) 1744 { 1745 PetscFunctionBegin; 1746 PetscValidBoolPointer(set, 2); 1747 options = options ? options : defaultoptions; 1748 *set = options->help; 1749 PetscFunctionReturn(PETSC_SUCCESS); 1750 } 1751 1752 PetscErrorCode PetscOptionsHasHelpIntro_Internal(PetscOptions options, PetscBool *set) 1753 { 1754 PetscFunctionBegin; 1755 PetscValidBoolPointer(set, 2); 1756 options = options ? options : defaultoptions; 1757 *set = options->help_intro; 1758 PetscFunctionReturn(PETSC_SUCCESS); 1759 } 1760 1761 /*@C 1762 PetscOptionsHasName - Determines whether a certain option is given in the database. This returns true whether the option is a number, string or Boolean, even 1763 if its value is set to false. 1764 1765 Not Collective 1766 1767 Input Parameters: 1768 + options - options database, use `NULL` for default global database 1769 . pre - string to prepend to the name or `NULL` 1770 - name - the option one is seeking 1771 1772 Output Parameter: 1773 . set - `PETSC_TRUE` if found else `PETSC_FALSE`. 1774 1775 Level: beginner 1776 1777 Note: 1778 In many cases you probably want to use `PetscOptionsGetBool()` instead of calling this, to allowing toggling values. 1779 1780 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 1781 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 1782 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 1783 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 1784 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 1785 `PetscOptionsFList()`, `PetscOptionsEList()` 1786 @*/ 1787 PetscErrorCode PetscOptionsHasName(PetscOptions options, const char pre[], const char name[], PetscBool *set) 1788 { 1789 const char *value; 1790 PetscBool flag; 1791 1792 PetscFunctionBegin; 1793 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 1794 if (set) *set = flag; 1795 PetscFunctionReturn(PETSC_SUCCESS); 1796 } 1797 1798 /*@C 1799 PetscOptionsGetAll - Lists all the options the program was run with in a single string. 1800 1801 Not Collective 1802 1803 Input Parameter: 1804 . options - the options database, use `NULL` for the default global database 1805 1806 Output Parameter: 1807 . copts - pointer where string pointer is stored 1808 1809 Level: advanced 1810 1811 Notes: 1812 The array and each entry in the array should be freed with `PetscFree()` 1813 1814 Each process may have different values depending on how the options were inserted into the database 1815 1816 .seealso: `PetscOptionsAllUsed()`, `PetscOptionsView()`, `PetscOptionsPush()`, `PetscOptionsPop()` 1817 @*/ 1818 PetscErrorCode PetscOptionsGetAll(PetscOptions options, char *copts[]) 1819 { 1820 PetscInt i; 1821 size_t len = 1, lent = 0; 1822 char *coptions = NULL; 1823 1824 PetscFunctionBegin; 1825 PetscValidPointer(copts, 2); 1826 options = options ? options : defaultoptions; 1827 /* count the length of the required string */ 1828 for (i = 0; i < options->N; i++) { 1829 PetscCall(PetscStrlen(options->names[i], &lent)); 1830 len += 2 + lent; 1831 if (options->values[i]) { 1832 PetscCall(PetscStrlen(options->values[i], &lent)); 1833 len += 1 + lent; 1834 } 1835 } 1836 PetscCall(PetscMalloc1(len, &coptions)); 1837 coptions[0] = 0; 1838 for (i = 0; i < options->N; i++) { 1839 PetscCall(PetscStrlcat(coptions, "-", len)); 1840 PetscCall(PetscStrlcat(coptions, options->names[i], len)); 1841 PetscCall(PetscStrlcat(coptions, " ", len)); 1842 if (options->values[i]) { 1843 PetscCall(PetscStrlcat(coptions, options->values[i], len)); 1844 PetscCall(PetscStrlcat(coptions, " ", len)); 1845 } 1846 } 1847 *copts = coptions; 1848 PetscFunctionReturn(PETSC_SUCCESS); 1849 } 1850 1851 /*@C 1852 PetscOptionsUsed - Indicates if PETSc has used a particular option set in the database 1853 1854 Not Collective 1855 1856 Input Parameters: 1857 + options - options database, use `NULL` for default global database 1858 - name - string name of option 1859 1860 Output Parameter: 1861 . used - `PETSC_TRUE` if the option was used, otherwise false, including if option was not found in options database 1862 1863 Level: advanced 1864 1865 Note: 1866 The value returned may be different on each process and depends on which options have been processed 1867 on the given process 1868 1869 .seealso: `PetscOptionsView()`, `PetscOptionsLeft()`, `PetscOptionsAllUsed()` 1870 @*/ 1871 PetscErrorCode PetscOptionsUsed(PetscOptions options, const char *name, PetscBool *used) 1872 { 1873 PetscInt i; 1874 1875 PetscFunctionBegin; 1876 PetscValidCharPointer(name, 2); 1877 PetscValidBoolPointer(used, 3); 1878 options = options ? options : defaultoptions; 1879 *used = PETSC_FALSE; 1880 for (i = 0; i < options->N; i++) { 1881 PetscCall(PetscStrcasecmp(options->names[i], name, used)); 1882 if (*used) { 1883 *used = options->used[i]; 1884 break; 1885 } 1886 } 1887 PetscFunctionReturn(PETSC_SUCCESS); 1888 } 1889 1890 /*@ 1891 PetscOptionsAllUsed - Returns a count of the number of options in the 1892 database that have never been selected. 1893 1894 Not Collective 1895 1896 Input Parameter: 1897 . options - options database, use `NULL` for default global database 1898 1899 Output Parameter: 1900 . N - count of options not used 1901 1902 Level: advanced 1903 1904 Note: 1905 The value returned may be different on each process and depends on which options have been processed 1906 on the given process 1907 1908 .seealso: `PetscOptionsView()` 1909 @*/ 1910 PetscErrorCode PetscOptionsAllUsed(PetscOptions options, PetscInt *N) 1911 { 1912 PetscInt i, n = 0; 1913 1914 PetscFunctionBegin; 1915 PetscValidIntPointer(N, 2); 1916 options = options ? options : defaultoptions; 1917 for (i = 0; i < options->N; i++) { 1918 if (!options->used[i]) n++; 1919 } 1920 *N = n; 1921 PetscFunctionReturn(PETSC_SUCCESS); 1922 } 1923 1924 /*@ 1925 PetscOptionsLeft - Prints to screen any options that were set and never used. 1926 1927 Not Collective 1928 1929 Input Parameter: 1930 . options - options database; use `NULL` for default global database 1931 1932 Options Database Key: 1933 . -options_left - activates `PetscOptionsAllUsed()` within `PetscFinalize()` 1934 1935 Level: advanced 1936 1937 Notes: 1938 This is rarely used directly, it is called by `PetscFinalize()` in debug more or if -options_left 1939 is passed otherwise to help users determine possible mistakes in their usage of options. This 1940 only prints values on process zero of `PETSC_COMM_WORLD`. 1941 1942 Other processes depending the objects 1943 used may have different options that are left unused. 1944 1945 .seealso: `PetscOptionsAllUsed()` 1946 @*/ 1947 PetscErrorCode PetscOptionsLeft(PetscOptions options) 1948 { 1949 PetscInt i; 1950 PetscInt cnt = 0; 1951 PetscOptions toptions; 1952 1953 PetscFunctionBegin; 1954 toptions = options ? options : defaultoptions; 1955 for (i = 0; i < toptions->N; i++) { 1956 if (!toptions->used[i]) { 1957 if (PetscCIOption(toptions->names[i])) continue; 1958 if (toptions->values[i]) { 1959 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Option left: name:-%s value: %s source: %s\n", toptions->names[i], toptions->values[i], PetscOptionSources[toptions->source[i]])); 1960 } else { 1961 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Option left: name:-%s (no value) source: %s\n", toptions->names[i], PetscOptionSources[toptions->source[i]])); 1962 } 1963 } 1964 } 1965 if (!options) { 1966 toptions = defaultoptions; 1967 while (toptions->previous) { 1968 cnt++; 1969 toptions = toptions->previous; 1970 } 1971 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)); 1972 } 1973 PetscFunctionReturn(PETSC_SUCCESS); 1974 } 1975 1976 /*@C 1977 PetscOptionsLeftGet - Returns all options that were set and never used. 1978 1979 Not Collective 1980 1981 Input Parameter: 1982 . options - options database, use `NULL` for default global database 1983 1984 Output Parameters: 1985 + N - count of options not used 1986 . names - names of options not used 1987 - values - values of options not used 1988 1989 Level: advanced 1990 1991 Notes: 1992 Users should call `PetscOptionsLeftRestore()` to free the memory allocated in this routine 1993 1994 The value returned may be different on each process and depends on which options have been processed 1995 on the given process 1996 1997 .seealso: `PetscOptionsAllUsed()`, `PetscOptionsLeft()` 1998 @*/ 1999 PetscErrorCode PetscOptionsLeftGet(PetscOptions options, PetscInt *N, char **names[], char **values[]) 2000 { 2001 PetscInt i, n; 2002 2003 PetscFunctionBegin; 2004 if (N) PetscValidIntPointer(N, 2); 2005 if (names) PetscValidPointer(names, 3); 2006 if (values) PetscValidPointer(values, 4); 2007 options = options ? options : defaultoptions; 2008 2009 /* The number of unused PETSc options */ 2010 n = 0; 2011 for (i = 0; i < options->N; i++) { 2012 if (PetscCIOption(options->names[i])) continue; 2013 if (!options->used[i]) n++; 2014 } 2015 if (N) *N = n; 2016 if (names) PetscCall(PetscMalloc1(n, names)); 2017 if (values) PetscCall(PetscMalloc1(n, values)); 2018 2019 n = 0; 2020 if (names || values) { 2021 for (i = 0; i < options->N; i++) { 2022 if (!options->used[i]) { 2023 if (PetscCIOption(options->names[i])) continue; 2024 if (names) (*names)[n] = options->names[i]; 2025 if (values) (*values)[n] = options->values[i]; 2026 n++; 2027 } 2028 } 2029 } 2030 PetscFunctionReturn(PETSC_SUCCESS); 2031 } 2032 2033 /*@C 2034 PetscOptionsLeftRestore - Free memory for the unused PETSc options obtained using `PetscOptionsLeftGet()`. 2035 2036 Not Collective 2037 2038 Input Parameters: 2039 + options - options database, use `NULL` for default global database 2040 . names - names of options not used 2041 - values - values of options not used 2042 2043 Level: advanced 2044 2045 .seealso: `PetscOptionsAllUsed()`, `PetscOptionsLeft()`, `PetscOptionsLeftGet()` 2046 @*/ 2047 PetscErrorCode PetscOptionsLeftRestore(PetscOptions options, PetscInt *N, char **names[], char **values[]) 2048 { 2049 PetscFunctionBegin; 2050 if (N) PetscValidIntPointer(N, 2); 2051 if (names) PetscValidPointer(names, 3); 2052 if (values) PetscValidPointer(values, 4); 2053 if (N) *N = 0; 2054 if (names) PetscCall(PetscFree(*names)); 2055 if (values) PetscCall(PetscFree(*values)); 2056 PetscFunctionReturn(PETSC_SUCCESS); 2057 } 2058 2059 /*@C 2060 PetscOptionsMonitorDefault - Print all options set value events using the supplied `PetscViewer`. 2061 2062 Logically Collective 2063 2064 Input Parameters: 2065 + name - option name string 2066 . value - option value string 2067 . source - The source for the option 2068 - ctx - a `PETSCVIEWERASCII` or `NULL` 2069 2070 Level: intermediate 2071 2072 Notes: 2073 If ctx is `NULL`, `PetscPrintf()` is used. 2074 The first MPI rank in the `PetscViewer` viewer actually prints the values, other 2075 processes may have different values set 2076 2077 If `PetscCIEnabled` then do not print the test harness options 2078 2079 .seealso: `PetscOptionsMonitorSet()` 2080 @*/ 2081 PetscErrorCode PetscOptionsMonitorDefault(const char name[], const char value[], PetscOptionSource source, void *ctx) 2082 { 2083 PetscFunctionBegin; 2084 if (PetscCIOption(name)) PetscFunctionReturn(PETSC_SUCCESS); 2085 2086 if (ctx) { 2087 PetscViewer viewer = (PetscViewer)ctx; 2088 if (!value) { 2089 PetscCall(PetscViewerASCIIPrintf(viewer, "Removing option: %s\n", name)); 2090 } else if (!value[0]) { 2091 PetscCall(PetscViewerASCIIPrintf(viewer, "Setting option: %s (no value) (source: %s)\n", name, PetscOptionSources[source])); 2092 } else { 2093 PetscCall(PetscViewerASCIIPrintf(viewer, "Setting option: %s = %s (source: %s)\n", name, value, PetscOptionSources[source])); 2094 } 2095 } else { 2096 MPI_Comm comm = PETSC_COMM_WORLD; 2097 if (!value) { 2098 PetscCall(PetscPrintf(comm, "Removing option: %s\n", name)); 2099 } else if (!value[0]) { 2100 PetscCall(PetscPrintf(comm, "Setting option: %s (no value) (source: %s)\n", name, PetscOptionSources[source])); 2101 } else { 2102 PetscCall(PetscPrintf(comm, "Setting option: %s = %s (source: %s)\n", name, value, PetscOptionSources[source])); 2103 } 2104 } 2105 PetscFunctionReturn(PETSC_SUCCESS); 2106 } 2107 2108 /*@C 2109 PetscOptionsMonitorSet - Sets an ADDITIONAL function to be called at every method that 2110 modified the PETSc options database. 2111 2112 Not Collective 2113 2114 Input Parameters: 2115 + monitor - pointer to function (if this is `NULL`, it turns off monitoring 2116 . mctx - [optional] context for private data for the 2117 monitor routine (use `NULL` if no context is desired) 2118 - monitordestroy - [optional] routine that frees monitor context 2119 (may be `NULL`) 2120 2121 Calling Sequence of `monitor`: 2122 $ PetscErrorCode monitor(const char name[], const char value[], void *mctx) 2123 + name - option name string 2124 . value - option value string 2125 . source - option source 2126 - mctx - optional monitoring context, as set by `PetscOptionsMonitorSet()` 2127 2128 Calling Sequence of `monitordestroy`: 2129 $ PetscErrorCode monitordestroy(void *cctx) 2130 2131 Options Database Key: 2132 See `PetscInitialize()` for options related to option database monitoring. 2133 2134 Level: intermediate 2135 2136 Notes: 2137 The default is to do nothing. To print the name and value of options 2138 being inserted into the database, use `PetscOptionsMonitorDefault()` as the monitoring routine, 2139 with a null monitoring context. 2140 2141 Several different monitoring routines may be set by calling 2142 `PetscOptionsMonitorSet()` multiple times; all will be called in the 2143 order in which they were set. 2144 2145 .seealso: `PetscOptionsMonitorDefault()`, `PetscInitialize()` 2146 @*/ 2147 PetscErrorCode PetscOptionsMonitorSet(PetscErrorCode (*monitor)(const char name[], const char value[], PetscOptionSource, void *), void *mctx, PetscErrorCode (*monitordestroy)(void **)) 2148 { 2149 PetscOptions options = defaultoptions; 2150 2151 PetscFunctionBegin; 2152 if (options->monitorCancel) PetscFunctionReturn(PETSC_SUCCESS); 2153 PetscCheck(options->numbermonitors < MAXOPTIONSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many PetscOptions monitors set"); 2154 options->monitor[options->numbermonitors] = monitor; 2155 options->monitordestroy[options->numbermonitors] = monitordestroy; 2156 options->monitorcontext[options->numbermonitors++] = (void *)mctx; 2157 PetscFunctionReturn(PETSC_SUCCESS); 2158 } 2159 2160 /* 2161 PetscOptionsStringToBool - Converts string to PetscBool, handles cases like "yes", "no", "true", "false", "0", "1", "off", "on". 2162 Empty string is considered as true. 2163 */ 2164 PetscErrorCode PetscOptionsStringToBool(const char value[], PetscBool *a) 2165 { 2166 PetscBool istrue, isfalse; 2167 size_t len; 2168 2169 PetscFunctionBegin; 2170 /* PetscStrlen() returns 0 for NULL or "" */ 2171 PetscCall(PetscStrlen(value, &len)); 2172 if (!len) { 2173 *a = PETSC_TRUE; 2174 PetscFunctionReturn(PETSC_SUCCESS); 2175 } 2176 PetscCall(PetscStrcasecmp(value, "TRUE", &istrue)); 2177 if (istrue) { 2178 *a = PETSC_TRUE; 2179 PetscFunctionReturn(PETSC_SUCCESS); 2180 } 2181 PetscCall(PetscStrcasecmp(value, "YES", &istrue)); 2182 if (istrue) { 2183 *a = PETSC_TRUE; 2184 PetscFunctionReturn(PETSC_SUCCESS); 2185 } 2186 PetscCall(PetscStrcasecmp(value, "1", &istrue)); 2187 if (istrue) { 2188 *a = PETSC_TRUE; 2189 PetscFunctionReturn(PETSC_SUCCESS); 2190 } 2191 PetscCall(PetscStrcasecmp(value, "on", &istrue)); 2192 if (istrue) { 2193 *a = PETSC_TRUE; 2194 PetscFunctionReturn(PETSC_SUCCESS); 2195 } 2196 PetscCall(PetscStrcasecmp(value, "FALSE", &isfalse)); 2197 if (isfalse) { 2198 *a = PETSC_FALSE; 2199 PetscFunctionReturn(PETSC_SUCCESS); 2200 } 2201 PetscCall(PetscStrcasecmp(value, "NO", &isfalse)); 2202 if (isfalse) { 2203 *a = PETSC_FALSE; 2204 PetscFunctionReturn(PETSC_SUCCESS); 2205 } 2206 PetscCall(PetscStrcasecmp(value, "0", &isfalse)); 2207 if (isfalse) { 2208 *a = PETSC_FALSE; 2209 PetscFunctionReturn(PETSC_SUCCESS); 2210 } 2211 PetscCall(PetscStrcasecmp(value, "off", &isfalse)); 2212 if (isfalse) { 2213 *a = PETSC_FALSE; 2214 PetscFunctionReturn(PETSC_SUCCESS); 2215 } 2216 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Unknown logical value: %s", value); 2217 } 2218 2219 /* 2220 PetscOptionsStringToInt - Converts a string to an integer value. Handles special cases such as "default" and "decide" 2221 */ 2222 PetscErrorCode PetscOptionsStringToInt(const char name[], PetscInt *a) 2223 { 2224 size_t len; 2225 PetscBool decide, tdefault, mouse; 2226 2227 PetscFunctionBegin; 2228 PetscCall(PetscStrlen(name, &len)); 2229 PetscCheck(len, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "character string of length zero has no numerical value"); 2230 2231 PetscCall(PetscStrcasecmp(name, "PETSC_DEFAULT", &tdefault)); 2232 if (!tdefault) PetscCall(PetscStrcasecmp(name, "DEFAULT", &tdefault)); 2233 PetscCall(PetscStrcasecmp(name, "PETSC_DECIDE", &decide)); 2234 if (!decide) PetscCall(PetscStrcasecmp(name, "DECIDE", &decide)); 2235 PetscCall(PetscStrcasecmp(name, "mouse", &mouse)); 2236 2237 if (tdefault) *a = PETSC_DEFAULT; 2238 else if (decide) *a = PETSC_DECIDE; 2239 else if (mouse) *a = -1; 2240 else { 2241 char *endptr; 2242 long strtolval; 2243 2244 strtolval = strtol(name, &endptr, 10); 2245 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); 2246 2247 #if defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE_ATOLL) 2248 (void)strtolval; 2249 *a = atoll(name); 2250 #elif defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE___INT64) 2251 (void)strtolval; 2252 *a = _atoi64(name); 2253 #else 2254 *a = (PetscInt)strtolval; 2255 #endif 2256 } 2257 PetscFunctionReturn(PETSC_SUCCESS); 2258 } 2259 2260 #if defined(PETSC_USE_REAL___FLOAT128) 2261 #include <quadmath.h> 2262 #endif 2263 2264 static PetscErrorCode PetscStrtod(const char name[], PetscReal *a, char **endptr) 2265 { 2266 PetscFunctionBegin; 2267 #if defined(PETSC_USE_REAL___FLOAT128) 2268 *a = strtoflt128(name, endptr); 2269 #else 2270 *a = (PetscReal)strtod(name, endptr); 2271 #endif 2272 PetscFunctionReturn(PETSC_SUCCESS); 2273 } 2274 2275 static PetscErrorCode PetscStrtoz(const char name[], PetscScalar *a, char **endptr, PetscBool *isImaginary) 2276 { 2277 PetscBool hasi = PETSC_FALSE; 2278 char *ptr; 2279 PetscReal strtoval; 2280 2281 PetscFunctionBegin; 2282 PetscCall(PetscStrtod(name, &strtoval, &ptr)); 2283 if (ptr == name) { 2284 strtoval = 1.; 2285 hasi = PETSC_TRUE; 2286 if (name[0] == 'i') { 2287 ptr++; 2288 } else if (name[0] == '+' && name[1] == 'i') { 2289 ptr += 2; 2290 } else if (name[0] == '-' && name[1] == 'i') { 2291 strtoval = -1.; 2292 ptr += 2; 2293 } 2294 } else if (*ptr == 'i') { 2295 hasi = PETSC_TRUE; 2296 ptr++; 2297 } 2298 *endptr = ptr; 2299 *isImaginary = hasi; 2300 if (hasi) { 2301 #if !defined(PETSC_USE_COMPLEX) 2302 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s contains imaginary but complex not supported ", name); 2303 #else 2304 *a = PetscCMPLX(0., strtoval); 2305 #endif 2306 } else { 2307 *a = strtoval; 2308 } 2309 PetscFunctionReturn(PETSC_SUCCESS); 2310 } 2311 2312 /* 2313 Converts a string to PetscReal value. Handles special cases like "default" and "decide" 2314 */ 2315 PetscErrorCode PetscOptionsStringToReal(const char name[], PetscReal *a) 2316 { 2317 size_t len; 2318 PetscBool match; 2319 char *endptr; 2320 2321 PetscFunctionBegin; 2322 PetscCall(PetscStrlen(name, &len)); 2323 PetscCheck(len, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "String of length zero has no numerical value"); 2324 2325 PetscCall(PetscStrcasecmp(name, "PETSC_DEFAULT", &match)); 2326 if (!match) PetscCall(PetscStrcasecmp(name, "DEFAULT", &match)); 2327 if (match) { 2328 *a = PETSC_DEFAULT; 2329 PetscFunctionReturn(PETSC_SUCCESS); 2330 } 2331 2332 PetscCall(PetscStrcasecmp(name, "PETSC_DECIDE", &match)); 2333 if (!match) PetscCall(PetscStrcasecmp(name, "DECIDE", &match)); 2334 if (match) { 2335 *a = PETSC_DECIDE; 2336 PetscFunctionReturn(PETSC_SUCCESS); 2337 } 2338 2339 PetscCall(PetscStrtod(name, a, &endptr)); 2340 PetscCheck((size_t)(endptr - name) == len, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s has no numeric value", name); 2341 PetscFunctionReturn(PETSC_SUCCESS); 2342 } 2343 2344 PetscErrorCode PetscOptionsStringToScalar(const char name[], PetscScalar *a) 2345 { 2346 PetscBool imag1; 2347 size_t len; 2348 PetscScalar val = 0.; 2349 char *ptr = NULL; 2350 2351 PetscFunctionBegin; 2352 PetscCall(PetscStrlen(name, &len)); 2353 PetscCheck(len, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "character string of length zero has no numerical value"); 2354 PetscCall(PetscStrtoz(name, &val, &ptr, &imag1)); 2355 #if defined(PETSC_USE_COMPLEX) 2356 if ((size_t)(ptr - name) < len) { 2357 PetscBool imag2; 2358 PetscScalar val2; 2359 2360 PetscCall(PetscStrtoz(ptr, &val2, &ptr, &imag2)); 2361 if (imag1) PetscCheck(imag2, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s: must specify imaginary component second", name); 2362 val = PetscCMPLX(PetscRealPart(val), PetscImaginaryPart(val2)); 2363 } 2364 #endif 2365 PetscCheck((size_t)(ptr - name) == len, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Input string %s has no numeric value ", name); 2366 *a = val; 2367 PetscFunctionReturn(PETSC_SUCCESS); 2368 } 2369 2370 /*@C 2371 PetscOptionsGetBool - Gets the Logical (true or false) value for a particular 2372 option in the database. 2373 2374 Not Collective 2375 2376 Input Parameters: 2377 + options - options database, use `NULL` for default global database 2378 . pre - the string to prepend to the name or `NULL` 2379 - name - the option one is seeking 2380 2381 Output Parameters: 2382 + ivalue - the logical value to return 2383 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2384 2385 Level: beginner 2386 2387 Notes: 2388 TRUE, true, YES, yes, nostring, and 1 all translate to `PETSC_TRUE` 2389 FALSE, false, NO, no, and 0 all translate to `PETSC_FALSE` 2390 2391 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 2392 is equivalent to -requested_bool true 2393 2394 If the user does not supply the option at all ivalue is NOT changed. Thus 2395 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2396 2397 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, 2398 `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsGetInt()`, `PetscOptionsBool()`, 2399 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2400 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2401 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2402 `PetscOptionsFList()`, `PetscOptionsEList()` 2403 @*/ 2404 PetscErrorCode PetscOptionsGetBool(PetscOptions options, const char pre[], const char name[], PetscBool *ivalue, PetscBool *set) 2405 { 2406 const char *value; 2407 PetscBool flag; 2408 2409 PetscFunctionBegin; 2410 PetscValidCharPointer(name, 3); 2411 if (ivalue) PetscValidBoolPointer(ivalue, 4); 2412 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2413 if (flag) { 2414 if (set) *set = PETSC_TRUE; 2415 PetscCall(PetscOptionsStringToBool(value, &flag)); 2416 if (ivalue) *ivalue = flag; 2417 } else { 2418 if (set) *set = PETSC_FALSE; 2419 } 2420 PetscFunctionReturn(PETSC_SUCCESS); 2421 } 2422 2423 /*@C 2424 PetscOptionsGetEList - Puts a list of option values that a single one may be selected from 2425 2426 Not Collective 2427 2428 Input Parameters: 2429 + options - options database, use `NULL` for default global database 2430 . pre - the string to prepend to the name or `NULL` 2431 . opt - option name 2432 . list - the possible choices (one of these must be selected, anything else is invalid) 2433 - ntext - number of choices 2434 2435 Output Parameters: 2436 + value - the index of the value to return (defaults to zero if the option name is given but no choice is listed) 2437 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2438 2439 Level: intermediate 2440 2441 Notes: 2442 If the user does not supply the option value is NOT changed. Thus 2443 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2444 2445 See `PetscOptionsFList()` for when the choices are given in a `PetscFunctionList` 2446 2447 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 2448 `PetscOptionsHasName()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2449 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2450 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2451 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2452 `PetscOptionsFList()`, `PetscOptionsEList()` 2453 @*/ 2454 PetscErrorCode PetscOptionsGetEList(PetscOptions options, const char pre[], const char opt[], const char *const *list, PetscInt ntext, PetscInt *value, PetscBool *set) 2455 { 2456 size_t alen, len = 0, tlen = 0; 2457 char *svalue; 2458 PetscBool aset, flg = PETSC_FALSE; 2459 PetscInt i; 2460 2461 PetscFunctionBegin; 2462 PetscValidCharPointer(opt, 3); 2463 for (i = 0; i < ntext; i++) { 2464 PetscCall(PetscStrlen(list[i], &alen)); 2465 if (alen > len) len = alen; 2466 tlen += len + 1; 2467 } 2468 len += 5; /* a little extra space for user mistypes */ 2469 PetscCall(PetscMalloc1(len, &svalue)); 2470 PetscCall(PetscOptionsGetString(options, pre, opt, svalue, len, &aset)); 2471 if (aset) { 2472 PetscCall(PetscEListFind(ntext, list, svalue, value, &flg)); 2473 if (!flg) { 2474 char *avail; 2475 2476 PetscCall(PetscMalloc1(tlen, &avail)); 2477 avail[0] = '\0'; 2478 for (i = 0; i < ntext; i++) { 2479 PetscCall(PetscStrlcat(avail, list[i], tlen)); 2480 PetscCall(PetscStrlcat(avail, " ", tlen)); 2481 } 2482 PetscCall(PetscStrtolower(avail)); 2483 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_USER, "Unknown option %s for -%s%s. Available options: %s", svalue, pre ? pre : "", opt + 1, avail); 2484 } 2485 if (set) *set = PETSC_TRUE; 2486 } else if (set) *set = PETSC_FALSE; 2487 PetscCall(PetscFree(svalue)); 2488 PetscFunctionReturn(PETSC_SUCCESS); 2489 } 2490 2491 /*@C 2492 PetscOptionsGetEnum - Gets the enum value for a particular option in the database. 2493 2494 Not Collective 2495 2496 Input Parameters: 2497 + options - options database, use `NULL` for default global database 2498 . pre - option prefix or `NULL` 2499 . opt - option name 2500 - list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null 2501 2502 Output Parameters: 2503 + value - the value to return 2504 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2505 2506 Level: beginner 2507 2508 Notes: 2509 If the user does not supply the option value is NOT changed. Thus 2510 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2511 2512 List is usually something like `PCASMTypes` or some other predefined list of enum names 2513 2514 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, `PetscOptionsGetInt()`, 2515 `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()` 2516 `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`, 2517 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2518 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2519 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2520 `PetscOptionsFList()`, `PetscOptionsEList()`, `PetscOptionsGetEList()`, `PetscOptionsEnum()` 2521 @*/ 2522 PetscErrorCode PetscOptionsGetEnum(PetscOptions options, const char pre[], const char opt[], const char *const *list, PetscEnum *value, PetscBool *set) 2523 { 2524 PetscInt ntext = 0, tval; 2525 PetscBool fset; 2526 2527 PetscFunctionBegin; 2528 PetscValidCharPointer(opt, 3); 2529 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"); 2530 PetscCheck(ntext >= 3, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "List argument must have at least two entries: typename and type prefix"); 2531 ntext -= 3; 2532 PetscCall(PetscOptionsGetEList(options, pre, opt, list, ntext, &tval, &fset)); 2533 /* with PETSC_USE_64BIT_INDICES sizeof(PetscInt) != sizeof(PetscEnum) */ 2534 if (fset) *value = (PetscEnum)tval; 2535 if (set) *set = fset; 2536 PetscFunctionReturn(PETSC_SUCCESS); 2537 } 2538 2539 /*@C 2540 PetscOptionsGetInt - Gets the integer value for a particular option in the database. 2541 2542 Not Collective 2543 2544 Input Parameters: 2545 + options - options database, use `NULL` for default global database 2546 . pre - the string to prepend to the name or `NULL` 2547 - name - the option one is seeking 2548 2549 Output Parameters: 2550 + ivalue - the integer value to return 2551 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2552 2553 Level: beginner 2554 2555 Notes: 2556 If the user does not supply the option ivalue is NOT changed. Thus 2557 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2558 2559 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, 2560 `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()` 2561 `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`, 2562 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2563 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2564 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2565 `PetscOptionsFList()`, `PetscOptionsEList()` 2566 @*/ 2567 PetscErrorCode PetscOptionsGetInt(PetscOptions options, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set) 2568 { 2569 const char *value; 2570 PetscBool flag; 2571 2572 PetscFunctionBegin; 2573 PetscValidCharPointer(name, 3); 2574 PetscValidIntPointer(ivalue, 4); 2575 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2576 if (flag) { 2577 if (!value) { 2578 if (set) *set = PETSC_FALSE; 2579 } else { 2580 if (set) *set = PETSC_TRUE; 2581 PetscCall(PetscOptionsStringToInt(value, ivalue)); 2582 } 2583 } else { 2584 if (set) *set = PETSC_FALSE; 2585 } 2586 PetscFunctionReturn(PETSC_SUCCESS); 2587 } 2588 2589 /*@C 2590 PetscOptionsGetReal - Gets the double precision value for a particular 2591 option in the database. 2592 2593 Not Collective 2594 2595 Input Parameters: 2596 + options - options database, use `NULL` for default global database 2597 . pre - string to prepend to each name or `NULL` 2598 - name - the option one is seeking 2599 2600 Output Parameters: 2601 + dvalue - the double value to return 2602 - set - `PETSC_TRUE` if found, `PETSC_FALSE` if not found 2603 2604 Level: beginner 2605 2606 Note: 2607 If the user does not supply the option dvalue is NOT changed. Thus 2608 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2609 2610 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2611 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2612 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2613 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2614 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2615 `PetscOptionsFList()`, `PetscOptionsEList()` 2616 @*/ 2617 PetscErrorCode PetscOptionsGetReal(PetscOptions options, const char pre[], const char name[], PetscReal *dvalue, PetscBool *set) 2618 { 2619 const char *value; 2620 PetscBool flag; 2621 2622 PetscFunctionBegin; 2623 PetscValidCharPointer(name, 3); 2624 PetscValidRealPointer(dvalue, 4); 2625 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2626 if (flag) { 2627 if (!value) { 2628 if (set) *set = PETSC_FALSE; 2629 } else { 2630 if (set) *set = PETSC_TRUE; 2631 PetscCall(PetscOptionsStringToReal(value, dvalue)); 2632 } 2633 } else { 2634 if (set) *set = PETSC_FALSE; 2635 } 2636 PetscFunctionReturn(PETSC_SUCCESS); 2637 } 2638 2639 /*@C 2640 PetscOptionsGetScalar - Gets the scalar value for a particular 2641 option in the database. 2642 2643 Not Collective 2644 2645 Input Parameters: 2646 + options - options database, use `NULL` for default global database 2647 . pre - string to prepend to each name or `NULL` 2648 - name - the option one is seeking 2649 2650 Output Parameters: 2651 + dvalue - the double value to return 2652 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2653 2654 Level: beginner 2655 2656 Usage: 2657 A complex number 2+3i must be specified with NO spaces 2658 2659 Note: 2660 If the user does not supply the option dvalue is NOT changed. Thus 2661 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2662 2663 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2664 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2665 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2666 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2667 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2668 `PetscOptionsFList()`, `PetscOptionsEList()` 2669 @*/ 2670 PetscErrorCode PetscOptionsGetScalar(PetscOptions options, const char pre[], const char name[], PetscScalar *dvalue, PetscBool *set) 2671 { 2672 const char *value; 2673 PetscBool flag; 2674 2675 PetscFunctionBegin; 2676 PetscValidCharPointer(name, 3); 2677 PetscValidScalarPointer(dvalue, 4); 2678 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2679 if (flag) { 2680 if (!value) { 2681 if (set) *set = PETSC_FALSE; 2682 } else { 2683 #if !defined(PETSC_USE_COMPLEX) 2684 PetscCall(PetscOptionsStringToReal(value, dvalue)); 2685 #else 2686 PetscCall(PetscOptionsStringToScalar(value, dvalue)); 2687 #endif 2688 if (set) *set = PETSC_TRUE; 2689 } 2690 } else { /* flag */ 2691 if (set) *set = PETSC_FALSE; 2692 } 2693 PetscFunctionReturn(PETSC_SUCCESS); 2694 } 2695 2696 /*@C 2697 PetscOptionsGetString - Gets the string value for a particular option in 2698 the database. 2699 2700 Not Collective 2701 2702 Input Parameters: 2703 + options - options database, use `NULL` for default global database 2704 . pre - string to prepend to name or `NULL` 2705 . name - the option one is seeking 2706 - len - maximum length of the string including null termination 2707 2708 Output Parameters: 2709 + string - location to copy string 2710 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2711 2712 Level: beginner 2713 2714 Note: 2715 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` 2716 2717 If the user does not use the option then the string is not changed. Thus 2718 you should ALWAYS initialize the string if you access it without first checking if the set flag is true. 2719 2720 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). 2721 2722 Fortran Note: 2723 The Fortran interface is slightly different from the C/C++ 2724 interface (len is not used). Sample usage in Fortran follows 2725 .vb 2726 character *20 string 2727 PetscErrorCode ierr 2728 PetscBool set 2729 call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s',string,set,ierr) 2730 .ve 2731 2732 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 2733 `PetscOptionsHasName()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2734 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2735 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2736 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2737 `PetscOptionsFList()`, `PetscOptionsEList()` 2738 @*/ 2739 PetscErrorCode PetscOptionsGetString(PetscOptions options, const char pre[], const char name[], char string[], size_t len, PetscBool *set) 2740 { 2741 const char *value; 2742 PetscBool flag; 2743 2744 PetscFunctionBegin; 2745 PetscValidCharPointer(name, 3); 2746 PetscValidCharPointer(string, 4); 2747 PetscCall(PetscOptionsFindPair(options, pre, name, &value, &flag)); 2748 if (!flag) { 2749 if (set) *set = PETSC_FALSE; 2750 } else { 2751 if (set) *set = PETSC_TRUE; 2752 if (value) PetscCall(PetscStrncpy(string, value, len)); 2753 else PetscCall(PetscArrayzero(string, len)); 2754 } 2755 PetscFunctionReturn(PETSC_SUCCESS); 2756 } 2757 2758 char *PetscOptionsGetStringMatlab(PetscOptions options, const char pre[], const char name[]) 2759 { 2760 const char *value; 2761 PetscBool flag; 2762 2763 PetscFunctionBegin; 2764 if (PetscOptionsFindPair(options, pre, name, &value, &flag)) PetscFunctionReturn(NULL); 2765 if (flag) PetscFunctionReturn((char *)value); 2766 PetscFunctionReturn(NULL); 2767 } 2768 2769 /*@C 2770 PetscOptionsGetBoolArray - Gets an array of Logical (true or false) values for a particular 2771 option in the database. The values must be separated with commas with no intervening spaces. 2772 2773 Not Collective 2774 2775 Input Parameters: 2776 + options - options database, use `NULL` for default global database 2777 . pre - string to prepend to each name or `NULL` 2778 - name - the option one is seeking 2779 2780 Output Parameters: 2781 + dvalue - the integer values to return 2782 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2783 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2784 2785 Level: beginner 2786 2787 Note: 2788 TRUE, true, YES, yes, nostring, and 1 all translate to `PETSC_TRUE`. FALSE, false, NO, no, and 0 all translate to `PETSC_FALSE` 2789 2790 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2791 `PetscOptionsGetString()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2792 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2793 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2794 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2795 `PetscOptionsFList()`, `PetscOptionsEList()` 2796 @*/ 2797 PetscErrorCode PetscOptionsGetBoolArray(PetscOptions options, const char pre[], const char name[], PetscBool dvalue[], PetscInt *nmax, PetscBool *set) 2798 { 2799 const char *svalue; 2800 char *value; 2801 PetscInt n = 0; 2802 PetscBool flag; 2803 PetscToken token; 2804 2805 PetscFunctionBegin; 2806 PetscValidCharPointer(name, 3); 2807 PetscValidBoolPointer(dvalue, 4); 2808 PetscValidIntPointer(nmax, 5); 2809 2810 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2811 if (!flag || !svalue) { 2812 if (set) *set = PETSC_FALSE; 2813 *nmax = 0; 2814 PetscFunctionReturn(PETSC_SUCCESS); 2815 } 2816 if (set) *set = PETSC_TRUE; 2817 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2818 PetscCall(PetscTokenFind(token, &value)); 2819 while (value && n < *nmax) { 2820 PetscCall(PetscOptionsStringToBool(value, dvalue)); 2821 PetscCall(PetscTokenFind(token, &value)); 2822 dvalue++; 2823 n++; 2824 } 2825 PetscCall(PetscTokenDestroy(&token)); 2826 *nmax = n; 2827 PetscFunctionReturn(PETSC_SUCCESS); 2828 } 2829 2830 /*@C 2831 PetscOptionsGetEnumArray - Gets an array of enum values for a particular option in the database. 2832 2833 Not Collective 2834 2835 Input Parameters: 2836 + options - options database, use `NULL` for default global database 2837 . pre - option prefix or `NULL` 2838 . name - option name 2839 - list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null 2840 2841 Output Parameters: 2842 + ivalue - the enum values to return 2843 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2844 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2845 2846 Level: beginner 2847 2848 Notes: 2849 The array must be passed as a comma separated list. 2850 2851 There must be no intervening spaces between the values. 2852 2853 list is usually something like `PCASMTypes` or some other predefined list of enum names. 2854 2855 .seealso: `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`, `PetscOptionsGetInt()`, 2856 `PetscOptionsGetEnum()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()` 2857 `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`, `PetscOptionsName()`, 2858 `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, 2859 `PetscOptionsScalar()`, `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2860 `PetscOptionsFList()`, `PetscOptionsEList()`, `PetscOptionsGetEList()`, `PetscOptionsEnum()` 2861 @*/ 2862 PetscErrorCode PetscOptionsGetEnumArray(PetscOptions options, const char pre[], const char name[], const char *const *list, PetscEnum ivalue[], PetscInt *nmax, PetscBool *set) 2863 { 2864 const char *svalue; 2865 char *value; 2866 PetscInt n = 0; 2867 PetscEnum evalue; 2868 PetscBool flag; 2869 PetscToken token; 2870 2871 PetscFunctionBegin; 2872 PetscValidCharPointer(name, 3); 2873 PetscValidPointer(list, 4); 2874 PetscValidPointer(ivalue, 5); 2875 PetscValidIntPointer(nmax, 6); 2876 2877 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2878 if (!flag || !svalue) { 2879 if (set) *set = PETSC_FALSE; 2880 *nmax = 0; 2881 PetscFunctionReturn(PETSC_SUCCESS); 2882 } 2883 if (set) *set = PETSC_TRUE; 2884 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2885 PetscCall(PetscTokenFind(token, &value)); 2886 while (value && n < *nmax) { 2887 PetscCall(PetscEnumFind(list, value, &evalue, &flag)); 2888 PetscCheck(flag, PETSC_COMM_SELF, PETSC_ERR_USER, "Unknown enum value '%s' for -%s%s", svalue, pre ? pre : "", name + 1); 2889 ivalue[n++] = evalue; 2890 PetscCall(PetscTokenFind(token, &value)); 2891 } 2892 PetscCall(PetscTokenDestroy(&token)); 2893 *nmax = n; 2894 PetscFunctionReturn(PETSC_SUCCESS); 2895 } 2896 2897 /*@C 2898 PetscOptionsGetIntArray - Gets an array of integer values for a particular option in the database. 2899 2900 Not Collective 2901 2902 Input Parameters: 2903 + options - options database, use `NULL` for default global database 2904 . pre - string to prepend to each name or `NULL` 2905 - name - the option one is seeking 2906 2907 Output Parameters: 2908 + ivalue - the integer values to return 2909 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 2910 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 2911 2912 Level: beginner 2913 2914 Notes: 2915 The array can be passed as 2916 + a comma separated list - 0,1,2,3,4,5,6,7 2917 . a range (start\-end+1) - 0-8 2918 . a range with given increment (start\-end+1:inc) - 0-7:2 2919 - a combination of values and ranges separated by commas - 0,1-8,8-15:2 2920 2921 There must be no intervening spaces between the values. 2922 2923 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 2924 `PetscOptionsGetString()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 2925 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 2926 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 2927 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 2928 `PetscOptionsFList()`, `PetscOptionsEList()` 2929 @*/ 2930 PetscErrorCode PetscOptionsGetIntArray(PetscOptions options, const char pre[], const char name[], PetscInt ivalue[], PetscInt *nmax, PetscBool *set) 2931 { 2932 const char *svalue; 2933 char *value; 2934 PetscInt n = 0, i, j, start, end, inc, nvalues; 2935 size_t len; 2936 PetscBool flag, foundrange; 2937 PetscToken token; 2938 2939 PetscFunctionBegin; 2940 PetscValidCharPointer(name, 3); 2941 PetscValidIntPointer(ivalue, 4); 2942 PetscValidIntPointer(nmax, 5); 2943 2944 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 2945 if (!flag || !svalue) { 2946 if (set) *set = PETSC_FALSE; 2947 *nmax = 0; 2948 PetscFunctionReturn(PETSC_SUCCESS); 2949 } 2950 if (set) *set = PETSC_TRUE; 2951 PetscCall(PetscTokenCreate(svalue, ',', &token)); 2952 PetscCall(PetscTokenFind(token, &value)); 2953 while (value && n < *nmax) { 2954 /* look for form d-D where d and D are integers */ 2955 foundrange = PETSC_FALSE; 2956 PetscCall(PetscStrlen(value, &len)); 2957 if (value[0] == '-') i = 2; 2958 else i = 1; 2959 for (; i < (int)len; i++) { 2960 if (value[i] == '-') { 2961 PetscCheck(i != (int)len - 1, PETSC_COMM_SELF, PETSC_ERR_USER, "Error in %" PetscInt_FMT "-th array entry %s", n, value); 2962 value[i] = 0; 2963 2964 PetscCall(PetscOptionsStringToInt(value, &start)); 2965 inc = 1; 2966 j = i + 1; 2967 for (; j < (int)len; j++) { 2968 if (value[j] == ':') { 2969 value[j] = 0; 2970 2971 PetscCall(PetscOptionsStringToInt(value + j + 1, &inc)); 2972 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); 2973 break; 2974 } 2975 } 2976 PetscCall(PetscOptionsStringToInt(value + i + 1, &end)); 2977 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); 2978 nvalues = (end - start) / inc + (end - start) % inc; 2979 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); 2980 for (; start < end; start += inc) { 2981 *ivalue = start; 2982 ivalue++; 2983 n++; 2984 } 2985 foundrange = PETSC_TRUE; 2986 break; 2987 } 2988 } 2989 if (!foundrange) { 2990 PetscCall(PetscOptionsStringToInt(value, ivalue)); 2991 ivalue++; 2992 n++; 2993 } 2994 PetscCall(PetscTokenFind(token, &value)); 2995 } 2996 PetscCall(PetscTokenDestroy(&token)); 2997 *nmax = n; 2998 PetscFunctionReturn(PETSC_SUCCESS); 2999 } 3000 3001 /*@C 3002 PetscOptionsGetRealArray - Gets an array of double precision values for a 3003 particular option in the database. The values must be separated with commas with no intervening spaces. 3004 3005 Not Collective 3006 3007 Input Parameters: 3008 + options - options database, use `NULL` for default global database 3009 . pre - string to prepend to each name or `NULL` 3010 - name - the option one is seeking 3011 3012 Output Parameters: 3013 + dvalue - the double values to return 3014 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 3015 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 3016 3017 Level: beginner 3018 3019 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 3020 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsBool()`, 3021 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 3022 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 3023 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 3024 `PetscOptionsFList()`, `PetscOptionsEList()` 3025 @*/ 3026 PetscErrorCode PetscOptionsGetRealArray(PetscOptions options, const char pre[], const char name[], PetscReal dvalue[], PetscInt *nmax, PetscBool *set) 3027 { 3028 const char *svalue; 3029 char *value; 3030 PetscInt n = 0; 3031 PetscBool flag; 3032 PetscToken token; 3033 3034 PetscFunctionBegin; 3035 PetscValidCharPointer(name, 3); 3036 PetscValidRealPointer(dvalue, 4); 3037 PetscValidIntPointer(nmax, 5); 3038 3039 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 3040 if (!flag || !svalue) { 3041 if (set) *set = PETSC_FALSE; 3042 *nmax = 0; 3043 PetscFunctionReturn(PETSC_SUCCESS); 3044 } 3045 if (set) *set = PETSC_TRUE; 3046 PetscCall(PetscTokenCreate(svalue, ',', &token)); 3047 PetscCall(PetscTokenFind(token, &value)); 3048 while (value && n < *nmax) { 3049 PetscCall(PetscOptionsStringToReal(value, dvalue++)); 3050 PetscCall(PetscTokenFind(token, &value)); 3051 n++; 3052 } 3053 PetscCall(PetscTokenDestroy(&token)); 3054 *nmax = n; 3055 PetscFunctionReturn(PETSC_SUCCESS); 3056 } 3057 3058 /*@C 3059 PetscOptionsGetScalarArray - Gets an array of scalars for a 3060 particular option in the database. The values must be separated with commas with no intervening spaces. 3061 3062 Not Collective 3063 3064 Input Parameters: 3065 + options - options database, use `NULL` for default global database 3066 . pre - string to prepend to each name or `NULL` 3067 - name - the option one is seeking 3068 3069 Output Parameters: 3070 + dvalue - the scalar values to return 3071 . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved 3072 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 3073 3074 Level: beginner 3075 3076 .seealso: `PetscOptionsGetInt()`, `PetscOptionsHasName()`, 3077 `PetscOptionsGetString()`, `PetscOptionsGetIntArray()`, `PetscOptionsBool()`, 3078 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 3079 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 3080 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 3081 `PetscOptionsFList()`, `PetscOptionsEList()` 3082 @*/ 3083 PetscErrorCode PetscOptionsGetScalarArray(PetscOptions options, const char pre[], const char name[], PetscScalar dvalue[], PetscInt *nmax, PetscBool *set) 3084 { 3085 const char *svalue; 3086 char *value; 3087 PetscInt n = 0; 3088 PetscBool flag; 3089 PetscToken token; 3090 3091 PetscFunctionBegin; 3092 PetscValidCharPointer(name, 3); 3093 PetscValidScalarPointer(dvalue, 4); 3094 PetscValidIntPointer(nmax, 5); 3095 3096 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 3097 if (!flag || !svalue) { 3098 if (set) *set = PETSC_FALSE; 3099 *nmax = 0; 3100 PetscFunctionReturn(PETSC_SUCCESS); 3101 } 3102 if (set) *set = PETSC_TRUE; 3103 PetscCall(PetscTokenCreate(svalue, ',', &token)); 3104 PetscCall(PetscTokenFind(token, &value)); 3105 while (value && n < *nmax) { 3106 PetscCall(PetscOptionsStringToScalar(value, dvalue++)); 3107 PetscCall(PetscTokenFind(token, &value)); 3108 n++; 3109 } 3110 PetscCall(PetscTokenDestroy(&token)); 3111 *nmax = n; 3112 PetscFunctionReturn(PETSC_SUCCESS); 3113 } 3114 3115 /*@C 3116 PetscOptionsGetStringArray - Gets an array of string values for a particular 3117 option in the database. The values must be separated with commas with no intervening spaces. 3118 3119 Not Collective; No Fortran Support 3120 3121 Input Parameters: 3122 + options - options database, use `NULL` for default global database 3123 . pre - string to prepend to name or `NULL` 3124 - name - the option one is seeking 3125 3126 Output Parameters: 3127 + strings - location to copy strings 3128 . nmax - On input maximum number of strings, on output the actual number of strings found 3129 - set - `PETSC_TRUE` if found, else `PETSC_FALSE` 3130 3131 Level: beginner 3132 3133 Notes: 3134 The nmax parameter is used for both input and output. 3135 3136 The user should pass in an array of pointers to char, to hold all the 3137 strings returned by this function. 3138 3139 The user is responsible for deallocating the strings that are 3140 returned. 3141 3142 .seealso: `PetscOptionsGetInt()`, `PetscOptionsGetReal()`, 3143 `PetscOptionsHasName()`, `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`, 3144 `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`, 3145 `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`, 3146 `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`, 3147 `PetscOptionsFList()`, `PetscOptionsEList()` 3148 @*/ 3149 PetscErrorCode PetscOptionsGetStringArray(PetscOptions options, const char pre[], const char name[], char *strings[], PetscInt *nmax, PetscBool *set) 3150 { 3151 const char *svalue; 3152 char *value; 3153 PetscInt n = 0; 3154 PetscBool flag; 3155 PetscToken token; 3156 3157 PetscFunctionBegin; 3158 PetscValidCharPointer(name, 3); 3159 PetscValidPointer(strings, 4); 3160 PetscValidIntPointer(nmax, 5); 3161 3162 PetscCall(PetscOptionsFindPair(options, pre, name, &svalue, &flag)); 3163 if (!flag || !svalue) { 3164 if (set) *set = PETSC_FALSE; 3165 *nmax = 0; 3166 PetscFunctionReturn(PETSC_SUCCESS); 3167 } 3168 if (set) *set = PETSC_TRUE; 3169 PetscCall(PetscTokenCreate(svalue, ',', &token)); 3170 PetscCall(PetscTokenFind(token, &value)); 3171 while (value && n < *nmax) { 3172 PetscCall(PetscStrallocpy(value, &strings[n])); 3173 PetscCall(PetscTokenFind(token, &value)); 3174 n++; 3175 } 3176 PetscCall(PetscTokenDestroy(&token)); 3177 *nmax = n; 3178 PetscFunctionReturn(PETSC_SUCCESS); 3179 } 3180 3181 /*@C 3182 PetscOptionsDeprecated - mark an option as deprecated, optionally replacing it with `newname` 3183 3184 Prints a deprecation warning, unless an option is supplied to suppress. 3185 3186 Logically Collective 3187 3188 Input Parameters: 3189 + pre - string to prepend to name or `NULL` 3190 . oldname - the old, deprecated option 3191 . newname - the new option, or `NULL` if option is purely removed 3192 . version - a string describing the version of first deprecation, e.g. "3.9" 3193 - info - additional information string, or `NULL`. 3194 3195 Options Database Key: 3196 . -options_suppress_deprecated_warnings - do not print deprecation warnings 3197 3198 Level: developer 3199 3200 Notes: 3201 If `newname` is provided then the options database will automatically check the database for `oldname`. 3202 3203 The old call `PetscOptionsXXX`(`oldname`) should be removed from the source code when both (1) the call to `PetscOptionsDeprecated()` occurs before the 3204 new call to `PetscOptionsXXX`(`newname`) and (2) the argument handling of the new call to `PetscOptionsXXX`(`newname`) is identical to the previous call. 3205 See `PTScotch_PartGraph_Seq()` for an example of when (1) fails and `SNESTestJacobian()` where an example of (2) fails. 3206 3207 Must be called between `PetscOptionsBegin()` (or `PetscObjectOptionsBegin()`) and `PetscOptionsEnd()`. 3208 Only the process of rank zero that owns the `PetscOptionsItems` are argument (managed by `PetscOptionsBegin()` or 3209 `PetscObjectOptionsBegin()` prints the information 3210 If newname is provided, the old option is replaced. Otherwise, it remains 3211 in the options database. 3212 If an option is not replaced, the info argument should be used to advise the user 3213 on how to proceed. 3214 There is a limit on the length of the warning printed, so very long strings 3215 provided as info may be truncated. 3216 3217 .seealso: `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsScalar()`, `PetscOptionsBool()`, `PetscOptionsString()`, `PetscOptionsSetValue()` 3218 @*/ 3219 PetscErrorCode PetscOptionsDeprecated_Private(PetscOptionItems *PetscOptionsObject, const char oldname[], const char newname[], const char version[], const char info[]) 3220 { 3221 PetscBool found, quiet; 3222 const char *value; 3223 const char *const quietopt = "-options_suppress_deprecated_warnings"; 3224 char msg[4096]; 3225 char *prefix = NULL; 3226 PetscOptions options = NULL; 3227 MPI_Comm comm = PETSC_COMM_SELF; 3228 3229 PetscFunctionBegin; 3230 PetscValidCharPointer(oldname, 2); 3231 PetscValidCharPointer(version, 4); 3232 if (PetscOptionsObject) { 3233 prefix = PetscOptionsObject->prefix; 3234 options = PetscOptionsObject->options; 3235 comm = PetscOptionsObject->comm; 3236 } 3237 PetscCall(PetscOptionsFindPair(options, prefix, oldname, &value, &found)); 3238 if (found) { 3239 if (newname) { 3240 if (prefix) PetscCall(PetscOptionsPrefixPush(options, prefix)); 3241 PetscCall(PetscOptionsSetValue(options, newname, value)); 3242 if (prefix) PetscCall(PetscOptionsPrefixPop(options)); 3243 PetscCall(PetscOptionsClearValue(options, oldname)); 3244 } 3245 quiet = PETSC_FALSE; 3246 PetscCall(PetscOptionsGetBool(options, NULL, quietopt, &quiet, NULL)); 3247 if (!quiet) { 3248 PetscCall(PetscStrncpy(msg, "** PETSc DEPRECATION WARNING ** : the option ", sizeof(msg))); 3249 PetscCall(PetscStrlcat(msg, oldname, sizeof(msg))); 3250 PetscCall(PetscStrlcat(msg, " is deprecated as of version ", sizeof(msg))); 3251 PetscCall(PetscStrlcat(msg, version, sizeof(msg))); 3252 PetscCall(PetscStrlcat(msg, " and will be removed in a future release.\n", sizeof(msg))); 3253 if (newname) { 3254 PetscCall(PetscStrlcat(msg, " Use the option ", sizeof(msg))); 3255 PetscCall(PetscStrlcat(msg, newname, sizeof(msg))); 3256 PetscCall(PetscStrlcat(msg, " instead.", sizeof(msg))); 3257 } 3258 if (info) { 3259 PetscCall(PetscStrlcat(msg, " ", sizeof(msg))); 3260 PetscCall(PetscStrlcat(msg, info, sizeof(msg))); 3261 } 3262 PetscCall(PetscStrlcat(msg, " (Silence this warning with ", sizeof(msg))); 3263 PetscCall(PetscStrlcat(msg, quietopt, sizeof(msg))); 3264 PetscCall(PetscStrlcat(msg, ")\n", sizeof(msg))); 3265 PetscCall(PetscPrintf(comm, "%s", msg)); 3266 } 3267 } 3268 PetscFunctionReturn(PETSC_SUCCESS); 3269 } 3270