Actual source code: options.c
petsc-master 2021-01-18
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
21: #if defined(PETSC_HAVE_YAML)
22: #include <yaml.h>
23: #endif
25: #if defined(PETSC_HAVE_STRCASECMP)
26: #define PetscOptNameCmp(a,b) strcasecmp(a,b)
27: #elif defined(PETSC_HAVE_STRICMP)
28: #define PetscOptNameCmp(a,b) stricmp(a,b)
29: #else
30: #define PetscOptNameCmp(a,b) Error_strcasecmp_not_found
31: #endif
33: #include <petsc/private/hashtable.h>
35: /* This assumes ASCII encoding and ignores locale settings */
36: /* Using tolower() is about 2X slower in microbenchmarks */
37: PETSC_STATIC_INLINE int PetscToLower(int c)
38: {
39: return ((c >= 'A') & (c <= 'Z')) ? c + 'a' - 'A' : c;
40: }
42: /* Bob Jenkins's one at a time hash function (case-insensitive) */
43: PETSC_STATIC_INLINE unsigned int PetscOptHash(const char key[])
44: {
45: unsigned int hash = 0;
46: while (*key) {
47: hash += PetscToLower(*key++);
48: hash += hash << 10;
49: hash ^= hash >> 6;
50: }
51: hash += hash << 3;
52: hash ^= hash >> 11;
53: hash += hash << 15;
54: return hash;
55: }
57: PETSC_STATIC_INLINE int PetscOptEqual(const char a[],const char b[])
58: {
59: return !PetscOptNameCmp(a,b);
60: }
62: KHASH_INIT(HO, kh_cstr_t, int, 1, PetscOptHash, PetscOptEqual)
64: /*
65: This table holds all the options set by the user. For simplicity, we use a static size database
66: */
67: #define MAXOPTNAME 512
68: #define MAXOPTIONS 512
69: #define MAXALIASES 25
70: #define MAXPREFIXES 25
71: #define MAXOPTIONSMONITORS 5
73: struct _n_PetscOptions {
74: PetscOptions previous;
75: int N; /* number of options */
76: char *names[MAXOPTIONS]; /* option names */
77: char *values[MAXOPTIONS]; /* option values */
78: PetscBool used[MAXOPTIONS]; /* flag option use */
79: PetscBool precedentProcessed;
81: /* Hash table */
82: khash_t(HO) *ht;
84: /* Prefixes */
85: int prefixind;
86: int prefixstack[MAXPREFIXES];
87: char prefix[MAXOPTNAME];
89: /* Aliases */
90: int Naliases; /* number or aliases */
91: char *aliases1[MAXALIASES]; /* aliased */
92: char *aliases2[MAXALIASES]; /* aliasee */
94: /* Help */
95: PetscBool help; /* flag whether "-help" is in the database */
96: PetscBool help_intro; /* flag whether "-help intro" is in the database */
98: /* Monitors */
99: PetscBool monitorFromOptions, monitorCancel;
100: PetscErrorCode (*monitor[MAXOPTIONSMONITORS])(const char[],const char[],void*); /* returns control to user after */
101: PetscErrorCode (*monitordestroy[MAXOPTIONSMONITORS])(void**); /* */
102: void *monitorcontext[MAXOPTIONSMONITORS]; /* to pass arbitrary user data into monitor */
103: PetscInt numbermonitors; /* to, for instance, detect options being set */
104: };
106: static PetscOptions defaultoptions = NULL; /* the options database routines query this object for options */
108: /* list of options which preceed others, i.e., are processed in PetscOptionsProcessPrecedentFlags() */
109: static const char *precedentOptions[] = {"-options_monitor","-options_monitor_cancel","-help","-skip_petscrc","-options_file_yaml","-options_string_yaml"};
110: enum PetscPrecedentOption {PO_OPTIONS_MONITOR,PO_OPTIONS_MONITOR_CANCEL,PO_HELP,PO_SKIP_PETSCRC,PO_OPTIONS_FILE_YAML,PO_OPTIONS_STRING_YAML,PO_NUM};
112: static PetscErrorCode PetscOptionsSetValue_Private(PetscOptions,const char[],const char[],int*);
114: /*
115: Options events monitor
116: */
117: static PetscErrorCode PetscOptionsMonitor(PetscOptions options,const char name[],const char value[])
118: {
119: PetscInt i;
122: if (!PetscErrorHandlingInitialized) return 0;
124: if (!value) value = "";
125: if (options->monitorFromOptions) {
126: PetscOptionsMonitorDefault(name,value,NULL);
127: }
128: for (i=0; i<options->numbermonitors; i++) {
129: (*options->monitor[i])(name,value,options->monitorcontext[i]);
130: }
131: return(0);
132: }
134: /*@
135: PetscOptionsCreate - Creates an empty options database.
137: Logically collective
139: Output Parameter:
140: . options - Options database object
142: Level: advanced
144: Developer Note: We may want eventually to pass a MPI_Comm to determine the ownership of the object
146: .seealso: PetscOptionsDestroy(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsInsert(), PetscOptionsSetValue()
147: @*/
148: PetscErrorCode PetscOptionsCreate(PetscOptions *options)
149: {
150: if (!options) return PETSC_ERR_ARG_NULL;
151: *options = (PetscOptions)calloc(1,sizeof(**options));
152: if (!*options) return PETSC_ERR_MEM;
153: return 0;
154: }
156: /*@
157: PetscOptionsDestroy - Destroys an option database.
159: Logically collective on whatever communicator was associated with the call to PetscOptionsCreate()
161: Input Parameter:
162: . options - the PetscOptions object
164: Level: advanced
166: .seealso: PetscOptionsInsert(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsInsert(), PetscOptionsSetValue()
167: @*/
168: PetscErrorCode PetscOptionsDestroy(PetscOptions *options)
169: {
172: if (!*options) return 0;
173: if ((*options)->previous) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"You are destroying an option that has been used with PetscOptionsPush() but does not have a corresponding PetscOptionsPop()");
174: PetscOptionsClear(*options);if (ierr) return ierr;
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: {
188: if (!defaultoptions) {
189: PetscOptionsCreate(&defaultoptions);if (ierr) return ierr;
190: }
191: return 0;
192: }
194: /*@
195: PetscOptionsPush - Push a new PetscOptions object as the default provider of options
196: Allows using different parts of a code to use different options databases
198: Logically Collective
200: Input Parameter:
201: . opt - the options obtained with PetscOptionsCreate()
203: Notes:
204: Use PetscOptionsPop() to return to the previous default options database
206: The collectivity of this routine is complex; only the MPI processes that call this routine will
207: have the affect of these options. If some processes that create objects call this routine and others do
208: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
209: on different ranks.
211: Level: advanced
213: .seealso: PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsInsert(), PetscOptionsSetValue(), PetscOptionsLeft()
215: @*/
216: PetscErrorCode PetscOptionsPush(PetscOptions opt)
217: {
221: PetscOptionsCreateDefault();
222: opt->previous = defaultoptions;
223: defaultoptions = opt;
224: return(0);
225: }
227: /*@
228: PetscOptionsPop - Pop the most recent PetscOptionsPush() to return to the previous default options
230: Logically collective on whatever communicator was associated with the call to PetscOptionsCreate()
232: Notes:
233: Use PetscOptionsPop() to return to the previous default options database
234: Allows using different parts of a code to use different options databases
236: Level: advanced
238: .seealso: PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsInsert(), PetscOptionsSetValue(), PetscOptionsLeft()
240: @*/
241: PetscErrorCode PetscOptionsPop(void)
242: {
243: PetscOptions current = defaultoptions;
246: if (!defaultoptions) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Missing default options");
247: if (!defaultoptions->previous) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscOptionsPop() called too many times");
248: defaultoptions = defaultoptions->previous;
249: current->previous = NULL;
250: return(0);
251: }
253: /*
254: PetscOptionsDestroyDefault - Destroys the default global options database
255: */
256: PetscErrorCode PetscOptionsDestroyDefault(void)
257: {
259: PetscOptions tmp;
261: /* Destroy any options that the user forgot to pop */
262: while (defaultoptions->previous) {
263: tmp = defaultoptions;
264: PetscOptionsPop();
265: PetscOptionsDestroy(&tmp);
266: }
267: PetscOptionsDestroy(&defaultoptions);if (ierr) return ierr;
268: return 0;
269: }
271: /*@C
272: PetscOptionsValidKey - PETSc Options database keys must begin with one or two dashes (-) followed by a letter.
274: Not collective
276: Input Parameter:
277: . key - string to check if valid
279: Output Parameter:
280: . valid - PETSC_TRUE if a valid key
282: Level: intermediate
283: @*/
284: PetscErrorCode PetscOptionsValidKey(const char key[],PetscBool *valid)
285: {
286: char *ptr;
291: *valid = PETSC_FALSE;
292: if (!key) return(0);
293: if (key[0] != '-') return(0);
294: if (key[1] == '-') key++;
295: if (!isalpha((int)key[1])) return(0);
296: (void) strtod(key,&ptr);
297: if (ptr != key && !(*ptr == '_' || isalnum((int)*ptr))) return(0);
298: *valid = PETSC_TRUE;
299: return(0);
300: }
302: /*@C
303: PetscOptionsInsertString - Inserts options into the database from a string
305: Logically Collective
307: Input Parameter:
308: + options - options object
309: - in_str - string that contains options separated by blanks
311: Level: intermediate
313: The collectivity of this routine is complex; only the MPI processes that call this routine will
314: have the affect of these options. If some processes that create objects call this routine and others do
315: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
316: on different ranks.
318: Contributed by Boyana Norris
320: .seealso: PetscOptionsSetValue(), PetscOptionsView(), PetscOptionsHasName(), PetscOptionsGetInt(),
321: PetscOptionsGetReal(), PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
322: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
323: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
324: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
325: PetscOptionsFList(), PetscOptionsEList(), PetscOptionsInsertFile()
326: @*/
327: PetscErrorCode PetscOptionsInsertString(PetscOptions options,const char in_str[])
328: {
329: char *first,*second;
331: PetscToken token;
332: PetscBool key,ispush,ispop,isopts;
335: PetscTokenCreate(in_str,' ',&token);
336: PetscTokenFind(token,&first);
337: while (first) {
338: PetscStrcasecmp(first,"-prefix_push",&ispush);
339: PetscStrcasecmp(first,"-prefix_pop",&ispop);
340: PetscStrcasecmp(first,"-options_file",&isopts);
341: PetscOptionsValidKey(first,&key);
342: if (ispush) {
343: PetscTokenFind(token,&second);
344: PetscOptionsPrefixPush(options,second);
345: PetscTokenFind(token,&first);
346: } else if (ispop) {
347: PetscOptionsPrefixPop(options);
348: PetscTokenFind(token,&first);
349: } else if (isopts) {
350: PetscTokenFind(token,&second);
351: PetscOptionsInsertFile(PETSC_COMM_SELF,options,second,PETSC_TRUE);
352: PetscTokenFind(token,&first);
353: } else if (key) {
354: PetscTokenFind(token,&second);
355: PetscOptionsValidKey(second,&key);
356: if (!key) {
357: PetscOptionsSetValue(options,first,second);
358: PetscTokenFind(token,&first);
359: } else {
360: PetscOptionsSetValue(options,first,NULL);
361: first = second;
362: }
363: } else {
364: PetscTokenFind(token,&first);
365: }
366: }
367: PetscTokenDestroy(&token);
368: return(0);
369: }
371: /*
372: Returns a line (ended by a \n, \r or null character of any length. Result should be freed with free()
373: */
374: static char *Petscgetline(FILE * f)
375: {
376: size_t size = 0;
377: size_t len = 0;
378: size_t last = 0;
379: char *buf = NULL;
381: if (feof(f)) return NULL;
382: do {
383: size += 1024; /* BUFSIZ is defined as "the optimal read size for this platform" */
384: buf = (char*)realloc((void*)buf,size); /* realloc(NULL,n) is the same as malloc(n) */
385: /* Actually do the read. Note that fgets puts a terminal '\0' on the
386: end of the string, so we make sure we overwrite this */
387: if (!fgets(buf+len,1024,f)) buf[len]=0;
388: PetscStrlen(buf,&len);
389: last = len - 1;
390: } while (!feof(f) && buf[last] != '\n' && buf[last] != '\r');
391: if (len) return buf;
392: free(buf);
393: return NULL;
394: }
396: /*@C
397: PetscOptionsInsertFile - Inserts options into the database from a file.
399: Collective
401: Input Parameter:
402: + comm - the processes that will share the options (usually PETSC_COMM_WORLD)
403: . options - options database, use NULL for default global database
404: . file - name of file
405: - require - if PETSC_TRUE will generate an error if the file does not exist
408: Notes:
409: Use # for lines that are comments and which should be ignored.
410: Usually, instead of using this command, one should list the file name in the call to PetscInitialize(), this insures that certain options
411: 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
412: calls to XXXSetFromOptions() it should not be used for options listed under PetscInitialize().
413: The collectivity of this routine is complex; only the MPI processes in comm will
414: have the affect of these options. If some processes that create objects call this routine and others do
415: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
416: on different ranks.
418: Level: developer
420: .seealso: PetscOptionsSetValue(), PetscOptionsView(), PetscOptionsHasName(), PetscOptionsGetInt(),
421: PetscOptionsGetReal(), PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
422: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
423: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
424: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
425: PetscOptionsFList(), PetscOptionsEList()
427: @*/
428: PetscErrorCode PetscOptionsInsertFile(MPI_Comm comm,PetscOptions options,const char file[],PetscBool require)
429: {
430: char *string,fname[PETSC_MAX_PATH_LEN],*vstring = NULL,*astring = NULL,*packed = NULL;
431: char *tokens[4];
433: size_t i,len,bytes;
434: FILE *fd;
435: PetscToken token=NULL;
436: int err;
437: char *cmatch;
438: const char cmt='#';
439: PetscInt line=1;
440: PetscMPIInt rank,cnt=0,acnt=0,counts[2];
441: PetscBool isdir,alias=PETSC_FALSE,valid;
444: MPI_Comm_rank(comm,&rank);CHKERRMPI(ierr);
445: PetscMemzero(tokens,sizeof(tokens));
446: if (!rank) {
447: cnt = 0;
448: acnt = 0;
450: PetscFixFilename(file,fname);
451: fd = fopen(fname,"r");
452: PetscTestDirectory(fname,'r',&isdir);
453: if (isdir && require) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Specified options file %s is a directory",fname);
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: PetscInfo1(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') {
470: string[i] = ' ';
471: }
472: }
473: PetscTokenCreate(string,' ',&token);
474: PetscTokenFind(token,&tokens[0]);
475: if (!tokens[0]) {
476: goto destroy;
477: } else if (!tokens[0][0]) { /* if token 0 is empty (string begins with spaces), redo */
478: PetscTokenFind(token,&tokens[0]);
479: }
480: for (i=1; i<4; i++) {
481: PetscTokenFind(token,&tokens[i]);
482: }
483: if (!tokens[0]) {
484: goto destroy;
485: } else if (tokens[0][0] == '-') {
486: PetscOptionsValidKey(tokens[0],&valid);
487: if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file %s line %D: invalid option %s",fname,line,tokens[0]);
488: PetscStrlen(tokens[0],&len);
489: PetscSegBufferGet(vseg,len+1,&vstring);
490: PetscArraycpy(vstring,tokens[0],len);
491: vstring[len] = ' ';
492: if (tokens[1]) {
493: PetscOptionsValidKey(tokens[1],&valid);
494: if (valid) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file %s line %D: cannot specify two options per line (%s %s)",fname,line,tokens[0],tokens[1]);
495: PetscStrlen(tokens[1],&len);
496: PetscSegBufferGet(vseg,len+3,&vstring);
497: vstring[0] = '"';
498: PetscArraycpy(vstring+1,tokens[1],len);
499: vstring[len+1] = '"';
500: vstring[len+2] = ' ';
501: }
502: } else {
503: PetscStrcasecmp(tokens[0],"alias",&alias);
504: if (alias) {
505: PetscOptionsValidKey(tokens[1],&valid);
506: if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file %s line %D: invalid aliased option %s",fname,line,tokens[1]);
507: if (!tokens[2]) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file %s line %D: alias missing for %s",fname,line,tokens[1]);
508: PetscOptionsValidKey(tokens[2],&valid);
509: if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file %s line %D: invalid aliasee option %s",fname,line,tokens[2]);
510: PetscStrlen(tokens[1],&len);
511: PetscSegBufferGet(aseg,len+1,&astring);
512: PetscArraycpy(astring,tokens[1],len);
513: astring[len] = ' ';
515: PetscStrlen(tokens[2],&len);
516: PetscSegBufferGet(aseg,len+1,&astring);
517: PetscArraycpy(astring,tokens[2],len);
518: astring[len] = ' ';
519: } else SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unknown first token in options file %s line %D: %s",fname,line,tokens[0]);
520: }
521: {
522: const char *extraToken = alias ? tokens[3] : tokens[2];
523: if (extraToken) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file %s line %D: extra token %s",fname,line,extraToken);
524: }
525: destroy:
526: free(string);
527: PetscTokenDestroy(&token);
528: alias = PETSC_FALSE;
529: line++;
530: }
531: err = fclose(fd);
532: if (err) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SYS,"fclose() failed on file %s",fname);
533: PetscSegBufferGetSize(aseg,&bytes); /* size without null termination */
534: PetscMPIIntCast(bytes,&acnt);
535: PetscSegBufferGet(aseg,1,&astring);
536: astring[0] = 0;
537: PetscSegBufferGetSize(vseg,&bytes); /* size without null termination */
538: PetscMPIIntCast(bytes,&cnt);
539: PetscSegBufferGet(vseg,1,&vstring);
540: vstring[0] = 0;
541: PetscMalloc1(2+acnt+cnt,&packed);
542: PetscSegBufferExtractTo(aseg,packed);
543: PetscSegBufferExtractTo(vseg,packed+acnt+1);
544: PetscSegBufferDestroy(&aseg);
545: PetscSegBufferDestroy(&vseg);
546: } else if (require) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Unable to open options file %s",fname);
547: }
549: counts[0] = acnt;
550: counts[1] = cnt;
551: err = MPI_Bcast(counts,2,MPI_INT,0,comm);
552: if (err) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in first MPI collective call, could be caused by using an incorrect mpiexec or a network problem, it can be caused by having VPN running: see https://www.mcs.anl.gov/petsc/documentation/faq.html");
553: acnt = counts[0];
554: cnt = counts[1];
555: if (rank) {
556: PetscMalloc1(2+acnt+cnt,&packed);
557: }
558: if (acnt || cnt) {
559: MPI_Bcast(packed,2+acnt+cnt,MPI_CHAR,0,comm);CHKERRMPI(ierr);
560: astring = packed;
561: vstring = packed + acnt + 1;
562: }
564: if (acnt) {
565: PetscTokenCreate(astring,' ',&token);
566: PetscTokenFind(token,&tokens[0]);
567: while (tokens[0]) {
568: PetscTokenFind(token,&tokens[1]);
569: PetscOptionsSetAlias(options,tokens[0],tokens[1]);
570: PetscTokenFind(token,&tokens[0]);
571: }
572: PetscTokenDestroy(&token);
573: }
575: if (cnt) {
576: PetscOptionsInsertString(options,vstring);
577: }
578: PetscFree(packed);
579: return(0);
580: }
582: static PetscErrorCode PetscOptionsInsertArgs(PetscOptions options,int argc,char *args[])
583: {
585: int left = argc - 1;
586: char **eargs = args + 1;
589: while (left) {
590: PetscBool isoptions_file,isprefixpush,isprefixpop,isp4,tisp4,isp4yourname,isp4rmrank,key;
591: PetscStrcasecmp(eargs[0],"-options_file",&isoptions_file);
592: PetscStrcasecmp(eargs[0],"-prefix_push",&isprefixpush);
593: PetscStrcasecmp(eargs[0],"-prefix_pop",&isprefixpop);
594: PetscStrcasecmp(eargs[0],"-p4pg",&isp4);
595: PetscStrcasecmp(eargs[0],"-p4yourname",&isp4yourname);
596: PetscStrcasecmp(eargs[0],"-p4rmrank",&isp4rmrank);
597: PetscStrcasecmp(eargs[0],"-p4wd",&tisp4);
598: isp4 = (PetscBool) (isp4 || tisp4);
599: PetscStrcasecmp(eargs[0],"-np",&tisp4);
600: isp4 = (PetscBool) (isp4 || tisp4);
601: PetscStrcasecmp(eargs[0],"-p4amslave",&tisp4);
602: PetscOptionsValidKey(eargs[0],&key);
604: if (!key) {
605: eargs++; left--;
606: } else if (isoptions_file) {
607: if (left <= 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing filename for -options_file filename option");
608: if (eargs[1][0] == '-') SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing filename for -options_file filename option");
609: PetscOptionsInsertFile(PETSC_COMM_WORLD,options,eargs[1],PETSC_TRUE);
610: eargs += 2; left -= 2;
611: } else if (isprefixpush) {
612: if (left <= 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing prefix for -prefix_push option");
613: if (eargs[1][0] == '-') SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing prefix for -prefix_push option (prefixes cannot start with '-')");
614: PetscOptionsPrefixPush(options,eargs[1]);
615: eargs += 2; left -= 2;
616: } else if (isprefixpop) {
617: PetscOptionsPrefixPop(options);
618: eargs++; left--;
620: /*
621: These are "bad" options that MPICH, etc put on the command line
622: we strip them out here.
623: */
624: } else if (tisp4 || isp4rmrank) {
625: eargs += 1; left -= 1;
626: } else if (isp4 || isp4yourname) {
627: eargs += 2; left -= 2;
628: } else {
629: PetscBool nextiskey = PETSC_FALSE;
630: if (left >= 2) {PetscOptionsValidKey(eargs[1],&nextiskey);}
631: if (left < 2 || nextiskey) {
632: PetscOptionsSetValue(options,eargs[0],NULL);
633: eargs++; left--;
634: } else {
635: PetscOptionsSetValue(options,eargs[0],eargs[1]);
636: eargs += 2; left -= 2;
637: }
638: }
639: }
640: return(0);
641: }
643: PETSC_STATIC_INLINE PetscErrorCode PetscOptionsStringToBoolIfSet_Private(enum PetscPrecedentOption opt,const char *val[],PetscBool set[],PetscBool *flg)
644: {
648: if (set[opt]) {
649: PetscOptionsStringToBool(val[opt],flg);
650: } else *flg = PETSC_FALSE;
651: return(0);
652: }
654: /* Process options with absolute precedence */
655: static PetscErrorCode PetscOptionsProcessPrecedentFlags(PetscOptions options,int argc,char *args[],PetscBool *skip_petscrc,PetscBool *skip_petscrc_set)
656: {
657: const char* const *opt = precedentOptions;
658: const size_t n = PO_NUM;
659: size_t o;
660: int a;
661: const char **val;
662: PetscBool *set;
663: PetscErrorCode ierr;
666: PetscCalloc2(n,&val,n,&set);
668: /* Look for options possibly set using PetscOptionsSetValue beforehand */
669: for (o=0; o<n; o++) {
670: PetscOptionsFindPair(options,NULL,opt[o],&val[o],&set[o]);
671: }
673: /* Loop through all args to collect last occuring value of each option */
674: for (a=1; a<argc; a++) {
675: PetscBool valid, eq;
677: PetscOptionsValidKey(args[a],&valid);
678: if (!valid) continue;
679: for (o=0; o<n; o++) {
680: PetscStrcasecmp(args[a],opt[o],&eq);
681: if (eq) {
682: set[o] = PETSC_TRUE;
683: if (a == argc-1 || !args[a+1] || !args[a+1][0] || args[a+1][0] == '-') val[o] = NULL;
684: else val[o] = args[a+1];
685: break;
686: }
687: }
688: }
690: /* Process flags */
691: PetscStrcasecmp(val[PO_HELP], "intro", &options->help_intro);
692: if (options->help_intro) options->help = PETSC_TRUE;
693: else {PetscOptionsStringToBoolIfSet_Private(PO_HELP, val,set,&options->help);}
694: PetscOptionsStringToBoolIfSet_Private(PO_OPTIONS_MONITOR_CANCEL,val,set,&options->monitorCancel);
695: PetscOptionsStringToBoolIfSet_Private(PO_OPTIONS_MONITOR, val,set,&options->monitorFromOptions);
696: PetscOptionsStringToBoolIfSet_Private(PO_SKIP_PETSCRC, val,set,skip_petscrc);
697: *skip_petscrc_set = set[PO_SKIP_PETSCRC];
699: /* Store precedent options in database and mark them as used */
700: for (o=0; o<n; o++) {
701: if (set[o]) {
702: int pos;
704: PetscOptionsSetValue_Private(options,opt[o],val[o],&pos);
705: options->used[pos] = PETSC_TRUE;
706: }
707: }
709: PetscFree2(val,set);
710: options->precedentProcessed = PETSC_TRUE;
711: return(0);
712: }
714: PETSC_STATIC_INLINE PetscErrorCode PetscOptionsSkipPrecedent(PetscOptions options,const char name[],PetscBool *flg)
715: {
716: int i;
719: *flg = PETSC_FALSE;
720: if (options->precedentProcessed) {
721: for (i=0; i<PO_NUM; i++) {
722: if (!PetscOptNameCmp(precedentOptions[i],name)) {
723: /* check if precedent option has been set already */
724: PetscOptionsFindPair(options,NULL,name,NULL,flg);
725: if (*flg) break;
726: }
727: }
728: }
729: return(0);
730: }
732: /*@C
733: PetscOptionsInsert - Inserts into the options database from the command line,
734: the environmental variable and a file.
736: Collective on PETSC_COMM_WORLD
738: Input Parameters:
739: + options - options database or NULL for the default global database
740: . argc - count of number of command line arguments
741: . args - the command line arguments
742: - file - [optional] PETSc database file, also checks ~/.petscrc, .petscrc and petscrc.
743: Use NULL to not check for code specific file.
744: Use -skip_petscrc in the code specific file (or command line) to skip ~/.petscrc, .petscrc and petscrc files.
746: Note:
747: Since PetscOptionsInsert() is automatically called by PetscInitialize(),
748: the user does not typically need to call this routine. PetscOptionsInsert()
749: can be called several times, adding additional entries into the database.
751: Options Database Keys:
752: . -options_file <filename> - read options from a file
754: See PetscInitialize() for options related to option database monitoring.
756: Level: advanced
758: .seealso: PetscOptionsDestroy(), PetscOptionsView(), PetscOptionsInsertString(), PetscOptionsInsertFile(),
759: PetscInitialize()
760: @*/
761: PetscErrorCode PetscOptionsInsert(PetscOptions options,int *argc,char ***args,const char file[])
762: {
764: PetscMPIInt rank;
765: char filename[PETSC_MAX_PATH_LEN];
766: PetscBool hasArgs = (argc && *argc) ? PETSC_TRUE : PETSC_FALSE;
767: PetscBool skipPetscrc = PETSC_FALSE, skipPetscrcSet = PETSC_FALSE;
771: if (hasArgs && !(args && *args)) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_NULL, "*argc > 1 but *args not given");
772: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr);
774: if (!options) {
775: PetscOptionsCreateDefault();
776: options = defaultoptions;
777: }
778: if (hasArgs) {
779: /* process options with absolute precedence */
780: PetscOptionsProcessPrecedentFlags(options,*argc,*args,&skipPetscrc,&skipPetscrcSet);
781: }
782: if (file && file[0]) {
783: PetscStrreplace(PETSC_COMM_WORLD,file,filename,PETSC_MAX_PATH_LEN);
784: PetscOptionsInsertFile(PETSC_COMM_WORLD,options,filename,PETSC_TRUE);
785: /* if -skip_petscrc has not been set from command line, check whether it has been set in the file */
786: if (!skipPetscrcSet) {PetscOptionsGetBool(options,NULL,"-skip_petscrc",&skipPetscrc,NULL);}
787: }
788: if (!skipPetscrc) {
789: PetscGetHomeDirectory(filename,PETSC_MAX_PATH_LEN-16);
790: /* PetscOptionsInsertFile() does a fopen() on rank0 only - so only rank0 HomeDir value is relavent */
791: if (filename[0]) { PetscStrcat(filename,"/.petscrc"); }
792: PetscOptionsInsertFile(PETSC_COMM_WORLD,options,filename,PETSC_FALSE);
793: PetscOptionsInsertFile(PETSC_COMM_WORLD,options,".petscrc",PETSC_FALSE);
794: PetscOptionsInsertFile(PETSC_COMM_WORLD,options,"petscrc",PETSC_FALSE);
795: }
797: /* insert environment options */
798: {
799: char *eoptions = NULL;
800: size_t len = 0;
801: if (!rank) {
802: eoptions = (char*)getenv("PETSC_OPTIONS");
803: PetscStrlen(eoptions,&len);
804: MPI_Bcast(&len,1,MPIU_SIZE_T,0,PETSC_COMM_WORLD);
805: } else {
806: MPI_Bcast(&len,1,MPIU_SIZE_T,0,PETSC_COMM_WORLD);CHKERRMPI(ierr);
807: if (len) {
808: PetscMalloc1(len+1,&eoptions);
809: }
810: }
811: if (len) {
812: MPI_Bcast(eoptions,len,MPI_CHAR,0,PETSC_COMM_WORLD);CHKERRMPI(ierr);
813: if (rank) eoptions[len] = 0;
814: PetscOptionsInsertString(options,eoptions);
815: if (rank) {PetscFree(eoptions);}
816: }
817: }
819: #if defined(PETSC_HAVE_YAML)
820: {
821: char *eoptions = NULL;
822: size_t len = 0;
823: if (!rank) {
824: eoptions = (char*)getenv("PETSC_OPTIONS_YAML");
825: PetscStrlen(eoptions,&len);
826: MPI_Bcast(&len,1,MPIU_SIZE_T,0,PETSC_COMM_WORLD);
827: } else {
828: MPI_Bcast(&len,1,MPIU_SIZE_T,0,PETSC_COMM_WORLD);CHKERRMPI(ierr);
829: if (len) {
830: PetscMalloc1(len+1,&eoptions);
831: }
832: }
833: if (len) {
834: MPI_Bcast(eoptions,len,MPI_CHAR,0,PETSC_COMM_WORLD);CHKERRMPI(ierr);
835: if (rank) eoptions[len] = 0;
836: PetscOptionsInsertStringYAML(options,eoptions);
837: if (rank) {PetscFree(eoptions);}
838: }
839: }
840: {
841: char yaml_file[PETSC_MAX_PATH_LEN];
842: char yaml_string[BUFSIZ];
843: PetscBool yaml_flg;
844: PetscOptionsGetString(NULL,NULL,"-options_file_yaml",yaml_file,sizeof(yaml_file),&yaml_flg);
845: if (yaml_flg) {
846: PetscOptionsInsertFileYAML(PETSC_COMM_WORLD,yaml_file,PETSC_TRUE);
847: }
848: PetscOptionsGetString(NULL,NULL,"-options_string_yaml",yaml_string,sizeof(yaml_string),&yaml_flg);
849: if (yaml_flg) {
850: PetscOptionsInsertStringYAML(NULL,yaml_string);
851: }
852: }
853: #endif
855: /* insert command line options here because they take precedence over arguments in petscrc/environment */
856: if (hasArgs) {PetscOptionsInsertArgs(options,*argc,*args);}
857: return(0);
858: }
860: /*@C
861: PetscOptionsView - Prints the options that have been loaded. This is
862: useful for debugging purposes.
864: Logically Collective on PetscViewer
866: Input Parameter:
867: + options - options database, use NULL for default global database
868: - viewer - must be an PETSCVIEWERASCII viewer
870: Options Database Key:
871: . -options_view - Activates PetscOptionsView() within PetscFinalize()
873: Notes:
874: Only the rank zero process of MPI_Comm used to create view prints the option values. Other processes
875: may have different values but they are not printed.
877: Level: advanced
879: .seealso: PetscOptionsAllUsed()
880: @*/
881: PetscErrorCode PetscOptionsView(PetscOptions options,PetscViewer viewer)
882: {
884: PetscInt i;
885: PetscBool isascii;
889: options = options ? options : defaultoptions;
890: if (!viewer) viewer = PETSC_VIEWER_STDOUT_WORLD;
891: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
892: if (!isascii) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Only supports ASCII viewer");
894: if (!options->N) {
895: PetscViewerASCIIPrintf(viewer,"#No PETSc Option Table entries\n");
896: return(0);
897: }
899: PetscViewerASCIIPrintf(viewer,"#PETSc Option Table entries:\n");
900: for (i=0; i<options->N; i++) {
901: if (options->values[i]) {
902: PetscViewerASCIIPrintf(viewer,"-%s %s\n",options->names[i],options->values[i]);
903: } else {
904: PetscViewerASCIIPrintf(viewer,"-%s\n",options->names[i]);
905: }
906: }
907: PetscViewerASCIIPrintf(viewer,"#End of PETSc Option Table entries\n");
908: return(0);
909: }
911: /*
912: Called by error handlers to print options used in run
913: */
914: PETSC_EXTERN PetscErrorCode PetscOptionsViewError(void)
915: {
916: PetscInt i;
917: PetscOptions options = defaultoptions;
920: if (options->N) {
921: (*PetscErrorPrintf)("PETSc Option Table entries:\n");
922: } else {
923: (*PetscErrorPrintf)("No PETSc Option Table entries\n");
924: }
925: for (i=0; i<options->N; i++) {
926: if (options->values[i]) {
927: (*PetscErrorPrintf)("-%s %s\n",options->names[i],options->values[i]);
928: } else {
929: (*PetscErrorPrintf)("-%s\n",options->names[i]);
930: }
931: }
932: return(0);
933: }
935: /*@C
936: PetscOptionsPrefixPush - Designate a prefix to be used by all options insertions to follow.
938: Logically Collective
940: Input Parameter:
941: + options - options database, or NULL for the default global database
942: - prefix - The string to append to the existing prefix
944: Options Database Keys:
945: + -prefix_push <some_prefix_> - push the given prefix
946: - -prefix_pop - pop the last prefix
948: Notes:
949: It is common to use this in conjunction with -options_file as in
951: $ -prefix_push system1_ -options_file system1rc -prefix_pop -prefix_push system2_ -options_file system2rc -prefix_pop
953: where the files no longer require all options to be prefixed with -system2_.
955: The collectivity of this routine is complex; only the MPI processes that call this routine will
956: have the affect of these options. If some processes that create objects call this routine and others do
957: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
958: on different ranks.
960: Level: advanced
962: .seealso: PetscOptionsPrefixPop(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsSetValue()
963: @*/
964: PetscErrorCode PetscOptionsPrefixPush(PetscOptions options,const char prefix[])
965: {
967: size_t n;
968: PetscInt start;
969: char key[MAXOPTNAME+1];
970: PetscBool valid;
974: options = options ? options : defaultoptions;
975: if (options->prefixind >= MAXPREFIXES) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Maximum depth of prefix stack %d exceeded, recompile \n src/sys/objects/options.c with larger value for MAXPREFIXES",MAXPREFIXES);
976: key[0] = '-'; /* keys must start with '-' */
977: PetscStrncpy(key+1,prefix,sizeof(key)-1);
978: PetscOptionsValidKey(key,&valid);
979: if (!valid) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Given prefix \"%s\" not valid (the first character must be a letter, do not include leading '-')",prefix);
980: start = options->prefixind ? options->prefixstack[options->prefixind-1] : 0;
981: PetscStrlen(prefix,&n);
982: if (n+1 > sizeof(options->prefix)-start) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Maximum prefix length %d exceeded",sizeof(options->prefix));
983: PetscArraycpy(options->prefix+start,prefix,n+1);
984: options->prefixstack[options->prefixind++] = start+n;
985: return(0);
986: }
988: /*@C
989: PetscOptionsPrefixPop - Remove the latest options prefix, see PetscOptionsPrefixPush() for details
991: Logically Collective on the MPI_Comm that called PetscOptionsPrefixPush()
993: Input Parameters:
994: . options - options database, or NULL for the default global database
996: Level: advanced
998: .seealso: PetscOptionsPrefixPush(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsSetValue()
999: @*/
1000: PetscErrorCode PetscOptionsPrefixPop(PetscOptions options)
1001: {
1002: PetscInt offset;
1005: options = options ? options : defaultoptions;
1006: if (options->prefixind < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"More prefixes popped than pushed");
1007: options->prefixind--;
1008: offset = options->prefixind ? options->prefixstack[options->prefixind-1] : 0;
1009: options->prefix[offset] = 0;
1010: return(0);
1011: }
1013: /*@C
1014: PetscOptionsClear - Removes all options form the database leaving it empty.
1016: Logically Collective
1018: Input Parameters:
1019: . options - options database, use NULL for the default global database
1021: The collectivity of this routine is complex; only the MPI processes that call this routine will
1022: have the affect of these options. If some processes that create objects call this routine and others do
1023: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1024: on different ranks.
1026: Level: developer
1028: .seealso: PetscOptionsInsert()
1029: @*/
1030: PetscErrorCode PetscOptionsClear(PetscOptions options)
1031: {
1032: PetscInt i;
1034: options = options ? options : defaultoptions;
1035: if (!options) return 0;
1037: for (i=0; i<options->N; i++) {
1038: if (options->names[i]) free(options->names[i]);
1039: if (options->values[i]) free(options->values[i]);
1040: }
1041: options->N = 0;
1043: for (i=0; i<options->Naliases; i++) {
1044: free(options->aliases1[i]);
1045: free(options->aliases2[i]);
1046: }
1047: options->Naliases = 0;
1049: /* destroy hash table */
1050: kh_destroy(HO,options->ht);
1051: options->ht = NULL;
1053: options->prefixind = 0;
1054: options->prefix[0] = 0;
1055: options->help = PETSC_FALSE;
1056: return 0;
1057: }
1059: /*@C
1060: PetscOptionsSetAlias - Makes a key and alias for another key
1062: Logically Collective
1064: Input Parameters:
1065: + options - options database, or NULL for default global database
1066: . newname - the alias
1067: - oldname - the name that alias will refer to
1069: Level: advanced
1071: The collectivity of this routine is complex; only the MPI processes that call this routine will
1072: have the affect of these options. If some processes that create objects call this routine and others do
1073: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1074: on different ranks.
1076: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),OptionsHasName(),
1077: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(),PetscOptionsBool(),
1078: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
1079: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
1080: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
1081: PetscOptionsFList(), PetscOptionsEList()
1082: @*/
1083: PetscErrorCode PetscOptionsSetAlias(PetscOptions options,const char newname[],const char oldname[])
1084: {
1085: PetscInt n;
1086: size_t len;
1087: PetscBool valid;
1093: options = options ? options : defaultoptions;
1094: PetscOptionsValidKey(newname,&valid);
1095: if (!valid) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid aliased option %s",newname);
1096: PetscOptionsValidKey(oldname,&valid);
1097: if (!valid) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid aliasee option %s",oldname);
1099: n = options->Naliases;
1100: if (n >= MAXALIASES) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_MEM,"You have defined to many PETSc options aliases, limit %d recompile \n src/sys/objects/options.c with larger value for MAXALIASES",MAXALIASES);
1102: newname++; oldname++;
1103: PetscStrlen(newname,&len);
1104: options->aliases1[n] = (char*)malloc((len+1)*sizeof(char));
1105: PetscStrcpy(options->aliases1[n],newname);
1106: PetscStrlen(oldname,&len);
1107: options->aliases2[n] = (char*)malloc((len+1)*sizeof(char));
1108: PetscStrcpy(options->aliases2[n],oldname);
1109: options->Naliases++;
1110: return(0);
1111: }
1113: /*@C
1114: PetscOptionsSetValue - Sets an option name-value pair in the options
1115: database, overriding whatever is already present.
1117: Logically Collective
1119: Input Parameters:
1120: + options - options database, use NULL for the default global database
1121: . name - name of option, this SHOULD have the - prepended
1122: - value - the option value (not used for all options, so can be NULL)
1124: Level: intermediate
1126: Note:
1127: This function can be called BEFORE PetscInitialize()
1129: The collectivity of this routine is complex; only the MPI processes that call this routine will
1130: have the affect of these options. If some processes that create objects call this routine and others do
1131: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1132: on different ranks.
1134: Developers Note: Uses malloc() directly because PETSc may not be initialized yet.
1136: .seealso: PetscOptionsInsert(), PetscOptionsClearValue()
1137: @*/
1138: PetscErrorCode PetscOptionsSetValue(PetscOptions options,const char name[],const char value[])
1139: {
1140: return PetscOptionsSetValue_Private(options,name,value,NULL);
1141: }
1143: static PetscErrorCode PetscOptionsSetValue_Private(PetscOptions options,const char name[],const char value[],int *pos)
1144: {
1145: size_t len;
1146: int N,n,i;
1147: char **names;
1148: char fullname[MAXOPTNAME] = "";
1149: PetscBool flg;
1152: if (!options) {
1153: PetscOptionsCreateDefault();if (ierr) return ierr;
1154: options = defaultoptions;
1155: }
1157: if (name[0] != '-') return PETSC_ERR_ARG_OUTOFRANGE;
1159: PetscOptionsSkipPrecedent(options,name,&flg);
1160: if (flg) return 0;
1162: name++; /* skip starting dash */
1164: if (options->prefixind > 0) {
1165: strncpy(fullname,options->prefix,sizeof(fullname));
1166: fullname[sizeof(fullname)-1] = 0;
1167: strncat(fullname,name,sizeof(fullname)-strlen(fullname)-1);
1168: fullname[sizeof(fullname)-1] = 0;
1169: name = fullname;
1170: }
1172: /* check against aliases */
1173: N = options->Naliases;
1174: for (i=0; i<N; i++) {
1175: int result = PetscOptNameCmp(options->aliases1[i],name);
1176: if (!result) { name = options->aliases2[i]; break; }
1177: }
1179: /* slow search */
1180: N = n = options->N;
1181: names = options->names;
1182: for (i=0; i<N; i++) {
1183: int result = PetscOptNameCmp(names[i],name);
1184: if (!result) {
1185: n = i; goto setvalue;
1186: } else if (result > 0) {
1187: n = i; break;
1188: }
1189: }
1190: if (N >= MAXOPTIONS) return PETSC_ERR_MEM;
1191: /* shift remaining values up 1 */
1192: for (i=N; i>n; i--) {
1193: options->names[i] = options->names[i-1];
1194: options->values[i] = options->values[i-1];
1195: options->used[i] = options->used[i-1];
1196: }
1197: options->names[n] = NULL;
1198: options->values[n] = NULL;
1199: options->used[n] = PETSC_FALSE;
1200: options->N++;
1202: /* destroy hash table */
1203: kh_destroy(HO,options->ht);
1204: options->ht = NULL;
1206: /* set new name */
1207: len = strlen(name);
1208: options->names[n] = (char*)malloc((len+1)*sizeof(char));
1209: if (!options->names[n]) return PETSC_ERR_MEM;
1210: strcpy(options->names[n],name);
1212: setvalue:
1213: /* set new value */
1214: if (options->values[n]) free(options->values[n]);
1215: len = value ? strlen(value) : 0;
1216: if (len) {
1217: options->values[n] = (char*)malloc((len+1)*sizeof(char));
1218: if (!options->values[n]) return PETSC_ERR_MEM;
1219: strcpy(options->values[n],value);
1220: } else {
1221: options->values[n] = NULL;
1222: }
1224: /* handle -help so that it can be set from anywhere */
1225: if (!PetscOptNameCmp(name,"help")) {
1226: options->help = PETSC_TRUE;
1227: PetscStrcasecmp(value, "intro", &options->help_intro);
1228: options->used[n] = PETSC_TRUE;
1229: }
1231: if (PetscErrorHandlingInitialized) {
1232: PetscOptionsMonitor(options,name,value);
1233: }
1234: if (pos) *pos = n;
1235: return 0;
1236: }
1238: /*@C
1239: PetscOptionsClearValue - Clears an option name-value pair in the options
1240: database, overriding whatever is already present.
1242: Logically Collective
1244: Input Parameter:
1245: + options - options database, use NULL for the default global database
1246: - name - name of option, this SHOULD have the - prepended
1248: Level: intermediate
1250: The collectivity of this routine is complex; only the MPI processes that call this routine will
1251: have the affect of these options. If some processes that create objects call this routine and others do
1252: not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1253: on different ranks.
1255: .seealso: PetscOptionsInsert()
1256: @*/
1257: PetscErrorCode PetscOptionsClearValue(PetscOptions options,const char name[])
1258: {
1259: int N,n,i;
1260: char **names;
1264: options = options ? options : defaultoptions;
1265: if (name[0] != '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Name must begin with '-': Instead %s",name);
1267: if (!PetscOptNameCmp(name,"-help")) options->help = PETSC_FALSE;
1269: name++; /* skip starting dash */
1271: /* slow search */
1272: N = n = options->N;
1273: names = options->names;
1274: for (i=0; i<N; i++) {
1275: int result = PetscOptNameCmp(names[i],name);
1276: if (!result) {
1277: n = i; break;
1278: } else if (result > 0) {
1279: n = N; break;
1280: }
1281: }
1282: if (n == N) return(0); /* it was not present */
1284: /* remove name and value */
1285: if (options->names[n]) free(options->names[n]);
1286: if (options->values[n]) free(options->values[n]);
1287: /* shift remaining values down 1 */
1288: for (i=n; i<N-1; i++) {
1289: options->names[i] = options->names[i+1];
1290: options->values[i] = options->values[i+1];
1291: options->used[i] = options->used[i+1];
1292: }
1293: options->N--;
1295: /* destroy hash table */
1296: kh_destroy(HO,options->ht);
1297: options->ht = NULL;
1299: PetscOptionsMonitor(options,name,NULL);
1300: return(0);
1301: }
1303: /*@C
1304: PetscOptionsFindPair - Gets an option name-value pair from the options database.
1306: Not Collective
1308: Input Parameters:
1309: + options - options database, use NULL for the default global database
1310: . pre - the string to prepend to the name or NULL, this SHOULD NOT have the "-" prepended
1311: - name - name of option, this SHOULD have the "-" prepended
1313: Output Parameters:
1314: + value - the option value (optional, not used for all options)
1315: - set - whether the option is set (optional)
1317: Notes:
1318: Each process may find different values or no value depending on how options were inserted into the database
1320: Level: developer
1322: .seealso: PetscOptionsSetValue(), PetscOptionsClearValue()
1323: @*/
1324: PetscErrorCode PetscOptionsFindPair(PetscOptions options,const char pre[],const char name[],const char *value[],PetscBool *set)
1325: {
1326: char buf[MAXOPTNAME];
1327: PetscBool usehashtable = PETSC_TRUE;
1328: PetscBool matchnumbers = PETSC_TRUE;
1332: options = options ? options : defaultoptions;
1333: if (pre && PetscUnlikely(pre[0] == '-')) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Prefix cannot begin with '-': Instead %s",pre);
1334: if (PetscUnlikely(name[0] != '-')) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Name must begin with '-': Instead %s",name);
1336: name++; /* skip starting dash */
1338: /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */
1339: if (pre && pre[0]) {
1340: char *ptr = buf;
1341: if (name[0] == '-') { *ptr++ = '-'; name++; }
1342: PetscStrncpy(ptr,pre,buf+sizeof(buf)-ptr);
1343: PetscStrlcat(buf,name,sizeof(buf));
1344: name = buf;
1345: }
1347: if (PetscDefined(USE_DEBUG)) {
1348: PetscBool valid;
1349: char key[MAXOPTNAME+1] = "-";
1350: PetscStrncpy(key+1,name,sizeof(key)-1);
1351: PetscOptionsValidKey(key,&valid);
1352: if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid option '%s' obtained from pre='%s' and name='%s'",key,pre?pre:"",name);
1353: }
1355: if (!options->ht && usehashtable) {
1356: int i,ret;
1357: khiter_t it;
1358: khash_t(HO) *ht;
1359: ht = kh_init(HO);
1360: if (PetscUnlikely(!ht)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_MEM,"Hash table allocation failed");
1361: ret = kh_resize(HO,ht,options->N*2); /* twice the required size to reduce risk of collisions */
1362: if (PetscUnlikely(ret)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_MEM,"Hash table allocation failed");
1363: for (i=0; i<options->N; i++) {
1364: it = kh_put(HO,ht,options->names[i],&ret);
1365: if (PetscUnlikely(ret != 1)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_MEM,"Hash table allocation failed");
1366: kh_val(ht,it) = i;
1367: }
1368: options->ht = ht;
1369: }
1371: if (usehashtable)
1372: { /* fast search */
1373: khash_t(HO) *ht = options->ht;
1374: khiter_t it = kh_get(HO,ht,name);
1375: if (it != kh_end(ht)) {
1376: int i = kh_val(ht,it);
1377: options->used[i] = PETSC_TRUE;
1378: if (value) *value = options->values[i];
1379: if (set) *set = PETSC_TRUE;
1380: return(0);
1381: }
1382: } else
1383: { /* slow search */
1384: int i, N = options->N;
1385: for (i=0; i<N; i++) {
1386: int result = PetscOptNameCmp(options->names[i],name);
1387: if (!result) {
1388: options->used[i] = PETSC_TRUE;
1389: if (value) *value = options->values[i];
1390: if (set) *set = PETSC_TRUE;
1391: return(0);
1392: } else if (result > 0) {
1393: break;
1394: }
1395: }
1396: }
1398: /*
1399: The following block slows down all lookups in the most frequent path (most lookups are unsuccessful).
1400: Maybe this special lookup mode should be enabled on request with a push/pop API.
1401: The feature of matching _%d_ used sparingly in the codebase.
1402: */
1403: if (matchnumbers) {
1404: int i,j,cnt = 0,locs[16],loce[16];
1405: /* determine the location and number of all _%d_ in the key */
1406: for (i=0; name[i]; i++) {
1407: if (name[i] == '_') {
1408: for (j=i+1; name[j]; j++) {
1409: if (name[j] >= '0' && name[j] <= '9') continue;
1410: if (name[j] == '_' && j > i+1) { /* found a number */
1411: locs[cnt] = i+1;
1412: loce[cnt++] = j+1;
1413: }
1414: i = j-1;
1415: break;
1416: }
1417: }
1418: }
1419: for (i=0; i<cnt; i++) {
1420: PetscBool found;
1421: char opt[MAXOPTNAME+1] = "-", tmp[MAXOPTNAME];
1422: PetscStrncpy(tmp,name,PetscMin((size_t)(locs[i]+1),sizeof(tmp)));
1423: PetscStrlcat(opt,tmp,sizeof(opt));
1424: PetscStrlcat(opt,name+loce[i],sizeof(opt));
1425: PetscOptionsFindPair(options,NULL,opt,value,&found);
1426: if (found) {if (set) *set = PETSC_TRUE; return(0);}
1427: }
1428: }
1430: if (set) *set = PETSC_FALSE;
1431: return(0);
1432: }
1434: /* Check whether any option begins with pre+name */
1435: PETSC_EXTERN PetscErrorCode PetscOptionsFindPairPrefix_Private(PetscOptions options,const char pre[], const char name[],const char *value[],PetscBool *set)
1436: {
1437: char buf[MAXOPTNAME];
1438: int numCnt = 0, locs[16],loce[16];
1442: options = options ? options : defaultoptions;
1443: if (pre && pre[0] == '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Prefix cannot begin with '-': Instead %s",pre);
1444: if (name[0] != '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Name must begin with '-': Instead %s",name);
1446: name++; /* skip starting dash */
1448: /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */
1449: if (pre && pre[0]) {
1450: char *ptr = buf;
1451: if (name[0] == '-') { *ptr++ = '-'; name++; }
1452: PetscStrncpy(ptr,pre,sizeof(buf)+(size_t)(ptr-buf));
1453: PetscStrlcat(buf,name,sizeof(buf));
1454: name = buf;
1455: }
1457: if (PetscDefined(USE_DEBUG)) {
1458: PetscBool valid;
1459: char key[MAXOPTNAME+1] = "-";
1460: PetscStrncpy(key+1,name,sizeof(key)-1);
1461: PetscOptionsValidKey(key,&valid);
1462: if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid option '%s' obtained from pre='%s' and name='%s'",key,pre?pre:"",name);
1463: }
1465: /* determine the location and number of all _%d_ in the key */
1466: {
1467: int i,j;
1468: for (i=0; name[i]; i++) {
1469: if (name[i] == '_') {
1470: for (j=i+1; name[j]; j++) {
1471: if (name[j] >= '0' && name[j] <= '9') continue;
1472: if (name[j] == '_' && j > i+1) { /* found a number */
1473: locs[numCnt] = i+1;
1474: loce[numCnt++] = j+1;
1475: }
1476: i = j-1;
1477: break;
1478: }
1479: }
1480: }
1481: }
1483: { /* slow search */
1484: int c, i;
1485: size_t len;
1486: PetscBool match;
1488: for (c = -1; c < numCnt; ++c) {
1489: char opt[MAXOPTNAME+1] = "", tmp[MAXOPTNAME];
1491: if (c < 0) {
1492: PetscStrcpy(opt,name);
1493: } else {
1494: PetscStrncpy(tmp,name,PetscMin((size_t)(locs[c]+1),sizeof(tmp)));
1495: PetscStrlcat(opt,tmp,sizeof(opt));
1496: PetscStrlcat(opt,name+loce[c],sizeof(opt));
1497: }
1498: PetscStrlen(opt,&len);
1499: for (i=0; i<options->N; i++) {
1500: PetscStrncmp(options->names[i],opt,len,&match);
1501: if (match) {
1502: options->used[i] = PETSC_TRUE;
1503: if (value) *value = options->values[i];
1504: if (set) *set = PETSC_TRUE;
1505: return(0);
1506: }
1507: }
1508: }
1509: }
1511: if (set) *set = PETSC_FALSE;
1512: return(0);
1513: }
1515: /*@C
1516: PetscOptionsReject - Generates an error if a certain option is given.
1518: Not Collective
1520: Input Parameters:
1521: + options - options database, use NULL for default global database
1522: . pre - the option prefix (may be NULL)
1523: . name - the option name one is seeking
1524: - mess - error message (may be NULL)
1526: Level: advanced
1528: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),OptionsHasName(),
1529: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
1530: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
1531: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
1532: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
1533: PetscOptionsFList(), PetscOptionsEList()
1534: @*/
1535: PetscErrorCode PetscOptionsReject(PetscOptions options,const char pre[],const char name[],const char mess[])
1536: {
1538: PetscBool flag = PETSC_FALSE;
1541: PetscOptionsHasName(options,pre,name,&flag);
1542: if (flag) {
1543: if (mess && mess[0]) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Program has disabled option: -%s%s with %s",pre?pre:"",name+1,mess);
1544: else SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Program has disabled option: -%s%s",pre?pre:"",name+1);
1545: }
1546: return(0);
1547: }
1549: /*@C
1550: PetscOptionsHasHelp - Determines whether the "-help" option is in the database.
1552: Not Collective
1554: Input Parameters:
1555: . options - options database, use NULL for default global database
1557: Output Parameters:
1558: . set - PETSC_TRUE if found else PETSC_FALSE.
1560: Level: advanced
1562: .seealso: PetscOptionsHasName()
1563: @*/
1564: PetscErrorCode PetscOptionsHasHelp(PetscOptions options,PetscBool *set)
1565: {
1568: options = options ? options : defaultoptions;
1569: *set = options->help;
1570: return(0);
1571: }
1573: PetscErrorCode PetscOptionsHasHelpIntro_Internal(PetscOptions options,PetscBool *set)
1574: {
1577: options = options ? options : defaultoptions;
1578: *set = options->help_intro;
1579: return(0);
1580: }
1582: /*@C
1583: PetscOptionsHasName - Determines whether a certain option is given in the database. This returns true whether the option is a number, string or boolean, even
1584: its value is set to false.
1586: Not Collective
1588: Input Parameters:
1589: + options - options database, use NULL for default global database
1590: . pre - string to prepend to the name or NULL
1591: - name - the option one is seeking
1593: Output Parameters:
1594: . set - PETSC_TRUE if found else PETSC_FALSE.
1596: Level: beginner
1598: Notes:
1599: In many cases you probably want to use PetscOptionsGetBool() instead of calling this, to allowing toggling values.
1601: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
1602: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
1603: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
1604: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
1605: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
1606: PetscOptionsFList(), PetscOptionsEList()
1607: @*/
1608: PetscErrorCode PetscOptionsHasName(PetscOptions options,const char pre[],const char name[],PetscBool *set)
1609: {
1610: const char *value;
1612: PetscBool flag;
1615: PetscOptionsFindPair(options,pre,name,&value,&flag);
1616: if (set) *set = flag;
1617: return(0);
1618: }
1620: /*@C
1621: PetscOptionsGetAll - Lists all the options the program was run with in a single string.
1623: Not Collective
1625: Input Parameter:
1626: . options - the options database, use NULL for the default global database
1628: Output Parameter:
1629: . copts - pointer where string pointer is stored
1631: Notes:
1632: The array and each entry in the array should be freed with PetscFree()
1633: Each process may have different values depending on how the options were inserted into the database
1635: Level: advanced
1637: .seealso: PetscOptionsAllUsed(), PetscOptionsView(), PetscOptionsPush(), PetscOptionsPop()
1638: @*/
1639: PetscErrorCode PetscOptionsGetAll(PetscOptions options,char *copts[])
1640: {
1642: PetscInt i;
1643: size_t len = 1,lent = 0;
1644: char *coptions = NULL;
1648: options = options ? options : defaultoptions;
1649: /* count the length of the required string */
1650: for (i=0; i<options->N; i++) {
1651: PetscStrlen(options->names[i],&lent);
1652: len += 2 + lent;
1653: if (options->values[i]) {
1654: PetscStrlen(options->values[i],&lent);
1655: len += 1 + lent;
1656: }
1657: }
1658: PetscMalloc1(len,&coptions);
1659: coptions[0] = 0;
1660: for (i=0; i<options->N; i++) {
1661: PetscStrcat(coptions,"-");
1662: PetscStrcat(coptions,options->names[i]);
1663: PetscStrcat(coptions," ");
1664: if (options->values[i]) {
1665: PetscStrcat(coptions,options->values[i]);
1666: PetscStrcat(coptions," ");
1667: }
1668: }
1669: *copts = coptions;
1670: return(0);
1671: }
1673: /*@C
1674: PetscOptionsUsed - Indicates if PETSc has used a particular option set in the database
1676: Not Collective
1678: Input Parameter:
1679: + options - options database, use NULL for default global database
1680: - name - string name of option
1682: Output Parameter:
1683: . used - PETSC_TRUE if the option was used, otherwise false, including if option was not found in options database
1685: Level: advanced
1687: Notes:
1688: The value returned may be different on each process and depends on which options have been processed
1689: on the given process
1691: .seealso: PetscOptionsView(), PetscOptionsLeft(), PetscOptionsAllUsed()
1692: @*/
1693: PetscErrorCode PetscOptionsUsed(PetscOptions options,const char *name,PetscBool *used)
1694: {
1695: PetscInt i;
1701: options = options ? options : defaultoptions;
1702: *used = PETSC_FALSE;
1703: for (i=0; i<options->N; i++) {
1704: PetscStrcmp(options->names[i],name,used);
1705: if (*used) {
1706: *used = options->used[i];
1707: break;
1708: }
1709: }
1710: return(0);
1711: }
1713: /*@
1714: PetscOptionsAllUsed - Returns a count of the number of options in the
1715: database that have never been selected.
1717: Not Collective
1719: Input Parameter:
1720: . options - options database, use NULL for default global database
1722: Output Parameter:
1723: . N - count of options not used
1725: Level: advanced
1727: Notes:
1728: The value returned may be different on each process and depends on which options have been processed
1729: on the given process
1731: .seealso: PetscOptionsView()
1732: @*/
1733: PetscErrorCode PetscOptionsAllUsed(PetscOptions options,PetscInt *N)
1734: {
1735: PetscInt i,n = 0;
1739: options = options ? options : defaultoptions;
1740: for (i=0; i<options->N; i++) {
1741: if (!options->used[i]) n++;
1742: }
1743: *N = n;
1744: return(0);
1745: }
1747: /*@
1748: PetscOptionsLeft - Prints to screen any options that were set and never used.
1750: Not Collective
1752: Input Parameter:
1753: . options - options database; use NULL for default global database
1755: Options Database Key:
1756: . -options_left - activates PetscOptionsAllUsed() within PetscFinalize()
1758: Notes:
1759: This is rarely used directly, it is called by PetscFinalize() in debug more or if -options_left
1760: is passed otherwise to help users determine possible mistakes in their usage of options. This
1761: only prints values on process zero of PETSC_COMM_WORLD. Other processes depending the objects
1762: used may have different options that are left unused.
1764: Level: advanced
1766: .seealso: PetscOptionsAllUsed()
1767: @*/
1768: PetscErrorCode PetscOptionsLeft(PetscOptions options)
1769: {
1771: PetscInt i;
1772: PetscInt cnt = 0;
1773: PetscOptions toptions;
1776: toptions = options ? options : defaultoptions;
1777: for (i=0; i<toptions->N; i++) {
1778: if (!toptions->used[i]) {
1779: if (toptions->values[i]) {
1780: PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s value: %s\n",toptions->names[i],toptions->values[i]);
1781: } else {
1782: PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s (no value)\n",toptions->names[i]);
1783: }
1784: }
1785: }
1786: if (!options) {
1787: toptions = defaultoptions;
1788: while (toptions->previous) {
1789: cnt++;
1790: toptions = toptions->previous;
1791: }
1792: if (cnt) {
1793: PetscPrintf(PETSC_COMM_WORLD,"Option left: You may have forgotten some calls to PetscOptionsPop(),\n PetscOptionsPop() has been called %D less times than PetscOptionsPush()\n",cnt);
1794: }
1795: }
1796: return(0);
1797: }
1799: /*@C
1800: PetscOptionsLeftGet - Returns all options that were set and never used.
1802: Not Collective
1804: Input Parameter:
1805: . options - options database, use NULL for default global database
1807: Output Parameter:
1808: + N - count of options not used
1809: . names - names of options not used
1810: - values - values of options not used
1812: Level: advanced
1814: Notes:
1815: Users should call PetscOptionsLeftRestore() to free the memory allocated in this routine
1816: Notes: The value returned may be different on each process and depends on which options have been processed
1817: on the given process
1819: .seealso: PetscOptionsAllUsed(), PetscOptionsLeft()
1820: @*/
1821: PetscErrorCode PetscOptionsLeftGet(PetscOptions options,PetscInt *N,char **names[],char **values[])
1822: {
1824: PetscInt i,n;
1830: options = options ? options : defaultoptions;
1832: /* The number of unused PETSc options */
1833: n = 0;
1834: for (i=0; i<options->N; i++) {
1835: if (!options->used[i]) n++;
1836: }
1837: if (N) { *N = n; }
1838: if (names) { PetscMalloc1(n,names); }
1839: if (values) { PetscMalloc1(n,values); }
1841: n = 0;
1842: if (names || values) {
1843: for (i=0; i<options->N; i++) {
1844: if (!options->used[i]) {
1845: if (names) (*names)[n] = options->names[i];
1846: if (values) (*values)[n] = options->values[i];
1847: n++;
1848: }
1849: }
1850: }
1851: return(0);
1852: }
1854: /*@C
1855: PetscOptionsLeftRestore - Free memory for the unused PETSc options obtained using PetscOptionsLeftGet.
1857: Not Collective
1859: Input Parameter:
1860: + options - options database, use NULL for default global database
1861: . names - names of options not used
1862: - values - values of options not used
1864: Level: advanced
1866: .seealso: PetscOptionsAllUsed(), PetscOptionsLeft(), PetscOptionsLeftGet()
1867: @*/
1868: PetscErrorCode PetscOptionsLeftRestore(PetscOptions options,PetscInt *N,char **names[],char **values[])
1869: {
1876: if (N) { *N = 0; }
1877: if (names) { PetscFree(*names); }
1878: if (values) { PetscFree(*values); }
1879: return(0);
1880: }
1882: /*@C
1883: PetscOptionsMonitorDefault - Print all options set value events using the supplied PetscViewer.
1885: Logically Collective on ctx
1887: Input Parameters:
1888: + name - option name string
1889: . value - option value string
1890: - ctx - an ASCII viewer or NULL
1892: Level: intermediate
1894: Notes:
1895: If ctx=NULL, PetscPrintf() is used.
1896: The first MPI rank in the PetscViewer viewer actually prints the values, other
1897: processes may have different values set
1899: .seealso: PetscOptionsMonitorSet()
1900: @*/
1901: PetscErrorCode PetscOptionsMonitorDefault(const char name[],const char value[],void *ctx)
1902: {
1906: if (ctx) {
1907: PetscViewer viewer = (PetscViewer)ctx;
1908: if (!value) {
1909: PetscViewerASCIIPrintf(viewer,"Removing option: %s\n",name,value);
1910: } else if (!value[0]) {
1911: PetscViewerASCIIPrintf(viewer,"Setting option: %s (no value)\n",name);
1912: } else {
1913: PetscViewerASCIIPrintf(viewer,"Setting option: %s = %s\n",name,value);
1914: }
1915: } else {
1916: MPI_Comm comm = PETSC_COMM_WORLD;
1917: if (!value) {
1918: PetscPrintf(comm,"Removing option: %s\n",name,value);
1919: } else if (!value[0]) {
1920: PetscPrintf(comm,"Setting option: %s (no value)\n",name);
1921: } else {
1922: PetscPrintf(comm,"Setting option: %s = %s\n",name,value);
1923: }
1924: }
1925: return(0);
1926: }
1928: /*@C
1929: PetscOptionsMonitorSet - Sets an ADDITIONAL function to be called at every method that
1930: modified the PETSc options database.
1932: Not Collective
1934: Input Parameters:
1935: + monitor - pointer to function (if this is NULL, it turns off monitoring
1936: . mctx - [optional] context for private data for the
1937: monitor routine (use NULL if no context is desired)
1938: - monitordestroy - [optional] routine that frees monitor context
1939: (may be NULL)
1941: Calling Sequence of monitor:
1942: $ monitor (const char name[], const char value[], void *mctx)
1944: + name - option name string
1945: . value - option value string
1946: - mctx - optional monitoring context, as set by PetscOptionsMonitorSet()
1948: Options Database Keys:
1949: See PetscInitialize() for options related to option database monitoring.
1951: Notes:
1952: The default is to do nothing. To print the name and value of options
1953: being inserted into the database, use PetscOptionsMonitorDefault() as the monitoring routine,
1954: with a null monitoring context.
1956: Several different monitoring routines may be set by calling
1957: PetscOptionsMonitorSet() multiple times; all will be called in the
1958: order in which they were set.
1960: Level: intermediate
1962: .seealso: PetscOptionsMonitorDefault(), PetscInitialize()
1963: @*/
1964: PetscErrorCode PetscOptionsMonitorSet(PetscErrorCode (*monitor)(const char name[], const char value[], void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
1965: {
1966: PetscOptions options = defaultoptions;
1969: if (options->monitorCancel) return(0);
1970: if (options->numbermonitors >= MAXOPTIONSMONITORS) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many PetscOptions monitors set");
1971: options->monitor[options->numbermonitors] = monitor;
1972: options->monitordestroy[options->numbermonitors] = monitordestroy;
1973: options->monitorcontext[options->numbermonitors++] = (void*)mctx;
1974: return(0);
1975: }
1977: /*
1978: PetscOptionsStringToBool - Converts string to PetscBool, handles cases like "yes", "no", "true", "false", "0", "1", "off", "on".
1979: Empty string is considered as true.
1980: */
1981: PetscErrorCode PetscOptionsStringToBool(const char value[],PetscBool *a)
1982: {
1983: PetscBool istrue,isfalse;
1984: size_t len;
1988: /* PetscStrlen() returns 0 for NULL or "" */
1989: PetscStrlen(value,&len);
1990: if (!len) {*a = PETSC_TRUE; return(0);}
1991: PetscStrcasecmp(value,"TRUE",&istrue);
1992: if (istrue) {*a = PETSC_TRUE; return(0);}
1993: PetscStrcasecmp(value,"YES",&istrue);
1994: if (istrue) {*a = PETSC_TRUE; return(0);}
1995: PetscStrcasecmp(value,"1",&istrue);
1996: if (istrue) {*a = PETSC_TRUE; return(0);}
1997: PetscStrcasecmp(value,"on",&istrue);
1998: if (istrue) {*a = PETSC_TRUE; return(0);}
1999: PetscStrcasecmp(value,"FALSE",&isfalse);
2000: if (isfalse) {*a = PETSC_FALSE; return(0);}
2001: PetscStrcasecmp(value,"NO",&isfalse);
2002: if (isfalse) {*a = PETSC_FALSE; return(0);}
2003: PetscStrcasecmp(value,"0",&isfalse);
2004: if (isfalse) {*a = PETSC_FALSE; return(0);}
2005: PetscStrcasecmp(value,"off",&isfalse);
2006: if (isfalse) {*a = PETSC_FALSE; return(0);}
2007: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unknown logical value: %s",value);
2008: }
2010: /*
2011: PetscOptionsStringToInt - Converts a string to an integer value. Handles special cases such as "default" and "decide"
2012: */
2013: PetscErrorCode PetscOptionsStringToInt(const char name[],PetscInt *a)
2014: {
2016: size_t len;
2017: PetscBool decide,tdefault,mouse;
2020: PetscStrlen(name,&len);
2021: if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"character string of length zero has no numerical value");
2023: PetscStrcasecmp(name,"PETSC_DEFAULT",&tdefault);
2024: if (!tdefault) {
2025: PetscStrcasecmp(name,"DEFAULT",&tdefault);
2026: }
2027: PetscStrcasecmp(name,"PETSC_DECIDE",&decide);
2028: if (!decide) {
2029: PetscStrcasecmp(name,"DECIDE",&decide);
2030: }
2031: PetscStrcasecmp(name,"mouse",&mouse);
2033: if (tdefault) *a = PETSC_DEFAULT;
2034: else if (decide) *a = PETSC_DECIDE;
2035: else if (mouse) *a = -1;
2036: else {
2037: char *endptr;
2038: long strtolval;
2040: strtolval = strtol(name,&endptr,10);
2041: if ((size_t) (endptr - name) != len) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s has no integer value (do not include . in it)",name);
2043: #if defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE_ATOLL)
2044: (void) strtolval;
2045: *a = atoll(name);
2046: #elif defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE___INT64)
2047: (void) strtolval;
2048: *a = _atoi64(name);
2049: #else
2050: *a = (PetscInt)strtolval;
2051: #endif
2052: }
2053: return(0);
2054: }
2056: #if defined(PETSC_USE_REAL___FLOAT128)
2057: #include <quadmath.h>
2058: #endif
2060: static PetscErrorCode PetscStrtod(const char name[],PetscReal *a,char **endptr)
2061: {
2063: #if defined(PETSC_USE_REAL___FLOAT128)
2064: *a = strtoflt128(name,endptr);
2065: #else
2066: *a = (PetscReal)strtod(name,endptr);
2067: #endif
2068: return(0);
2069: }
2071: static PetscErrorCode PetscStrtoz(const char name[],PetscScalar *a,char **endptr,PetscBool *isImaginary)
2072: {
2073: PetscBool hasi = PETSC_FALSE;
2074: char *ptr;
2075: PetscReal strtoval;
2079: PetscStrtod(name,&strtoval,&ptr);
2080: if (ptr == name) {
2081: strtoval = 1.;
2082: hasi = PETSC_TRUE;
2083: if (name[0] == 'i') {
2084: ptr++;
2085: } else if (name[0] == '+' && name[1] == 'i') {
2086: ptr += 2;
2087: } else if (name[0] == '-' && name[1] == 'i') {
2088: strtoval = -1.;
2089: ptr += 2;
2090: }
2091: } else if (*ptr == 'i') {
2092: hasi = PETSC_TRUE;
2093: ptr++;
2094: }
2095: *endptr = ptr;
2096: *isImaginary = hasi;
2097: if (hasi) {
2098: #if !defined(PETSC_USE_COMPLEX)
2099: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s contains imaginary but complex not supported ",name);
2100: #else
2101: *a = PetscCMPLX(0.,strtoval);
2102: #endif
2103: } else {
2104: *a = strtoval;
2105: }
2106: return(0);
2107: }
2109: /*
2110: Converts a string to PetscReal value. Handles special cases like "default" and "decide"
2111: */
2112: PetscErrorCode PetscOptionsStringToReal(const char name[],PetscReal *a)
2113: {
2114: size_t len;
2115: PetscBool match;
2116: char *endptr;
2120: PetscStrlen(name,&len);
2121: if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"String of length zero has no numerical value");
2123: PetscStrcasecmp(name,"PETSC_DEFAULT",&match);
2124: if (!match) {
2125: PetscStrcasecmp(name,"DEFAULT",&match);
2126: }
2127: if (match) {*a = PETSC_DEFAULT; return(0);}
2129: PetscStrcasecmp(name,"PETSC_DECIDE",&match);
2130: if (!match) {
2131: PetscStrcasecmp(name,"DECIDE",&match);
2132: }
2133: if (match) {*a = PETSC_DECIDE; return(0);}
2135: PetscStrtod(name,a,&endptr);
2136: if ((size_t) (endptr - name) != len) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s has no numeric value",name);
2137: return(0);
2138: }
2140: PetscErrorCode PetscOptionsStringToScalar(const char name[],PetscScalar *a)
2141: {
2142: PetscBool imag1;
2143: size_t len;
2144: PetscScalar val = 0.;
2145: char *ptr = NULL;
2149: PetscStrlen(name,&len);
2150: if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"character string of length zero has no numerical value");
2151: PetscStrtoz(name,&val,&ptr,&imag1);
2152: #if defined(PETSC_USE_COMPLEX)
2153: if ((size_t) (ptr - name) < len) {
2154: PetscBool imag2;
2155: PetscScalar val2;
2157: PetscStrtoz(ptr,&val2,&ptr,&imag2);
2158: if (imag1 || !imag2) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s: must specify imaginary component second",name);
2159: val = PetscCMPLX(PetscRealPart(val),PetscImaginaryPart(val2));
2160: }
2161: #endif
2162: if ((size_t) (ptr - name) != len) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s has no numeric value ",name);
2163: *a = val;
2164: return(0);
2165: }
2167: /*@C
2168: PetscOptionsGetBool - Gets the Logical (true or false) value for a particular
2169: option in the database.
2171: Not Collective
2173: Input Parameters:
2174: + options - options database, use NULL for default global database
2175: . pre - the string to prepend to the name or NULL
2176: - name - the option one is seeking
2178: Output Parameter:
2179: + ivalue - the logical value to return
2180: - set - PETSC_TRUE if found, else PETSC_FALSE
2182: Level: beginner
2184: Notes:
2185: TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE
2186: FALSE, false, NO, no, and 0 all translate to PETSC_FALSE
2188: 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
2189: is equivalent to -requested_bool true
2191: If the user does not supply the option at all ivalue is NOT changed. Thus
2192: you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.
2194: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
2195: PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsGetInt(), PetscOptionsBool(),
2196: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2197: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2198: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2199: PetscOptionsFList(), PetscOptionsEList()
2200: @*/
2201: PetscErrorCode PetscOptionsGetBool(PetscOptions options,const char pre[],const char name[],PetscBool *ivalue,PetscBool *set)
2202: {
2203: const char *value;
2204: PetscBool flag;
2210: PetscOptionsFindPair(options,pre,name,&value,&flag);
2211: if (flag) {
2212: if (set) *set = PETSC_TRUE;
2213: PetscOptionsStringToBool(value, &flag);
2214: if (ivalue) *ivalue = flag;
2215: } else {
2216: if (set) *set = PETSC_FALSE;
2217: }
2218: return(0);
2219: }
2221: /*@C
2222: PetscOptionsGetEList - Puts a list of option values that a single one may be selected from
2224: Not Collective
2226: Input Parameters:
2227: + options - options database, use NULL for default global database
2228: . pre - the string to prepend to the name or NULL
2229: . opt - option name
2230: . list - the possible choices (one of these must be selected, anything else is invalid)
2231: - ntext - number of choices
2233: Output Parameter:
2234: + value - the index of the value to return (defaults to zero if the option name is given but no choice is listed)
2235: - set - PETSC_TRUE if found, else PETSC_FALSE
2237: Level: intermediate
2239: Notes:
2240: If the user does not supply the option value is NOT changed. Thus
2241: you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.
2243: See PetscOptionsFList() for when the choices are given in a PetscFunctionList()
2245: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
2246: PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2247: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2248: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2249: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2250: PetscOptionsFList(), PetscOptionsEList()
2251: @*/
2252: PetscErrorCode PetscOptionsGetEList(PetscOptions options,const char pre[],const char opt[],const char * const *list,PetscInt ntext,PetscInt *value,PetscBool *set)
2253: {
2255: size_t alen,len = 0, tlen = 0;
2256: char *svalue;
2257: PetscBool aset,flg = PETSC_FALSE;
2258: PetscInt i;
2262: for (i=0; i<ntext; i++) {
2263: PetscStrlen(list[i],&alen);
2264: if (alen > len) len = alen;
2265: tlen += len + 1;
2266: }
2267: len += 5; /* a little extra space for user mistypes */
2268: PetscMalloc1(len,&svalue);
2269: PetscOptionsGetString(options,pre,opt,svalue,len,&aset);
2270: if (aset) {
2271: PetscEListFind(ntext,list,svalue,value,&flg);
2272: if (!flg) {
2273: char *avail,*pavl;
2275: PetscMalloc1(tlen,&avail);
2276: pavl = avail;
2277: for (i=0; i<ntext; i++) {
2278: PetscStrlen(list[i],&alen);
2279: PetscStrcpy(pavl,list[i]);
2280: pavl += alen;
2281: PetscStrcpy(pavl," ");
2282: pavl += 1;
2283: }
2284: PetscStrtolower(avail);
2285: SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_USER,"Unknown option %s for -%s%s. Available options: %s",svalue,pre ? pre : "",opt+1,avail);
2286: }
2287: if (set) *set = PETSC_TRUE;
2288: } else if (set) *set = PETSC_FALSE;
2289: PetscFree(svalue);
2290: return(0);
2291: }
2293: /*@C
2294: PetscOptionsGetEnum - Gets the enum value for a particular option in the database.
2296: Not Collective
2298: Input Parameters:
2299: + options - options database, use NULL for default global database
2300: . pre - option prefix or NULL
2301: . opt - option name
2302: . list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null
2303: - defaultv - the default (current) value
2305: Output Parameter:
2306: + value - the value to return
2307: - set - PETSC_TRUE if found, else PETSC_FALSE
2309: Level: beginner
2311: Notes:
2312: If the user does not supply the option value is NOT changed. Thus
2313: you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.
2315: List is usually something like PCASMTypes or some other predefined list of enum names
2317: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), PetscOptionsGetInt(),
2318: PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
2319: PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
2320: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2321: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2322: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2323: PetscOptionsFList(), PetscOptionsEList(), PetscOptionsGetEList(), PetscOptionsEnum()
2324: @*/
2325: PetscErrorCode PetscOptionsGetEnum(PetscOptions options,const char pre[],const char opt[],const char * const *list,PetscEnum *value,PetscBool *set)
2326: {
2328: PetscInt ntext = 0,tval;
2329: PetscBool fset;
2333: while (list[ntext++]) {
2334: if (ntext > 50) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"List argument appears to be wrong or have more than 50 entries");
2335: }
2336: if (ntext < 3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"List argument must have at least two entries: typename and type prefix");
2337: ntext -= 3;
2338: PetscOptionsGetEList(options,pre,opt,list,ntext,&tval,&fset);
2339: /* with PETSC_USE_64BIT_INDICES sizeof(PetscInt) != sizeof(PetscEnum) */
2340: if (fset) *value = (PetscEnum)tval;
2341: if (set) *set = fset;
2342: return(0);
2343: }
2345: /*@C
2346: PetscOptionsGetInt - Gets the integer value for a particular option in the database.
2348: Not Collective
2350: Input Parameters:
2351: + options - options database, use NULL for default global database
2352: . pre - the string to prepend to the name or NULL
2353: - name - the option one is seeking
2355: Output Parameter:
2356: + ivalue - the integer value to return
2357: - set - PETSC_TRUE if found, else PETSC_FALSE
2359: Level: beginner
2361: Notes:
2362: If the user does not supply the option ivalue is NOT changed. Thus
2363: you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.
2365: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
2366: PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
2367: PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
2368: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2369: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2370: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2371: PetscOptionsFList(), PetscOptionsEList()
2372: @*/
2373: PetscErrorCode PetscOptionsGetInt(PetscOptions options,const char pre[],const char name[],PetscInt *ivalue,PetscBool *set)
2374: {
2375: const char *value;
2377: PetscBool flag;
2382: PetscOptionsFindPair(options,pre,name,&value,&flag);
2383: if (flag) {
2384: if (!value) {
2385: if (set) *set = PETSC_FALSE;
2386: } else {
2387: if (set) *set = PETSC_TRUE;
2388: PetscOptionsStringToInt(value,ivalue);
2389: }
2390: } else {
2391: if (set) *set = PETSC_FALSE;
2392: }
2393: return(0);
2394: }
2396: /*@C
2397: PetscOptionsGetReal - Gets the double precision value for a particular
2398: option in the database.
2400: Not Collective
2402: Input Parameters:
2403: + options - options database, use NULL for default global database
2404: . pre - string to prepend to each name or NULL
2405: - name - the option one is seeking
2407: Output Parameter:
2408: + dvalue - the double value to return
2409: - set - PETSC_TRUE if found, PETSC_FALSE if not found
2411: Notes:
2412: If the user does not supply the option dvalue is NOT changed. Thus
2413: you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.
2415: Level: beginner
2417: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2418: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(),PetscOptionsBool(),
2419: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2420: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2421: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2422: PetscOptionsFList(), PetscOptionsEList()
2423: @*/
2424: PetscErrorCode PetscOptionsGetReal(PetscOptions options,const char pre[],const char name[],PetscReal *dvalue,PetscBool *set)
2425: {
2426: const char *value;
2427: PetscBool flag;
2433: PetscOptionsFindPair(options,pre,name,&value,&flag);
2434: if (flag) {
2435: if (!value) {
2436: if (set) *set = PETSC_FALSE;
2437: } else {
2438: if (set) *set = PETSC_TRUE;
2439: PetscOptionsStringToReal(value,dvalue);
2440: }
2441: } else {
2442: if (set) *set = PETSC_FALSE;
2443: }
2444: return(0);
2445: }
2447: /*@C
2448: PetscOptionsGetScalar - Gets the scalar value for a particular
2449: option in the database.
2451: Not Collective
2453: Input Parameters:
2454: + options - options database, use NULL for default global database
2455: . pre - string to prepend to each name or NULL
2456: - name - the option one is seeking
2458: Output Parameter:
2459: + dvalue - the double value to return
2460: - set - PETSC_TRUE if found, else PETSC_FALSE
2462: Level: beginner
2464: Usage:
2465: A complex number 2+3i must be specified with NO spaces
2467: Notes:
2468: If the user does not supply the option dvalue is NOT changed. Thus
2469: you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.
2471: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2472: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2473: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2474: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2475: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2476: PetscOptionsFList(), PetscOptionsEList()
2477: @*/
2478: PetscErrorCode PetscOptionsGetScalar(PetscOptions options,const char pre[],const char name[],PetscScalar *dvalue,PetscBool *set)
2479: {
2480: const char *value;
2481: PetscBool flag;
2487: PetscOptionsFindPair(options,pre,name,&value,&flag);
2488: if (flag) {
2489: if (!value) {
2490: if (set) *set = PETSC_FALSE;
2491: } else {
2492: #if !defined(PETSC_USE_COMPLEX)
2493: PetscOptionsStringToReal(value,dvalue);
2494: #else
2495: PetscOptionsStringToScalar(value,dvalue);
2496: #endif
2497: if (set) *set = PETSC_TRUE;
2498: }
2499: } else { /* flag */
2500: if (set) *set = PETSC_FALSE;
2501: }
2502: return(0);
2503: }
2505: /*@C
2506: PetscOptionsGetString - Gets the string value for a particular option in
2507: the database.
2509: Not Collective
2511: Input Parameters:
2512: + options - options database, use NULL for default global database
2513: . pre - string to prepend to name or NULL
2514: . name - the option one is seeking
2515: - len - maximum length of the string including null termination
2517: Output Parameters:
2518: + string - location to copy string
2519: - set - PETSC_TRUE if found, else PETSC_FALSE
2521: Level: beginner
2523: Fortran Note:
2524: The Fortran interface is slightly different from the C/C++
2525: interface (len is not used). Sample usage in Fortran follows
2526: .vb
2527: character *20 string
2528: PetscErrorCode ierr
2529: PetscBool set
2530: call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s',string,set,ierr)
2531: .ve
2533: Notes:
2534: 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
2536: If the user does not use the option then the string is not changed. Thus
2537: you should ALWAYS initialize the string if you access it without first checking if the set flag is true.
2539: Note:
2540: 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).
2542: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
2543: PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2544: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2545: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2546: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2547: PetscOptionsFList(), PetscOptionsEList()
2548: @*/
2549: PetscErrorCode PetscOptionsGetString(PetscOptions options,const char pre[],const char name[],char string[],size_t len,PetscBool *set)
2550: {
2551: const char *value;
2552: PetscBool flag;
2558: PetscOptionsFindPair(options,pre,name,&value,&flag);
2559: if (!flag) {
2560: if (set) *set = PETSC_FALSE;
2561: } else {
2562: if (set) *set = PETSC_TRUE;
2563: if (value) {
2564: PetscStrncpy(string,value,len);
2565: } else {
2566: PetscArrayzero(string,len);
2567: }
2568: }
2569: return(0);
2570: }
2572: char *PetscOptionsGetStringMatlab(PetscOptions options,const char pre[],const char name[])
2573: {
2574: const char *value;
2575: PetscBool flag;
2579: PetscOptionsFindPair(options,pre,name,&value,&flag);if (ierr) PetscFunctionReturn(NULL);
2580: if (flag) PetscFunctionReturn((char*)value);
2581: else PetscFunctionReturn(NULL);
2582: }
2584: /*@C
2585: PetscOptionsGetBoolArray - Gets an array of Logical (true or false) values for a particular
2586: option in the database. The values must be separated with commas with
2587: no intervening spaces.
2589: Not Collective
2591: Input Parameters:
2592: + options - options database, use NULL for default global database
2593: . pre - string to prepend to each name or NULL
2594: . name - the option one is seeking
2595: - nmax - maximum number of values to retrieve
2597: Output Parameter:
2598: + dvalue - the integer values to return
2599: . nmax - actual number of values retreived
2600: - set - PETSC_TRUE if found, else PETSC_FALSE
2602: Level: beginner
2604: Notes:
2605: TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE
2606: FALSE, false, NO, no, and 0 all translate to PETSC_FALSE
2608: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2609: PetscOptionsGetString(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2610: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2611: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2612: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2613: PetscOptionsFList(), PetscOptionsEList()
2614: @*/
2615: PetscErrorCode PetscOptionsGetBoolArray(PetscOptions options,const char pre[],const char name[],PetscBool dvalue[],PetscInt *nmax,PetscBool *set)
2616: {
2617: const char *svalue;
2618: char *value;
2620: PetscInt n = 0;
2621: PetscBool flag;
2622: PetscToken token;
2629: PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2630: if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; return(0);}
2631: if (set) *set = PETSC_TRUE;
2632: PetscTokenCreate(svalue,',',&token);
2633: PetscTokenFind(token,&value);
2634: while (value && n < *nmax) {
2635: PetscOptionsStringToBool(value,dvalue);
2636: PetscTokenFind(token,&value);
2637: dvalue++;
2638: n++;
2639: }
2640: PetscTokenDestroy(&token);
2641: *nmax = n;
2642: return(0);
2643: }
2645: /*@C
2646: PetscOptionsGetEnumArray - Gets an array of enum values for a particular option in the database.
2648: Not Collective
2650: Input Parameters:
2651: + options - options database, use NULL for default global database
2652: . pre - option prefix or NULL
2653: . name - option name
2654: . list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null
2655: - nmax - maximum number of values to retrieve
2657: Output Parameters:
2658: + ivalue - the enum values to return
2659: . nmax - actual number of values retreived
2660: - set - PETSC_TRUE if found, else PETSC_FALSE
2662: Level: beginner
2664: Notes:
2665: The array must be passed as a comma separated list.
2667: There must be no intervening spaces between the values.
2669: list is usually something like PCASMTypes or some other predefined list of enum names.
2671: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), PetscOptionsGetInt(),
2672: PetscOptionsGetEnum(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
2673: PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(), PetscOptionsName(),
2674: PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), PetscOptionsStringArray(),PetscOptionsRealArray(),
2675: PetscOptionsScalar(), PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2676: PetscOptionsFList(), PetscOptionsEList(), PetscOptionsGetEList(), PetscOptionsEnum()
2677: @*/
2678: PetscErrorCode PetscOptionsGetEnumArray(PetscOptions options,const char pre[],const char name[],const char *const *list,PetscEnum ivalue[],PetscInt *nmax,PetscBool *set)
2679: {
2680: const char *svalue;
2681: char *value;
2682: PetscInt n = 0;
2683: PetscEnum evalue;
2684: PetscBool flag;
2685: PetscToken token;
2694: PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2695: if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; return(0);}
2696: if (set) *set = PETSC_TRUE;
2697: PetscTokenCreate(svalue,',',&token);
2698: PetscTokenFind(token,&value);
2699: while (value && n < *nmax) {
2700: PetscEnumFind(list,value,&evalue,&flag);
2701: if (!flag) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_USER,"Unknown enum value '%s' for -%s%s",svalue,pre ? pre : "",name+1);
2702: ivalue[n++] = evalue;
2703: PetscTokenFind(token,&value);
2704: }
2705: PetscTokenDestroy(&token);
2706: *nmax = n;
2707: return(0);
2708: }
2710: /*@C
2711: PetscOptionsGetIntArray - Gets an array of integer values for a particular
2712: option in the database.
2714: Not Collective
2716: Input Parameters:
2717: + options - options database, use NULL for default global database
2718: . pre - string to prepend to each name or NULL
2719: . name - the option one is seeking
2720: - nmax - maximum number of values to retrieve
2722: Output Parameter:
2723: + ivalue - the integer values to return
2724: . nmax - actual number of values retreived
2725: - set - PETSC_TRUE if found, else PETSC_FALSE
2727: Level: beginner
2729: Notes:
2730: The array can be passed as
2731: a comma separated list: 0,1,2,3,4,5,6,7
2732: a range (start-end+1): 0-8
2733: a range with given increment (start-end+1:inc): 0-7:2
2734: a combination of values and ranges separated by commas: 0,1-8,8-15:2
2736: There must be no intervening spaces between the values.
2738: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2739: PetscOptionsGetString(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2740: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2741: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2742: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2743: PetscOptionsFList(), PetscOptionsEList()
2744: @*/
2745: PetscErrorCode PetscOptionsGetIntArray(PetscOptions options,const char pre[],const char name[],PetscInt ivalue[],PetscInt *nmax,PetscBool *set)
2746: {
2747: const char *svalue;
2748: char *value;
2750: PetscInt n = 0,i,j,start,end,inc,nvalues;
2751: size_t len;
2752: PetscBool flag,foundrange;
2753: PetscToken token;
2760: PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2761: if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; return(0);}
2762: if (set) *set = PETSC_TRUE;
2763: PetscTokenCreate(svalue,',',&token);
2764: PetscTokenFind(token,&value);
2765: while (value && n < *nmax) {
2766: /* look for form d-D where d and D are integers */
2767: foundrange = PETSC_FALSE;
2768: PetscStrlen(value,&len);
2769: if (value[0] == '-') i=2;
2770: else i=1;
2771: for (;i<(int)len; i++) {
2772: if (value[i] == '-') {
2773: if (i == (int)len-1) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry %s\n",n,value);
2774: value[i] = 0;
2776: PetscOptionsStringToInt(value,&start);
2777: inc = 1;
2778: j = i+1;
2779: for (;j<(int)len; j++) {
2780: if (value[j] == ':') {
2781: value[j] = 0;
2783: PetscOptionsStringToInt(value+j+1,&inc);
2784: if (inc <= 0) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry,%s cannot have negative increment",n,value+j+1);
2785: break;
2786: }
2787: }
2788: PetscOptionsStringToInt(value+i+1,&end);
2789: if (end <= start) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry, %s-%s cannot have decreasing list",n,value,value+i+1);
2790: nvalues = (end-start)/inc + (end-start)%inc;
2791: if (n + nvalues > *nmax) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry, not enough space left in array (%D) to contain entire range from %D to %D",n,*nmax-n,start,end);
2792: for (;start<end; start+=inc) {
2793: *ivalue = start; ivalue++;n++;
2794: }
2795: foundrange = PETSC_TRUE;
2796: break;
2797: }
2798: }
2799: if (!foundrange) {
2800: PetscOptionsStringToInt(value,ivalue);
2801: ivalue++;
2802: n++;
2803: }
2804: PetscTokenFind(token,&value);
2805: }
2806: PetscTokenDestroy(&token);
2807: *nmax = n;
2808: return(0);
2809: }
2811: /*@C
2812: PetscOptionsGetRealArray - Gets an array of double precision values for a
2813: particular option in the database. The values must be separated with
2814: commas with no intervening spaces.
2816: Not Collective
2818: Input Parameters:
2819: + options - options database, use NULL for default global database
2820: . pre - string to prepend to each name or NULL
2821: . name - the option one is seeking
2822: - nmax - maximum number of values to retrieve
2824: Output Parameters:
2825: + dvalue - the double values to return
2826: . nmax - actual number of values retreived
2827: - set - PETSC_TRUE if found, else PETSC_FALSE
2829: Level: beginner
2831: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2832: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
2833: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2834: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2835: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2836: PetscOptionsFList(), PetscOptionsEList()
2837: @*/
2838: PetscErrorCode PetscOptionsGetRealArray(PetscOptions options,const char pre[],const char name[],PetscReal dvalue[],PetscInt *nmax,PetscBool *set)
2839: {
2840: const char *svalue;
2841: char *value;
2843: PetscInt n = 0;
2844: PetscBool flag;
2845: PetscToken token;
2852: PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2853: if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; return(0);}
2854: if (set) *set = PETSC_TRUE;
2855: PetscTokenCreate(svalue,',',&token);
2856: PetscTokenFind(token,&value);
2857: while (value && n < *nmax) {
2858: PetscOptionsStringToReal(value,dvalue++);
2859: PetscTokenFind(token,&value);
2860: n++;
2861: }
2862: PetscTokenDestroy(&token);
2863: *nmax = n;
2864: return(0);
2865: }
2867: /*@C
2868: PetscOptionsGetScalarArray - Gets an array of scalars for a
2869: particular option in the database. The values must be separated with
2870: commas with no intervening spaces.
2872: Not Collective
2874: Input Parameters:
2875: + options - options database, use NULL for default global database
2876: . pre - string to prepend to each name or NULL
2877: . name - the option one is seeking
2878: - nmax - maximum number of values to retrieve
2880: Output Parameters:
2881: + dvalue - the scalar values to return
2882: . nmax - actual number of values retreived
2883: - set - PETSC_TRUE if found, else PETSC_FALSE
2885: Level: beginner
2887: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2888: PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
2889: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2890: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2891: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2892: PetscOptionsFList(), PetscOptionsEList()
2893: @*/
2894: PetscErrorCode PetscOptionsGetScalarArray(PetscOptions options,const char pre[],const char name[],PetscScalar dvalue[],PetscInt *nmax,PetscBool *set)
2895: {
2896: const char *svalue;
2897: char *value;
2899: PetscInt n = 0;
2900: PetscBool flag;
2901: PetscToken token;
2908: PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2909: if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; return(0);}
2910: if (set) *set = PETSC_TRUE;
2911: PetscTokenCreate(svalue,',',&token);
2912: PetscTokenFind(token,&value);
2913: while (value && n < *nmax) {
2914: PetscOptionsStringToScalar(value,dvalue++);
2915: PetscTokenFind(token,&value);
2916: n++;
2917: }
2918: PetscTokenDestroy(&token);
2919: *nmax = n;
2920: return(0);
2921: }
2923: /*@C
2924: PetscOptionsGetStringArray - Gets an array of string values for a particular
2925: option in the database. The values must be separated with commas with
2926: no intervening spaces.
2928: Not Collective
2930: Input Parameters:
2931: + options - options database, use NULL for default global database
2932: . pre - string to prepend to name or NULL
2933: . name - the option one is seeking
2934: - nmax - maximum number of strings
2936: Output Parameters:
2937: + strings - location to copy strings
2938: . nmax - the number of strings found
2939: - set - PETSC_TRUE if found, else PETSC_FALSE
2941: Level: beginner
2943: Notes:
2944: The nmax parameter is used for both input and output.
2946: The user should pass in an array of pointers to char, to hold all the
2947: strings returned by this function.
2949: The user is responsible for deallocating the strings that are
2950: returned. The Fortran interface for this routine is not supported.
2952: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
2953: PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2954: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2955: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2956: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2957: PetscOptionsFList(), PetscOptionsEList()
2958: @*/
2959: PetscErrorCode PetscOptionsGetStringArray(PetscOptions options,const char pre[],const char name[],char *strings[],PetscInt *nmax,PetscBool *set)
2960: {
2961: const char *svalue;
2962: char *value;
2964: PetscInt n = 0;
2965: PetscBool flag;
2966: PetscToken token;
2973: PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2974: if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; return(0);}
2975: if (set) *set = PETSC_TRUE;
2976: PetscTokenCreate(svalue,',',&token);
2977: PetscTokenFind(token,&value);
2978: while (value && n < *nmax) {
2979: PetscStrallocpy(value,&strings[n]);
2980: PetscTokenFind(token,&value);
2981: n++;
2982: }
2983: PetscTokenDestroy(&token);
2984: *nmax = n;
2985: return(0);
2986: }
2988: /*@C
2989: PetscOptionsDeprecated - mark an option as deprecated, optionally replacing it with a new one
2991: Prints a deprecation warning, unless an option is supplied to suppress.
2993: Logically Collective
2995: Input Parameters:
2996: + pre - string to prepend to name or NULL
2997: . oldname - the old, deprecated option
2998: . newname - the new option, or NULL if option is purely removed
2999: . version - a string describing the version of first deprecation, e.g. "3.9"
3000: - info - additional information string, or NULL.
3002: Options Database Keys:
3003: . -options_suppress_deprecated_warnings - do not print deprecation warnings
3005: Notes:
3006: Must be called between PetscOptionsBegin() (or PetscObjectOptionsBegin()) and PetscOptionsEnd().
3007: Only the proces of rank zero that owns the PetscOptionsItems are argument (managed by PetscOptionsBegin() or
3008: PetscObjectOptionsBegin() prints the information
3009: If newname is provided, the old option is replaced. Otherwise, it remains
3010: in the options database.
3011: If an option is not replaced, the info argument should be used to advise the user
3012: on how to proceed.
3013: There is a limit on the length of the warning printed, so very long strings
3014: provided as info may be truncated.
3016: Level: developer
3018: .seealso: PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsScalar(), PetscOptionsBool(), PetscOptionsString(), PetscOptionsSetValue()
3020: @*/
3021: PetscErrorCode PetscOptionsDeprecated_Private(PetscOptionItems *PetscOptionsObject,const char oldname[],const char newname[],const char version[],const char info[])
3022: {
3023: PetscErrorCode ierr;
3024: PetscBool found,quiet;
3025: const char *value;
3026: const char * const quietopt="-options_suppress_deprecated_warnings";
3027: char msg[4096];
3028: char *prefix = NULL;
3029: PetscOptions options = NULL;
3030: MPI_Comm comm = PETSC_COMM_SELF;
3035: if (PetscOptionsObject) {
3036: prefix = PetscOptionsObject->prefix;
3037: options = PetscOptionsObject->options;
3038: comm = PetscOptionsObject->comm;
3039: }
3040: PetscOptionsFindPair(options,prefix,oldname,&value,&found);
3041: if (found) {
3042: if (newname) {
3043: if (prefix) {
3044: PetscOptionsPrefixPush(options,prefix);
3045: }
3046: PetscOptionsSetValue(options,newname,value);
3047: if (prefix) {
3048: PetscOptionsPrefixPop(options);
3049: }
3050: PetscOptionsClearValue(options,oldname);
3051: }
3052: quiet = PETSC_FALSE;
3053: PetscOptionsGetBool(options,NULL,quietopt,&quiet,NULL);
3054: if (!quiet) {
3055: PetscStrcpy(msg,"** PETSc DEPRECATION WARNING ** : the option ");
3056: PetscStrcat(msg,oldname);
3057: PetscStrcat(msg," is deprecated as of version ");
3058: PetscStrcat(msg,version);
3059: PetscStrcat(msg," and will be removed in a future release.");
3060: if (newname) {
3061: PetscStrcat(msg," Please use the option ");
3062: PetscStrcat(msg,newname);
3063: PetscStrcat(msg," instead.");
3064: }
3065: if (info) {
3066: PetscStrcat(msg," ");
3067: PetscStrcat(msg,info);
3068: }
3069: PetscStrcat(msg," (Silence this warning with ");
3070: PetscStrcat(msg,quietopt);
3071: PetscStrcat(msg,")\n");
3072: PetscPrintf(comm,msg);
3073: }
3074: }
3075: return(0);
3076: }