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