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