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