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