Actual source code: options.c

  1: /* Define Feature test macros to make sure atoll is available (SVr4, POSIX.1-2001, 4.3BSD, C99), not in (C89 and POSIX.1-1996) */
  2: #define PETSC_DESIRE_FEATURE_TEST_MACROS /* for atoll() */

  4: /*
  5:    These routines simplify the use of command line, file options, etc., and are used to manipulate the options database.
  6:    This provides the low-level interface, the high level interface is in aoptions.c

  8:    Some routines use regular malloc and free because it cannot know  what malloc is requested with the
  9:    options database until it has already processed the input.
 10: */

 12: #include <petsc/private/petscimpl.h>
 13: #include <petscviewer.h>
 14: #include <ctype.h>
 15: #if defined(PETSC_HAVE_MALLOC_H)
 16: #include <malloc.h>
 17: #endif
 18: #if defined(PETSC_HAVE_STRINGS_H)
 19: #  include <strings.h>          /* strcasecmp */
 20: #endif

 22: #if defined(PETSC_HAVE_STRCASECMP)
 23: #define PetscOptNameCmp(a,b) strcasecmp(a,b)
 24: #elif defined(PETSC_HAVE_STRICMP)
 25: #define PetscOptNameCmp(a,b) stricmp(a,b)
 26: #else
 27: #define PetscOptNameCmp(a,b) Error_strcasecmp_not_found
 28: #endif

 30: #include <petsc/private/hashtable.h>

 32: /* This assumes ASCII encoding and ignores locale settings */
 33: /* Using tolower() is about 2X slower in microbenchmarks   */
 34: static inline int PetscToLower(int c)
 35: {
 36:   return ((c >= 'A') & (c <= 'Z')) ? c + 'a' - 'A' : c;
 37: }

 39: /* Bob Jenkins's one at a time hash function (case-insensitive) */
 40: static inline unsigned int PetscOptHash(const char key[])
 41: {
 42:   unsigned int hash = 0;
 43:   while (*key) {
 44:     hash += PetscToLower(*key++);
 45:     hash += hash << 10;
 46:     hash ^= hash >>  6;
 47:   }
 48:   hash += hash <<  3;
 49:   hash ^= hash >> 11;
 50:   hash += hash << 15;
 51:   return hash;
 52: }

 54: static inline int PetscOptEqual(const char a[],const char b[])
 55: {
 56:   return !PetscOptNameCmp(a,b);
 57: }

 59: KHASH_INIT(HO, kh_cstr_t, int, 1, PetscOptHash, PetscOptEqual)

 61: /*
 62:     This table holds all the options set by the user. For simplicity, we use a static size database
 63: */
 64: #define MAXOPTNAME PETSC_MAX_OPTION_NAME
 65: #define MAXOPTIONS 512
 66: #define MAXALIASES  25
 67: #define MAXPREFIXES 25
 68: #define MAXOPTIONSMONITORS 5

 70: struct  _n_PetscOptions {
 71:   PetscOptions   previous;
 72:   int            N;                    /* number of options */
 73:   char           *names[MAXOPTIONS];   /* option names */
 74:   char           *values[MAXOPTIONS];  /* option values */
 75:   PetscBool      used[MAXOPTIONS];     /* flag option use */
 76:   PetscBool      precedentProcessed;

 78:   /* Hash table */
 79:   khash_t(HO)    *ht;

 81:   /* Prefixes */
 82:   int            prefixind;
 83:   int            prefixstack[MAXPREFIXES];
 84:   char           prefix[MAXOPTNAME];

 86:   /* Aliases */
 87:   int            Naliases;                   /* number or aliases */
 88:   char           *aliases1[MAXALIASES];      /* aliased */
 89:   char           *aliases2[MAXALIASES];      /* aliasee */

 91:   /* Help */
 92:   PetscBool      help;       /* flag whether "-help" is in the database */
 93:   PetscBool      help_intro; /* flag whether "-help intro" is in the database */

 95:   /* Monitors */
 96:   PetscBool      monitorFromOptions, monitorCancel;
 97:   PetscErrorCode (*monitor[MAXOPTIONSMONITORS])(const char[],const char[],void*); /* returns control to user after */
 98:   PetscErrorCode (*monitordestroy[MAXOPTIONSMONITORS])(void**);         /* */
 99:   void           *monitorcontext[MAXOPTIONSMONITORS];                  /* to pass arbitrary user data into monitor */
100:   PetscInt       numbermonitors;                                       /* to, for instance, detect options being set */
101: };

103: static PetscOptions defaultoptions = NULL;  /* the options database routines query this object for options */

105: /* list of options which preceed others, i.e., are processed in PetscOptionsProcessPrecedentFlags() */
106: static const char *precedentOptions[] = {"-options_monitor","-options_monitor_cancel","-help","-skip_petscrc"};
107: enum PetscPrecedentOption {PO_OPTIONS_MONITOR,PO_OPTIONS_MONITOR_CANCEL,PO_HELP,PO_SKIP_PETSCRC,PO_NUM};

109: static PetscErrorCode PetscOptionsSetValue_Private(PetscOptions,const char[],const char[],int*);

111: /*
112:     Options events monitor
113: */
114: static PetscErrorCode PetscOptionsMonitor(PetscOptions options,const char name[],const char value[])
115: {
116:   if (!value) value = "";
117:   if (options->monitorFromOptions) PetscOptionsMonitorDefault(name,value,NULL);
118:   for (PetscInt i=0; i<options->numbermonitors; i++) (*options->monitor[i])(name,value,options->monitorcontext[i]);
119:   return 0;
120: }

122: /*@
123:    PetscOptionsCreate - Creates an empty options database.

125:    Logically collective

127:    Output Parameter:
128: .  options - Options database object

130:    Level: advanced

132:    Developer Note: We may want eventually to pass a MPI_Comm to determine the ownership of the object

134: .seealso: PetscOptionsDestroy(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsInsert(), PetscOptionsSetValue()
135: @*/
136: PetscErrorCode PetscOptionsCreate(PetscOptions *options)
137: {
139:   *options = (PetscOptions)calloc(1,sizeof(**options));
141:   return 0;
142: }

144: /*@
145:     PetscOptionsDestroy - Destroys an option database.

147:     Logically collective on whatever communicator was associated with the call to PetscOptionsCreate()

149:   Input Parameter:
150: .  options - the PetscOptions object

152:    Level: advanced

154: .seealso: PetscOptionsInsert(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsInsert(), PetscOptionsSetValue()
155: @*/
156: PetscErrorCode PetscOptionsDestroy(PetscOptions *options)
157: {
158:   if (!*options) return 0;
160:   PetscOptionsClear(*options);
161:   /* XXX what about monitors ? */
162:   free(*options);
163:   *options = NULL;
164:   return 0;
165: }

167: /*
168:     PetscOptionsCreateDefault - Creates the default global options database
169: */
170: PetscErrorCode PetscOptionsCreateDefault(void)
171: {
172:   if (PetscUnlikely(!defaultoptions)) PetscOptionsCreate(&defaultoptions);
173:   return 0;
174: }

176: /*@
177:       PetscOptionsPush - Push a new PetscOptions object as the default provider of options
178:                          Allows using different parts of a code to use different options databases

180:   Logically Collective

182:   Input Parameter:
183: .   opt - the options obtained with PetscOptionsCreate()

185:   Notes:
186:   Use PetscOptionsPop() to return to the previous default options database

188:   The collectivity of this routine is complex; only the MPI processes that call this routine will
189:   have the affect of these options. If some processes that create objects call this routine and others do
190:   not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
191:   on different ranks.

193:    Level: advanced

195: .seealso: PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsInsert(), PetscOptionsSetValue(), PetscOptionsLeft()

197: @*/
198: PetscErrorCode PetscOptionsPush(PetscOptions opt)
199: {
200:   PetscOptionsCreateDefault();
201:   opt->previous  = defaultoptions;
202:   defaultoptions = opt;
203:   return 0;
204: }

206: /*@
207:       PetscOptionsPop - Pop the most recent PetscOptionsPush() to return to the previous default options

209:       Logically collective on whatever communicator was associated with the call to PetscOptionsCreate()

211:   Notes:
212:   Use PetscOptionsPop() to return to the previous default options database
213:   Allows using different parts of a code to use different options databases

215:    Level: advanced

217: .seealso: PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsInsert(), PetscOptionsSetValue(), PetscOptionsLeft()

219: @*/
220: PetscErrorCode PetscOptionsPop(void)
221: {
222:   PetscOptions current = defaultoptions;

226:   defaultoptions    = defaultoptions->previous;
227:   current->previous = NULL;
228:   return 0;
229: }

231: /*
232:     PetscOptionsDestroyDefault - Destroys the default global options database
233: */
234: PetscErrorCode PetscOptionsDestroyDefault(void)
235: {
236:   if (!defaultoptions) return 0;
237:   /* Destroy any options that the user forgot to pop */
238:   while (defaultoptions->previous) {
239:     PetscOptions tmp = defaultoptions;

241:     PetscOptionsPop();
242:     PetscOptionsDestroy(&tmp);
243:   }
244:   PetscOptionsDestroy(&defaultoptions);
245:   return 0;
246: }

248: /*@C
249:    PetscOptionsValidKey - PETSc Options database keys must begin with one or two dashes (-) followed by a letter.

251:    Not collective

253:    Input Parameter:
254: .  key - string to check if valid

256:    Output Parameter:
257: .  valid - PETSC_TRUE if a valid key

259:    Level: intermediate
260: @*/
261: PetscErrorCode PetscOptionsValidKey(const char key[],PetscBool *valid)
262: {
263:   char *ptr;

267:   *valid = PETSC_FALSE;
268:   if (!key) return 0;
269:   if (key[0] != '-') return 0;
270:   if (key[1] == '-') key++;
271:   if (!isalpha((int)key[1])) return 0;
272:   (void) strtod(key,&ptr);
273:   if (ptr != key && !(*ptr == '_' || isalnum((int)*ptr))) return 0;
274:   *valid = PETSC_TRUE;
275:   return 0;
276: }

278: /*@C
279:    PetscOptionsInsertString - Inserts options into the database from a string

281:    Logically Collective

283:    Input Parameters:
284: +  options - options object
285: -  in_str - string that contains options separated by blanks

287:    Level: intermediate

289:   The collectivity of this routine is complex; only the MPI processes that call this routine will
290:   have the affect of these options. If some processes that create objects call this routine and others do
291:   not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
292:   on different ranks.

294:    Contributed by Boyana Norris

296: .seealso: PetscOptionsSetValue(), PetscOptionsView(), PetscOptionsHasName(), PetscOptionsGetInt(),
297:           PetscOptionsGetReal(), PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
298:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
299:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
300:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
301:           PetscOptionsFList(), PetscOptionsEList(), PetscOptionsInsertFile()
302: @*/
303: PetscErrorCode PetscOptionsInsertString(PetscOptions options,const char in_str[])
304: {
305:   MPI_Comm       comm = PETSC_COMM_SELF;
306:   char           *first,*second;
307:   PetscToken     token;

309:   PetscTokenCreate(in_str,' ',&token);
310:   PetscTokenFind(token,&first);
311:   while (first) {
312:     PetscBool isfile,isfileyaml,isstringyaml,ispush,ispop,key;
313:     PetscStrcasecmp(first,"-options_file",&isfile);
314:     PetscStrcasecmp(first,"-options_file_yaml",&isfileyaml);
315:     PetscStrcasecmp(first,"-options_string_yaml",&isstringyaml);
316:     PetscStrcasecmp(first,"-prefix_push",&ispush);
317:     PetscStrcasecmp(first,"-prefix_pop",&ispop);
318:     PetscOptionsValidKey(first,&key);
319:     if (!key) {
320:       PetscTokenFind(token,&first);
321:     } else if (isfile) {
322:       PetscTokenFind(token,&second);
323:       PetscOptionsInsertFile(comm,options,second,PETSC_TRUE);
324:       PetscTokenFind(token,&first);
325:     } else if (isfileyaml) {
326:       PetscTokenFind(token,&second);
327:       PetscOptionsInsertFileYAML(comm,options,second,PETSC_TRUE);
328:       PetscTokenFind(token,&first);
329:     } else if (isstringyaml) {
330:       PetscTokenFind(token,&second);
331:       PetscOptionsInsertStringYAML(options,second);
332:       PetscTokenFind(token,&first);
333:     } else if (ispush) {
334:       PetscTokenFind(token,&second);
335:       PetscOptionsPrefixPush(options,second);
336:       PetscTokenFind(token,&first);
337:     } else if (ispop) {
338:       PetscOptionsPrefixPop(options);
339:       PetscTokenFind(token,&first);
340:     } else {
341:       PetscTokenFind(token,&second);
342:       PetscOptionsValidKey(second,&key);
343:       if (!key) {
344:         PetscOptionsSetValue(options,first,second);
345:         PetscTokenFind(token,&first);
346:       } else {
347:         PetscOptionsSetValue(options,first,NULL);
348:         first = second;
349:       }
350:     }
351:   }
352:   PetscTokenDestroy(&token);
353:   return 0;
354: }

356: /*
357:     Returns a line (ended by a \n, \r or null character of any length. Result should be freed with free()
358: */
359: static char *Petscgetline(FILE * f)
360: {
361:   size_t size  = 0;
362:   size_t len   = 0;
363:   size_t last  = 0;
364:   char   *buf  = NULL;

366:   if (feof(f)) return NULL;
367:   do {
368:     size += 1024; /* BUFSIZ is defined as "the optimal read size for this platform" */
369:     buf   = (char*)realloc((void*)buf,size); /* realloc(NULL,n) is the same as malloc(n) */
370:     /* Actually do the read. Note that fgets puts a terminal '\0' on the
371:     end of the string, so we make sure we overwrite this */
372:     if (!fgets(buf+len,1024,f)) buf[len]=0;
373:     PetscStrlen(buf,&len);
374:     last = len - 1;
375:   } while (!feof(f) && buf[last] != '\n' && buf[last] != '\r');
376:   if (len) return buf;
377:   free(buf);
378:   return NULL;
379: }

381: static PetscErrorCode PetscOptionsFilename(MPI_Comm comm,const char file[],char filename[PETSC_MAX_PATH_LEN],PetscBool *yaml)
382: {
383:   char           fname[PETSC_MAX_PATH_LEN+8],path[PETSC_MAX_PATH_LEN+8],*tail;

385:   *yaml = PETSC_FALSE;
386:   PetscStrreplace(comm,file,fname,sizeof(fname));
387:   PetscFixFilename(fname,path);
388:   PetscStrendswith(path,":yaml",yaml);
389:   if (*yaml) {
390:     PetscStrrchr(path,':',&tail);
391:     tail[-1] = 0; /* remove ":yaml" suffix from path */
392:   }
393:   PetscStrncpy(filename,path,PETSC_MAX_PATH_LEN);
394:   /* check for standard YAML and JSON filename extensions */
395:   if (!*yaml) PetscStrendswith(filename,".yaml",yaml);
396:   if (!*yaml) PetscStrendswith(filename,".yml", yaml);
397:   if (!*yaml) PetscStrendswith(filename,".json",yaml);
398:   if (!*yaml) { /* check file contents */
399:     PetscMPIInt rank;
400:     MPI_Comm_rank(comm,&rank);
401:     if (rank == 0) {
402:       FILE *fh = fopen(filename,"r");
403:       if (fh) {
404:         char buf[6] = "";
405:         if (fread(buf,1,6,fh) > 0) {
406:           PetscStrncmp(buf,"%YAML ",6,yaml);  /* check for '%YAML' tag */
407:           if (!*yaml) PetscStrncmp(buf,"---",3,yaml);  /* check for document start */
408:         }
409:         (void)fclose(fh);
410:       }
411:     }
412:     MPI_Bcast(yaml,1,MPIU_BOOL,0,comm);
413:   }
414:   return 0;
415: }

417: static PetscErrorCode PetscOptionsInsertFilePetsc(MPI_Comm comm,PetscOptions options,const char file[],PetscBool require)
418: {
419:   char           *string,*vstring = NULL,*astring = NULL,*packed = NULL;
420:   char           *tokens[4];
421:   size_t         i,len,bytes;
422:   FILE           *fd;
423:   PetscToken     token=NULL;
424:   int            err;
425:   char           *cmatch;
426:   const char     cmt='#';
427:   PetscInt       line=1;
428:   PetscMPIInt    rank,cnt=0,acnt=0,counts[2];
429:   PetscBool      isdir,alias=PETSC_FALSE,valid;

431:   PetscMemzero(tokens,sizeof(tokens));
432:   MPI_Comm_rank(comm,&rank);
433:   if (rank == 0) {
434:     char fpath[PETSC_MAX_PATH_LEN];
435:     char fname[PETSC_MAX_PATH_LEN];

437:     PetscStrreplace(PETSC_COMM_SELF,file,fpath,sizeof(fpath));
438:     PetscFixFilename(fpath,fname);

440:     fd   = fopen(fname,"r");
441:     PetscTestDirectory(fname,'r',&isdir);
443:     if (fd && !isdir) {
444:       PetscSegBuffer vseg,aseg;
445:       PetscSegBufferCreate(1,4000,&vseg);
446:       PetscSegBufferCreate(1,2000,&aseg);

448:       /* the following line will not work when opening initial files (like .petscrc) since info is not yet set */
449:       PetscInfo(NULL,"Opened options file %s\n",file);

451:       while ((string = Petscgetline(fd))) {
452:         /* eliminate comments from each line */
453:         PetscStrchr(string,cmt,&cmatch);
454:         if (cmatch) *cmatch = 0;
455:         PetscStrlen(string,&len);
456:         /* replace tabs, ^M, \n with " " */
457:         for (i=0; i<len; i++) {
458:           if (string[i] == '\t' || string[i] == '\r' || string[i] == '\n') {
459:             string[i] = ' ';
460:           }
461:         }
462:         PetscTokenCreate(string,' ',&token);
463:         PetscTokenFind(token,&tokens[0]);
464:         if (!tokens[0]) {
465:           goto destroy;
466:         } else if (!tokens[0][0]) { /* if token 0 is empty (string begins with spaces), redo */
467:           PetscTokenFind(token,&tokens[0]);
468:         }
469:         for (i=1; i<4; i++) {
470:           PetscTokenFind(token,&tokens[i]);
471:         }
472:         if (!tokens[0]) {
473:           goto destroy;
474:         } else if (tokens[0][0] == '-') {
475:           PetscOptionsValidKey(tokens[0],&valid);
477:           PetscStrlen(tokens[0],&len);
478:           PetscSegBufferGet(vseg,len+1,&vstring);
479:           PetscArraycpy(vstring,tokens[0],len);
480:           vstring[len] = ' ';
481:           if (tokens[1]) {
482:             PetscOptionsValidKey(tokens[1],&valid);
484:             PetscStrlen(tokens[1],&len);
485:             PetscSegBufferGet(vseg,len+3,&vstring);
486:             vstring[0] = '"';
487:             PetscArraycpy(vstring+1,tokens[1],len);
488:             vstring[len+1] = '"';
489:             vstring[len+2] = ' ';
490:           }
491:         } else {
492:           PetscStrcasecmp(tokens[0],"alias",&alias);
493:           if (alias) {
494:             PetscOptionsValidKey(tokens[1],&valid);
497:             PetscOptionsValidKey(tokens[2],&valid);
499:             PetscStrlen(tokens[1],&len);
500:             PetscSegBufferGet(aseg,len+1,&astring);
501:             PetscArraycpy(astring,tokens[1],len);
502:             astring[len] = ' ';

504:             PetscStrlen(tokens[2],&len);
505:             PetscSegBufferGet(aseg,len+1,&astring);
506:             PetscArraycpy(astring,tokens[2],len);
507:             astring[len] = ' ';
508:           } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unknown first token in options file %s line %" PetscInt_FMT ": %s",fname,line,tokens[0]);
509:         }
510:         {
511:           const char *extraToken = alias ? tokens[3] : tokens[2];
513:         }
514: destroy:
515:         free(string);
516:         PetscTokenDestroy(&token);
517:         alias = PETSC_FALSE;
518:         line++;
519:       }
520:       err = fclose(fd);
522:       PetscSegBufferGetSize(aseg,&bytes); /* size without null termination */
523:       PetscMPIIntCast(bytes,&acnt);
524:       PetscSegBufferGet(aseg,1,&astring);
525:       astring[0] = 0;
526:       PetscSegBufferGetSize(vseg,&bytes); /* size without null termination */
527:       PetscMPIIntCast(bytes,&cnt);
528:       PetscSegBufferGet(vseg,1,&vstring);
529:       vstring[0] = 0;
530:       PetscMalloc1(2+acnt+cnt,&packed);
531:       PetscSegBufferExtractTo(aseg,packed);
532:       PetscSegBufferExtractTo(vseg,packed+acnt+1);
533:       PetscSegBufferDestroy(&aseg);
534:       PetscSegBufferDestroy(&vseg);
536:   }

538:   counts[0] = acnt;
539:   counts[1] = cnt;
540:   err = MPI_Bcast(counts,2,MPI_INT,0,comm);
542:   acnt = counts[0];
543:   cnt = counts[1];
544:   if (rank) {
545:     PetscMalloc1(2+acnt+cnt,&packed);
546:   }
547:   if (acnt || cnt) {
548:     MPI_Bcast(packed,2+acnt+cnt,MPI_CHAR,0,comm);
549:     astring = packed;
550:     vstring = packed + acnt + 1;
551:   }

553:   if (acnt) {
554:     PetscTokenCreate(astring,' ',&token);
555:     PetscTokenFind(token,&tokens[0]);
556:     while (tokens[0]) {
557:       PetscTokenFind(token,&tokens[1]);
558:       PetscOptionsSetAlias(options,tokens[0],tokens[1]);
559:       PetscTokenFind(token,&tokens[0]);
560:     }
561:     PetscTokenDestroy(&token);
562:   }

564:   if (cnt) {
565:     PetscOptionsInsertString(options,vstring);
566:   }
567:   PetscFree(packed);
568:   return 0;
569: }

571: /*@C
572:      PetscOptionsInsertFile - Inserts options into the database from a file.

574:      Collective

576:   Input Parameters:
577: +   comm - the processes that will share the options (usually PETSC_COMM_WORLD)
578: .   options - options database, use NULL for default global database
579: .   file - name of file,
580:            ".yml" and ".yaml" filename extensions are inserted as YAML options,
581:            append ":yaml" to filename to force YAML options.
582: -   require - if PETSC_TRUE will generate an error if the file does not exist

584:   Notes:
585:    Use  # for lines that are comments and which should be ignored.
586:    Usually, instead of using this command, one should list the file name in the call to PetscInitialize(), this insures that certain options
587:    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
588:    calls to XXXSetFromOptions() it should not be used for options listed under PetscInitialize().
589:    The collectivity of this routine is complex; only the MPI processes in comm will
590:    have the affect of these options. If some processes that create objects call this routine and others do
591:    not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
592:    on different ranks.

594:   Level: developer

596: .seealso: PetscOptionsSetValue(), PetscOptionsView(), PetscOptionsHasName(), PetscOptionsGetInt(),
597:           PetscOptionsGetReal(), PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
598:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
599:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
600:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
601:           PetscOptionsFList(), PetscOptionsEList()

603: @*/
604: PetscErrorCode PetscOptionsInsertFile(MPI_Comm comm,PetscOptions options,const char file[],PetscBool require)
605: {
606:   char           filename[PETSC_MAX_PATH_LEN];
607:   PetscBool      yaml;

609:   PetscOptionsFilename(comm,file,filename,&yaml);
610:   if (yaml) {
611:     PetscOptionsInsertFileYAML(comm,options,filename,require);
612:   } else {
613:     PetscOptionsInsertFilePetsc(comm,options,filename,require);
614:   }
615:   return 0;
616: }

618: /*@C
619:    PetscOptionsInsertArgs - Inserts options into the database from a array of strings

621:    Logically Collective

623:    Input Parameters:
624: +  options - options object
625: .  argc - the array length
626: -  args - the string array

628:    Level: intermediate

630: .seealso: PetscOptions, PetscOptionsInsertString(), PetscOptionsInsertFile()
631: @*/
632: PetscErrorCode PetscOptionsInsertArgs(PetscOptions options,int argc,char *args[])
633: {
634:   MPI_Comm       comm = PETSC_COMM_WORLD;
635:   int            left          = PetscMax(argc,0);
636:   char           *const *eargs = args;

638:   while (left) {
639:     PetscBool isfile,isfileyaml,isstringyaml,ispush,ispop,key;
640:     PetscStrcasecmp(eargs[0],"-options_file",&isfile);
641:     PetscStrcasecmp(eargs[0],"-options_file_yaml",&isfileyaml);
642:     PetscStrcasecmp(eargs[0],"-options_string_yaml",&isstringyaml);
643:     PetscStrcasecmp(eargs[0],"-prefix_push",&ispush);
644:     PetscStrcasecmp(eargs[0],"-prefix_pop",&ispop);
645:     PetscOptionsValidKey(eargs[0],&key);
646:     if (!key) {
647:       eargs++; left--;
648:     } else if (isfile) {
650:       PetscOptionsInsertFile(comm,options,eargs[1],PETSC_TRUE);
651:       eargs += 2; left -= 2;
652:     } else if (isfileyaml) {
654:       PetscOptionsInsertFileYAML(comm,options,eargs[1],PETSC_TRUE);
655:       eargs += 2; left -= 2;
656:     } else if (isstringyaml) {
658:       PetscOptionsInsertStringYAML(options,eargs[1]);
659:       eargs += 2; left -= 2;
660:     } else if (ispush) {
663:       PetscOptionsPrefixPush(options,eargs[1]);
664:       eargs += 2; left -= 2;
665:     } else if (ispop) {
666:       PetscOptionsPrefixPop(options);
667:       eargs++; left--;
668:     } else {
669:       PetscBool nextiskey = PETSC_FALSE;
670:       if (left >= 2) PetscOptionsValidKey(eargs[1],&nextiskey);
671:       if (left < 2 || nextiskey) {
672:         PetscOptionsSetValue(options,eargs[0],NULL);
673:         eargs++; left--;
674:       } else {
675:         PetscOptionsSetValue(options,eargs[0],eargs[1]);
676:         eargs += 2; left -= 2;
677:       }
678:     }
679:   }
680:   return 0;
681: }

683: static inline PetscErrorCode PetscOptionsStringToBoolIfSet_Private(enum PetscPrecedentOption opt,const char *val[],PetscBool set[],PetscBool *flg)
684: {
685:   if (set[opt]) {
686:     PetscOptionsStringToBool(val[opt],flg);
687:   } else *flg = PETSC_FALSE;
688:   return 0;
689: }

691: /* Process options with absolute precedence */
692: static PetscErrorCode PetscOptionsProcessPrecedentFlags(PetscOptions options,int argc,char *args[],PetscBool *skip_petscrc,PetscBool *skip_petscrc_set)
693: {
694:   const char* const *opt = precedentOptions;
695:   const size_t      n = PO_NUM;
696:   size_t            o;
697:   int               a;
698:   const char        **val;
699:   char              **cval;
700:   PetscBool         *set;

702:   PetscCalloc2(n,&cval,n,&set);
703:   val = (const char**) cval;

705:   /* Look for options possibly set using PetscOptionsSetValue beforehand */
706:   for (o=0; o<n; o++) {
707:     PetscOptionsFindPair(options,NULL,opt[o],&val[o],&set[o]);
708:   }

710:   /* Loop through all args to collect last occurring value of each option */
711:   for (a=1; a<argc; a++) {
712:     PetscBool valid, eq;

714:     PetscOptionsValidKey(args[a],&valid);
715:     if (!valid) continue;
716:     for (o=0; o<n; o++) {
717:       PetscStrcasecmp(args[a],opt[o],&eq);
718:       if (eq) {
719:         set[o] = PETSC_TRUE;
720:         if (a == argc-1 || !args[a+1] || !args[a+1][0] || args[a+1][0] == '-') val[o] = NULL;
721:         else val[o] = args[a+1];
722:         break;
723:       }
724:     }
725:   }

727:   /* Process flags */
728:   PetscStrcasecmp(val[PO_HELP], "intro", &options->help_intro);
729:   if (options->help_intro) options->help = PETSC_TRUE;
730:   else PetscOptionsStringToBoolIfSet_Private(PO_HELP,            val,set,&options->help);
731:   PetscOptionsStringToBoolIfSet_Private(PO_OPTIONS_MONITOR_CANCEL,val,set,&options->monitorCancel);
732:   PetscOptionsStringToBoolIfSet_Private(PO_OPTIONS_MONITOR,       val,set,&options->monitorFromOptions);
733:   PetscOptionsStringToBoolIfSet_Private(PO_SKIP_PETSCRC,          val,set,skip_petscrc);
734:   *skip_petscrc_set = set[PO_SKIP_PETSCRC];

736:   /* Store precedent options in database and mark them as used */
737:   for (o=0; o<n; o++) {
738:     if (set[o]) {
739:       PetscOptionsSetValue_Private(options,opt[o],val[o],&a);
740:       options->used[a] = PETSC_TRUE;
741:     }
742:   }
743:   PetscFree2(cval,set);
744:   options->precedentProcessed = PETSC_TRUE;
745:   return 0;
746: }

748: static inline PetscErrorCode PetscOptionsSkipPrecedent(PetscOptions options,const char name[],PetscBool *flg)
749: {
751:   *flg = PETSC_FALSE;
752:   if (options->precedentProcessed) {
753:     for (int i = 0; i < PO_NUM; ++i) {
754:       if (!PetscOptNameCmp(precedentOptions[i],name)) {
755:         /* check if precedent option has been set already */
756:         PetscOptionsFindPair(options,NULL,name,NULL,flg);
757:         if (*flg) break;
758:       }
759:     }
760:   }
761:   return 0;
762: }

764: /*@C
765:    PetscOptionsInsert - Inserts into the options database from the command line,
766:                         the environmental variable and a file.

768:    Collective on PETSC_COMM_WORLD

770:    Input Parameters:
771: +  options - options database or NULL for the default global database
772: .  argc - count of number of command line arguments
773: .  args - the command line arguments
774: -  file - [optional] PETSc database file, append ":yaml" to filename to specify YAML options format.
775:           Use NULL or empty string to not check for code specific file.
776:           Also checks ~/.petscrc, .petscrc and petscrc.
777:           Use -skip_petscrc in the code specific file (or command line) to skip ~/.petscrc, .petscrc and petscrc files.

779:    Note:
780:    Since PetscOptionsInsert() is automatically called by PetscInitialize(),
781:    the user does not typically need to call this routine. PetscOptionsInsert()
782:    can be called several times, adding additional entries into the database.

784:    Options Database Keys:
785: +   -options_file <filename> - read options from a file
786: -   -options_file_yaml <filename> - read options from a YAML file

788:    See PetscInitialize() for options related to option database monitoring.

790:    Level: advanced

792: .seealso: PetscOptionsDestroy(), PetscOptionsView(), PetscOptionsInsertString(), PetscOptionsInsertFile(),
793:           PetscInitialize()
794: @*/
795: PetscErrorCode PetscOptionsInsert(PetscOptions options,int *argc,char ***args,const char file[])
796: {
797:   MPI_Comm       comm = PETSC_COMM_WORLD;
798:   PetscMPIInt    rank;
799:   PetscBool      hasArgs = (argc && *argc) ? PETSC_TRUE : PETSC_FALSE;
800:   PetscBool      skipPetscrc = PETSC_FALSE, skipPetscrcSet = PETSC_FALSE;

803:   MPI_Comm_rank(comm,&rank);

805:   if (!options) {
806:     PetscOptionsCreateDefault();
807:     options = defaultoptions;
808:   }
809:   if (hasArgs) {
810:     /* process options with absolute precedence */
811:     PetscOptionsProcessPrecedentFlags(options,*argc,*args,&skipPetscrc,&skipPetscrcSet);
812:   }
813:   if (file && file[0]) {
814:     PetscOptionsInsertFile(comm,options,file,PETSC_TRUE);
815:     /* if -skip_petscrc has not been set from command line, check whether it has been set in the file */
816:     if (!skipPetscrcSet) PetscOptionsGetBool(options,NULL,"-skip_petscrc",&skipPetscrc,NULL);
817:   }
818:   if (!skipPetscrc) {
819:     char filename[PETSC_MAX_PATH_LEN];
820:     PetscGetHomeDirectory(filename,sizeof(filename));
821:     MPI_Bcast(filename,(int)sizeof(filename),MPI_CHAR,0,comm);
822:     if (filename[0]) PetscStrcat(filename,"/.petscrc");
823:     PetscOptionsInsertFile(comm,options,filename,PETSC_FALSE);
824:     PetscOptionsInsertFile(comm,options,".petscrc",PETSC_FALSE);
825:     PetscOptionsInsertFile(comm,options,"petscrc",PETSC_FALSE);
826:   }

828:   /* insert environment options */
829:   {
830:     char   *eoptions = NULL;
831:     size_t len       = 0;
832:     if (rank == 0) {
833:       eoptions = (char*)getenv("PETSC_OPTIONS");
834:       PetscStrlen(eoptions,&len);
835:     }
836:     MPI_Bcast(&len,1,MPIU_SIZE_T,0,comm);
837:     if (len) {
838:       if (rank) PetscMalloc1(len+1,&eoptions);
839:       MPI_Bcast(eoptions,len,MPI_CHAR,0,comm);
840:       if (rank) eoptions[len] = 0;
841:       PetscOptionsInsertString(options,eoptions);
842:       if (rank) PetscFree(eoptions);
843:     }
844:   }

846:   /* insert YAML environment options */
847:   {
848:     char   *eoptions = NULL;
849:     size_t len       = 0;
850:     if (rank == 0) {
851:       eoptions = (char*)getenv("PETSC_OPTIONS_YAML");
852:       PetscStrlen(eoptions,&len);
853:     }
854:     MPI_Bcast(&len,1,MPIU_SIZE_T,0,comm);
855:     if (len) {
856:       if (rank) PetscMalloc1(len+1,&eoptions);
857:       MPI_Bcast(eoptions,len,MPI_CHAR,0,comm);
858:       if (rank) eoptions[len] = 0;
859:       PetscOptionsInsertStringYAML(options,eoptions);
860:       if (rank) PetscFree(eoptions);
861:     }
862:   }

864:   /* insert command line options here because they take precedence over arguments in petscrc/environment */
865:   if (hasArgs) PetscOptionsInsertArgs(options,*argc-1,*args+1);
866:   return 0;
867: }

869: /*@C
870:    PetscOptionsView - Prints the options that have been loaded. This is
871:    useful for debugging purposes.

873:    Logically Collective on PetscViewer

875:    Input Parameters:
876: +  options - options database, use NULL for default global database
877: -  viewer - must be an PETSCVIEWERASCII viewer

879:    Options Database Key:
880: .  -options_view - Activates PetscOptionsView() within PetscFinalize()

882:    Notes:
883:    Only the rank zero process of MPI_Comm used to create view prints the option values. Other processes
884:    may have different values but they are not printed.

886:    Level: advanced

888: .seealso: PetscOptionsAllUsed()
889: @*/
890: PetscErrorCode PetscOptionsView(PetscOptions options,PetscViewer viewer)
891: {
892:   PetscInt       i;
893:   PetscBool      isascii;

896:   options = options ? options : defaultoptions;
897:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_WORLD;
898:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);

901:   if (!options->N) {
902:     PetscViewerASCIIPrintf(viewer,"#No PETSc Option Table entries\n");
903:     return 0;
904:   }

906:   PetscViewerASCIIPrintf(viewer,"#PETSc Option Table entries:\n");
907:   for (i=0; i<options->N; i++) {
908:     if (options->values[i]) {
909:       PetscViewerASCIIPrintf(viewer,"-%s %s\n",options->names[i],options->values[i]);
910:     } else {
911:       PetscViewerASCIIPrintf(viewer,"-%s\n",options->names[i]);
912:     }
913:   }
914:   PetscViewerASCIIPrintf(viewer,"#End of PETSc Option Table entries\n");
915:   return 0;
916: }

918: /*
919:    Called by error handlers to print options used in run
920: */
921: PETSC_EXTERN PetscErrorCode PetscOptionsViewError(void)
922: {
923:   PetscInt     i;
924:   PetscOptions options = defaultoptions;

926:   if (options->N) {
927:     (*PetscErrorPrintf)("PETSc Option Table entries:\n");
928:   } else {
929:     (*PetscErrorPrintf)("No PETSc Option Table entries\n");
930:   }
931:   for (i=0; i<options->N; i++) {
932:     if (options->values[i]) {
933:       (*PetscErrorPrintf)("-%s %s\n",options->names[i],options->values[i]);
934:     } else {
935:       (*PetscErrorPrintf)("-%s\n",options->names[i]);
936:     }
937:   }
938:   return 0;
939: }

941: /*@C
942:    PetscOptionsPrefixPush - Designate a prefix to be used by all options insertions to follow.

944:    Logically Collective

946:    Input Parameters:
947: +  options - options database, or NULL for the default global database
948: -  prefix - The string to append to the existing prefix

950:    Options Database Keys:
951: +   -prefix_push <some_prefix_> - push the given prefix
952: -   -prefix_pop - pop the last prefix

954:    Notes:
955:    It is common to use this in conjunction with -options_file as in

957: $ -prefix_push system1_ -options_file system1rc -prefix_pop -prefix_push system2_ -options_file system2rc -prefix_pop

959:    where the files no longer require all options to be prefixed with -system2_.

961:    The collectivity of this routine is complex; only the MPI processes that call this routine will
962:    have the affect of these options. If some processes that create objects call this routine and others do
963:    not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
964:    on different ranks.

966: Level: advanced

968: .seealso: PetscOptionsPrefixPop(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsSetValue()
969: @*/
970: PetscErrorCode PetscOptionsPrefixPush(PetscOptions options,const char prefix[])
971: {
972:   size_t         n;
973:   PetscInt       start;
974:   char           key[MAXOPTNAME+1];
975:   PetscBool      valid;

978:   options = options ? options : defaultoptions;
980:   key[0] = '-'; /* keys must start with '-' */
981:   PetscStrncpy(key+1,prefix,sizeof(key)-1);
982:   PetscOptionsValidKey(key,&valid);
983:   if (!valid && options->prefixind > 0 && isdigit((int)prefix[0])) valid = PETSC_TRUE; /* If the prefix stack is not empty, make numbers a valid prefix */
985:   start = options->prefixind ? options->prefixstack[options->prefixind-1] : 0;
986:   PetscStrlen(prefix,&n);
988:   PetscArraycpy(options->prefix+start,prefix,n+1);
989:   options->prefixstack[options->prefixind++] = start+n;
990:   return 0;
991: }

993: /*@C
994:    PetscOptionsPrefixPop - Remove the latest options prefix, see PetscOptionsPrefixPush() for details

996:    Logically Collective on the MPI_Comm that called PetscOptionsPrefixPush()

998:   Input Parameters:
999: .  options - options database, or NULL for the default global database

1001:    Level: advanced

1003: .seealso: PetscOptionsPrefixPush(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsSetValue()
1004: @*/
1005: PetscErrorCode PetscOptionsPrefixPop(PetscOptions options)
1006: {
1007:   PetscInt offset;

1009:   options = options ? options : defaultoptions;
1011:   options->prefixind--;
1012:   offset = options->prefixind ? options->prefixstack[options->prefixind-1] : 0;
1013:   options->prefix[offset] = 0;
1014:   return 0;
1015: }

1017: /*@C
1018:     PetscOptionsClear - Removes all options form the database leaving it empty.

1020:     Logically Collective

1022:   Input Parameters:
1023: .  options - options database, use NULL for the default global database

1025:    The collectivity of this routine is complex; only the MPI processes that call this routine will
1026:    have the affect of these options. If some processes that create objects call this routine and others do
1027:    not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1028:    on different ranks.

1030:    Level: developer

1032: .seealso: PetscOptionsInsert()
1033: @*/
1034: PetscErrorCode PetscOptionsClear(PetscOptions options)
1035: {
1036:   PetscInt i;

1038:   options = options ? options : defaultoptions;
1039:   if (!options) return 0;

1041:   for (i=0; i<options->N; i++) {
1042:     if (options->names[i])  free(options->names[i]);
1043:     if (options->values[i]) free(options->values[i]);
1044:   }
1045:   options->N = 0;

1047:   for (i=0; i<options->Naliases; i++) {
1048:     free(options->aliases1[i]);
1049:     free(options->aliases2[i]);
1050:   }
1051:   options->Naliases = 0;

1053:   /* destroy hash table */
1054:   kh_destroy(HO,options->ht);
1055:   options->ht = NULL;

1057:   options->prefixind = 0;
1058:   options->prefix[0] = 0;
1059:   options->help      = PETSC_FALSE;
1060:   return 0;
1061: }

1063: /*@C
1064:    PetscOptionsSetAlias - Makes a key and alias for another key

1066:    Logically Collective

1068:    Input Parameters:
1069: +  options - options database, or NULL for default global database
1070: .  newname - the alias
1071: -  oldname - the name that alias will refer to

1073:    Level: advanced

1075:    The collectivity of this routine is complex; only the MPI processes that call this routine will
1076:    have the affect of these options. If some processes that create objects call this routine and others do
1077:    not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1078:    on different ranks.

1080: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),OptionsHasName(),
1081:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(),PetscOptionsBool(),
1082:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
1083:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
1084:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
1085:           PetscOptionsFList(), PetscOptionsEList()
1086: @*/
1087: PetscErrorCode PetscOptionsSetAlias(PetscOptions options,const char newname[],const char oldname[])
1088: {
1089:   PetscInt       n;
1090:   size_t         len;
1091:   PetscBool      valid;

1095:   options = options ? options : defaultoptions;
1096:   PetscOptionsValidKey(newname,&valid);
1098:   PetscOptionsValidKey(oldname,&valid);

1101:   n = options->Naliases;

1104:   newname++; oldname++;
1105:   PetscStrlen(newname,&len);
1106:   options->aliases1[n] = (char*)malloc((len+1)*sizeof(char));
1107:   PetscStrcpy(options->aliases1[n],newname);
1108:   PetscStrlen(oldname,&len);
1109:   options->aliases2[n] = (char*)malloc((len+1)*sizeof(char));
1110:   PetscStrcpy(options->aliases2[n],oldname);
1111:   options->Naliases++;
1112:   return 0;
1113: }

1115: /*@C
1116:    PetscOptionsSetValue - Sets an option name-value pair in the options
1117:    database, overriding whatever is already present.

1119:    Logically Collective

1121:    Input Parameters:
1122: +  options - options database, use NULL for the default global database
1123: .  name - name of option, this SHOULD have the - prepended
1124: -  value - the option value (not used for all options, so can be NULL)

1126:    Level: intermediate

1128:    Note:
1129:    This function can be called BEFORE PetscInitialize()

1131:    The collectivity of this routine is complex; only the MPI processes that call this routine will
1132:    have the affect of these options. If some processes that create objects call this routine and others do
1133:    not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options
1134:    on different ranks.

1136:    Developers Note: Uses malloc() directly because PETSc may not be initialized yet.

1138: .seealso: PetscOptionsInsert(), PetscOptionsClearValue()
1139: @*/
1140: PetscErrorCode PetscOptionsSetValue(PetscOptions options,const char name[],const char value[])
1141: {
1142:   PetscOptionsSetValue_Private(options,name,value,NULL);
1143:   return 0;
1144: }

1146: static PetscErrorCode PetscOptionsSetValue_Private(PetscOptions options,const char name[],const char value[],int *pos)
1147: {
1148:   size_t      len;
1149:   int         N,n,i;
1150:   char      **names;
1151:   char        fullname[MAXOPTNAME] = "";
1152:   PetscBool   flg;

1154:   if (!options) {
1155:     PetscOptionsCreateDefault();
1156:     options = defaultoptions;
1157:   }

1160:   PetscOptionsSkipPrecedent(options,name,&flg);
1161:   if (flg) return 0;

1163:   name++; /* skip starting dash */

1165:   if (options->prefixind > 0) {
1166:     strncpy(fullname,options->prefix,sizeof(fullname));
1167:     fullname[sizeof(fullname)-1] = 0;
1168:     strncat(fullname,name,sizeof(fullname)-strlen(fullname)-1);
1169:     fullname[sizeof(fullname)-1] = 0;
1170:     name = fullname;
1171:   }

1173:   /* check against aliases */
1174:   N = options->Naliases;
1175:   for (i=0; i<N; i++) {
1176:     int result = PetscOptNameCmp(options->aliases1[i],name);
1177:     if (!result) { name = options->aliases2[i]; break; }
1178:   }

1180:   /* slow search */
1181:   N = n = options->N;
1182:   names = options->names;
1183:   for (i=0; i<N; i++) {
1184:     int result = PetscOptNameCmp(names[i],name);
1185:     if (!result) {
1186:       n = i; goto setvalue;
1187:     } else if (result > 0) {
1188:       n = i; break;
1189:     }
1190:   }

1193:   /* shift remaining values up 1 */
1194:   for (i=N; i>n; i--) {
1195:     options->names[i]  = options->names[i-1];
1196:     options->values[i] = options->values[i-1];
1197:     options->used[i]   = options->used[i-1];
1198:   }
1199:   options->names[n]  = NULL;
1200:   options->values[n] = NULL;
1201:   options->used[n]   = PETSC_FALSE;
1202:   options->N++;

1204:   /* destroy hash table */
1205:   kh_destroy(HO,options->ht);
1206:   options->ht = NULL;

1208:   /* set new name */
1209:   len = strlen(name);
1210:   options->names[n] = (char*)malloc((len+1)*sizeof(char));
1212:   strcpy(options->names[n],name);

1214: setvalue:
1215:   /* set new value */
1216:   if (options->values[n]) free(options->values[n]);
1217:   len = value ? strlen(value) : 0;
1218:   if (len) {
1219:     options->values[n] = (char*)malloc((len+1)*sizeof(char));
1220:     if (!options->values[n]) return PETSC_ERR_MEM;
1221:     strcpy(options->values[n],value);
1222:   } else {
1223:     options->values[n] = NULL;
1224:   }

1226:   /* handle -help so that it can be set from anywhere */
1227:   if (!PetscOptNameCmp(name,"help")) {
1228:     options->help = PETSC_TRUE;
1229:     options->help_intro = (value && !PetscOptNameCmp(value,"intro")) ? PETSC_TRUE : PETSC_FALSE;
1230:     options->used[n] = PETSC_TRUE;
1231:   }

1233:   PetscOptionsMonitor(options,name,value);
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 Parameters:
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;

1262:   options = options ? options : defaultoptions;
1264:   if (!PetscOptNameCmp(name,"-help")) options->help = options->help_intro = PETSC_FALSE;

1266:   name++; /* skip starting dash */

1268:   /* slow search */
1269:   N = n = options->N;
1270:   names = options->names;
1271:   for (i=0; i<N; i++) {
1272:     int result = PetscOptNameCmp(names[i],name);
1273:     if (!result) {
1274:       n = i; break;
1275:     } else if (result > 0) {
1276:       n = N; break;
1277:     }
1278:   }
1279:   if (n == N) return 0; /* it was not present */

1281:   /* remove name and value */
1282:   if (options->names[n])  free(options->names[n]);
1283:   if (options->values[n]) free(options->values[n]);
1284:   /* shift remaining values down 1 */
1285:   for (i=n; i<N-1; i++) {
1286:     options->names[i]  = options->names[i+1];
1287:     options->values[i] = options->values[i+1];
1288:     options->used[i]   = options->used[i+1];
1289:   }
1290:   options->N--;

1292:   /* destroy hash table */
1293:   kh_destroy(HO,options->ht);
1294:   options->ht = NULL;

1296:   PetscOptionsMonitor(options,name,NULL);
1297:   return 0;
1298: }

1300: /*@C
1301:    PetscOptionsFindPair - Gets an option name-value pair from the options database.

1303:    Not Collective

1305:    Input Parameters:
1306: +  options - options database, use NULL for the default global database
1307: .  pre - the string to prepend to the name or NULL, this SHOULD NOT have the "-" prepended
1308: -  name - name of option, this SHOULD have the "-" prepended

1310:    Output Parameters:
1311: +  value - the option value (optional, not used for all options)
1312: -  set - whether the option is set (optional)

1314:    Notes:
1315:    Each process may find different values or no value depending on how options were inserted into the database

1317:    Level: developer

1319: .seealso: PetscOptionsSetValue(), PetscOptionsClearValue()
1320: @*/
1321: PetscErrorCode PetscOptionsFindPair(PetscOptions options,const char pre[],const char name[],const char *value[],PetscBool *set)
1322: {
1323:   char           buf[MAXOPTNAME];
1324:   PetscBool      usehashtable = PETSC_TRUE;
1325:   PetscBool      matchnumbers = PETSC_TRUE;

1327:   options = options ? options : defaultoptions;

1331:   name++; /* skip starting dash */

1333:   /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */
1334:   if (pre && pre[0]) {
1335:     char *ptr = buf;
1336:     if (name[0] == '-') { *ptr++ = '-';  name++; }
1337:     PetscStrncpy(ptr,pre,buf+sizeof(buf)-ptr);
1338:     PetscStrlcat(buf,name,sizeof(buf));
1339:     name = buf;
1340:   }

1342:   if (PetscDefined(USE_DEBUG)) {
1343:     PetscBool valid;
1344:     char      key[MAXOPTNAME+1] = "-";
1345:     PetscStrncpy(key+1,name,sizeof(key)-1);
1346:     PetscOptionsValidKey(key,&valid);
1348:   }

1350:   if (!options->ht && usehashtable) {
1351:     int i,ret;
1352:     khiter_t it;
1353:     khash_t(HO) *ht;
1354:     ht = kh_init(HO);
1356:     ret = kh_resize(HO,ht,options->N*2); /* twice the required size to reduce risk of collisions */
1358:     for (i=0; i<options->N; i++) {
1359:       it = kh_put(HO,ht,options->names[i],&ret);
1361:       kh_val(ht,it) = i;
1362:     }
1363:     options->ht = ht;
1364:   }

1366:   if (usehashtable)
1367:   { /* fast search */
1368:     khash_t(HO) *ht = options->ht;
1369:     khiter_t it = kh_get(HO,ht,name);
1370:     if (it != kh_end(ht)) {
1371:       int i = kh_val(ht,it);
1372:       options->used[i]  = PETSC_TRUE;
1373:       if (value) *value = options->values[i];
1374:       if (set)   *set   = PETSC_TRUE;
1375:       return 0;
1376:     }
1377:   } else
1378:   { /* slow search */
1379:     int i, N = options->N;
1380:     for (i=0; i<N; i++) {
1381:       int result = PetscOptNameCmp(options->names[i],name);
1382:       if (!result) {
1383:         options->used[i]  = PETSC_TRUE;
1384:         if (value) *value = options->values[i];
1385:         if (set)   *set   = PETSC_TRUE;
1386:         return 0;
1387:       } else if (result > 0) {
1388:         break;
1389:       }
1390:     }
1391:   }

1393:   /*
1394:    The following block slows down all lookups in the most frequent path (most lookups are unsuccessful).
1395:    Maybe this special lookup mode should be enabled on request with a push/pop API.
1396:    The feature of matching _%d_ used sparingly in the codebase.
1397:    */
1398:   if (matchnumbers) {
1399:     int i,j,cnt = 0,locs[16],loce[16];
1400:     /* determine the location and number of all _%d_ in the key */
1401:     for (i=0; name[i]; i++) {
1402:       if (name[i] == '_') {
1403:         for (j=i+1; name[j]; j++) {
1404:           if (name[j] >= '0' && name[j] <= '9') continue;
1405:           if (name[j] == '_' && j > i+1) { /* found a number */
1406:             locs[cnt]   = i+1;
1407:             loce[cnt++] = j+1;
1408:           }
1409:           i = j-1;
1410:           break;
1411:         }
1412:       }
1413:     }
1414:     for (i=0; i<cnt; i++) {
1415:       PetscBool found;
1416:       char      opt[MAXOPTNAME+1] = "-", tmp[MAXOPTNAME];
1417:       PetscStrncpy(tmp,name,PetscMin((size_t)(locs[i]+1),sizeof(tmp)));
1418:       PetscStrlcat(opt,tmp,sizeof(opt));
1419:       PetscStrlcat(opt,name+loce[i],sizeof(opt));
1420:       PetscOptionsFindPair(options,NULL,opt,value,&found);
1421:       if (found) {if (set) *set = PETSC_TRUE; return 0;}
1422:     }
1423:   }

1425:   if (set) *set = PETSC_FALSE;
1426:   return 0;
1427: }

1429: /* Check whether any option begins with pre+name */
1430: PETSC_EXTERN PetscErrorCode PetscOptionsFindPairPrefix_Private(PetscOptions options,const char pre[], const char name[],const char *value[],PetscBool *set)
1431: {
1432:   char           buf[MAXOPTNAME];
1433:   int            numCnt = 0, locs[16],loce[16];

1435:   options = options ? options : defaultoptions;

1439:   name++; /* skip starting dash */

1441:   /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */
1442:   if (pre && pre[0]) {
1443:     char *ptr = buf;
1444:     if (name[0] == '-') { *ptr++ = '-';  name++; }
1445:     PetscStrncpy(ptr,pre,sizeof(buf)+(size_t)(ptr-buf));
1446:     PetscStrlcat(buf,name,sizeof(buf));
1447:     name = buf;
1448:   }

1450:   if (PetscDefined(USE_DEBUG)) {
1451:     PetscBool valid;
1452:     char      key[MAXOPTNAME+1] = "-";
1453:     PetscStrncpy(key+1,name,sizeof(key)-1);
1454:     PetscOptionsValidKey(key,&valid);
1456:   }

1458:   /* determine the location and number of all _%d_ in the key */
1459:   {
1460:     int i,j;
1461:     for (i=0; name[i]; i++) {
1462:       if (name[i] == '_') {
1463:         for (j=i+1; name[j]; j++) {
1464:           if (name[j] >= '0' && name[j] <= '9') continue;
1465:           if (name[j] == '_' && j > i+1) { /* found a number */
1466:             locs[numCnt]   = i+1;
1467:             loce[numCnt++] = j+1;
1468:           }
1469:           i = j-1;
1470:           break;
1471:         }
1472:       }
1473:     }
1474:   }

1476:   { /* slow search */
1477:     int       c, i;
1478:     size_t    len;
1479:     PetscBool match;

1481:     for (c = -1; c < numCnt; ++c) {
1482:       char opt[MAXOPTNAME+1] = "", tmp[MAXOPTNAME];

1484:       if (c < 0) {
1485:         PetscStrcpy(opt,name);
1486:       } else {
1487:         PetscStrncpy(tmp,name,PetscMin((size_t)(locs[c]+1),sizeof(tmp)));
1488:         PetscStrlcat(opt,tmp,sizeof(opt));
1489:         PetscStrlcat(opt,name+loce[c],sizeof(opt));
1490:       }
1491:       PetscStrlen(opt,&len);
1492:       for (i=0; i<options->N; i++) {
1493:         PetscStrncmp(options->names[i],opt,len,&match);
1494:         if (match) {
1495:           options->used[i]  = PETSC_TRUE;
1496:           if (value) *value = options->values[i];
1497:           if (set)   *set   = PETSC_TRUE;
1498:           return 0;
1499:         }
1500:       }
1501:     }
1502:   }

1504:   if (set) *set = PETSC_FALSE;
1505:   return 0;
1506: }

1508: /*@C
1509:    PetscOptionsReject - Generates an error if a certain option is given.

1511:    Not Collective

1513:    Input Parameters:
1514: +  options - options database, use NULL for default global database
1515: .  pre - the option prefix (may be NULL)
1516: .  name - the option name one is seeking
1517: -  mess - error message (may be NULL)

1519:    Level: advanced

1521: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),OptionsHasName(),
1522:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
1523:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
1524:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
1525:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
1526:           PetscOptionsFList(), PetscOptionsEList()
1527: @*/
1528: PetscErrorCode PetscOptionsReject(PetscOptions options,const char pre[],const char name[],const char mess[])
1529: {
1530:   PetscBool      flag = PETSC_FALSE;

1532:   PetscOptionsHasName(options,pre,name,&flag);
1533:   if (flag) {
1535:     else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Program has disabled option: -%s%s",pre?pre:"",name+1);
1536:   }
1537:   return 0;
1538: }

1540: /*@C
1541:    PetscOptionsHasHelp - Determines whether the "-help" option is in the database.

1543:    Not Collective

1545:    Input Parameters:
1546: .  options - options database, use NULL for default global database

1548:    Output Parameters:
1549: .  set - PETSC_TRUE if found else PETSC_FALSE.

1551:    Level: advanced

1553: .seealso: PetscOptionsHasName()
1554: @*/
1555: PetscErrorCode PetscOptionsHasHelp(PetscOptions options,PetscBool *set)
1556: {
1558:   options = options ? options : defaultoptions;
1559:   *set = options->help;
1560:   return 0;
1561: }

1563: PetscErrorCode PetscOptionsHasHelpIntro_Internal(PetscOptions options,PetscBool *set)
1564: {
1566:   options = options ? options : defaultoptions;
1567:   *set = options->help_intro;
1568:   return 0;
1569: }

1571: /*@C
1572:    PetscOptionsHasName - Determines whether a certain option is given in the database. This returns true whether the option is a number, string or Boolean, even
1573:                       if its value is set to false.

1575:    Not Collective

1577:    Input Parameters:
1578: +  options - options database, use NULL for default global database
1579: .  pre - string to prepend to the name or NULL
1580: -  name - the option one is seeking

1582:    Output Parameters:
1583: .  set - PETSC_TRUE if found else PETSC_FALSE.

1585:    Level: beginner

1587:    Notes:
1588:    In many cases you probably want to use PetscOptionsGetBool() instead of calling this, to allowing toggling values.

1590: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
1591:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
1592:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
1593:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
1594:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
1595:           PetscOptionsFList(), PetscOptionsEList()
1596: @*/
1597: PetscErrorCode PetscOptionsHasName(PetscOptions options,const char pre[],const char name[],PetscBool *set)
1598: {
1599:   const char     *value;
1600:   PetscBool      flag;

1602:   PetscOptionsFindPair(options,pre,name,&value,&flag);
1603:   if (set) *set = flag;
1604:   return 0;
1605: }

1607: /*@C
1608:    PetscOptionsGetAll - Lists all the options the program was run with in a single string.

1610:    Not Collective

1612:    Input Parameter:
1613: .  options - the options database, use NULL for the default global database

1615:    Output Parameter:
1616: .  copts - pointer where string pointer is stored

1618:    Notes:
1619:     The array and each entry in the array should be freed with PetscFree()
1620:     Each process may have different values depending on how the options were inserted into the database

1622:    Level: advanced

1624: .seealso: PetscOptionsAllUsed(), PetscOptionsView(), PetscOptionsPush(), PetscOptionsPop()
1625: @*/
1626: PetscErrorCode PetscOptionsGetAll(PetscOptions options,char *copts[])
1627: {
1628:   PetscInt       i;
1629:   size_t         len = 1,lent = 0;
1630:   char           *coptions = NULL;

1633:   options = options ? options : defaultoptions;
1634:   /* count the length of the required string */
1635:   for (i=0; i<options->N; i++) {
1636:     PetscStrlen(options->names[i],&lent);
1637:     len += 2 + lent;
1638:     if (options->values[i]) {
1639:       PetscStrlen(options->values[i],&lent);
1640:       len += 1 + lent;
1641:     }
1642:   }
1643:   PetscMalloc1(len,&coptions);
1644:   coptions[0] = 0;
1645:   for (i=0; i<options->N; i++) {
1646:     PetscStrcat(coptions,"-");
1647:     PetscStrcat(coptions,options->names[i]);
1648:     PetscStrcat(coptions," ");
1649:     if (options->values[i]) {
1650:       PetscStrcat(coptions,options->values[i]);
1651:       PetscStrcat(coptions," ");
1652:     }
1653:   }
1654:   *copts = coptions;
1655:   return 0;
1656: }

1658: /*@C
1659:    PetscOptionsUsed - Indicates if PETSc has used a particular option set in the database

1661:    Not Collective

1663:    Input Parameters:
1664: +  options - options database, use NULL for default global database
1665: -  name - string name of option

1667:    Output Parameter:
1668: .  used - PETSC_TRUE if the option was used, otherwise false, including if option was not found in options database

1670:    Level: advanced

1672:    Notes:
1673:    The value returned may be different on each process and depends on which options have been processed
1674:    on the given process

1676: .seealso: PetscOptionsView(), PetscOptionsLeft(), PetscOptionsAllUsed()
1677: @*/
1678: PetscErrorCode PetscOptionsUsed(PetscOptions options,const char *name,PetscBool *used)
1679: {
1680:   PetscInt       i;

1684:   options = options ? options : defaultoptions;
1685:   *used = PETSC_FALSE;
1686:   for (i=0; i<options->N; i++) {
1687:     PetscStrcasecmp(options->names[i],name,used);
1688:     if (*used) {
1689:       *used = options->used[i];
1690:       break;
1691:     }
1692:   }
1693:   return 0;
1694: }

1696: /*@
1697:    PetscOptionsAllUsed - Returns a count of the number of options in the
1698:    database that have never been selected.

1700:    Not Collective

1702:    Input Parameter:
1703: .  options - options database, use NULL for default global database

1705:    Output Parameter:
1706: .  N - count of options not used

1708:    Level: advanced

1710:    Notes:
1711:    The value returned may be different on each process and depends on which options have been processed
1712:    on the given process

1714: .seealso: PetscOptionsView()
1715: @*/
1716: PetscErrorCode PetscOptionsAllUsed(PetscOptions options,PetscInt *N)
1717: {
1718:   PetscInt     i,n = 0;

1721:   options = options ? options : defaultoptions;
1722:   for (i=0; i<options->N; i++) {
1723:     if (!options->used[i]) n++;
1724:   }
1725:   *N = n;
1726:   return 0;
1727: }

1729: /*@
1730:    PetscOptionsLeft - Prints to screen any options that were set and never used.

1732:    Not Collective

1734:    Input Parameter:
1735: .  options - options database; use NULL for default global database

1737:    Options Database Key:
1738: .  -options_left - activates PetscOptionsAllUsed() within PetscFinalize()

1740:    Notes:
1741:       This is rarely used directly, it is called by PetscFinalize() in debug more or if -options_left
1742:       is passed otherwise to help users determine possible mistakes in their usage of options. This
1743:       only prints values on process zero of PETSC_COMM_WORLD. Other processes depending the objects
1744:       used may have different options that are left unused.

1746:    Level: advanced

1748: .seealso: PetscOptionsAllUsed()
1749: @*/
1750: PetscErrorCode PetscOptionsLeft(PetscOptions options)
1751: {
1752:   PetscInt       i;
1753:   PetscInt       cnt = 0;
1754:   PetscOptions   toptions;

1756:   toptions = options ? options : defaultoptions;
1757:   for (i=0; i<toptions->N; i++) {
1758:     if (!toptions->used[i]) {
1759:       if (toptions->values[i]) {
1760:         PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s value: %s\n",toptions->names[i],toptions->values[i]);
1761:       } else {
1762:         PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s (no value)\n",toptions->names[i]);
1763:       }
1764:     }
1765:   }
1766:   if (!options) {
1767:     toptions = defaultoptions;
1768:     while (toptions->previous) {
1769:       cnt++;
1770:       toptions = toptions->previous;
1771:     }
1772:     if (cnt) {
1773:       PetscPrintf(PETSC_COMM_WORLD,"Option left: You may have forgotten some calls to PetscOptionsPop(),\n             PetscOptionsPop() has been called %" PetscInt_FMT " less times than PetscOptionsPush()\n",cnt);
1774:     }
1775:   }
1776:   return 0;
1777: }

1779: /*@C
1780:    PetscOptionsLeftGet - Returns all options that were set and never used.

1782:    Not Collective

1784:    Input Parameter:
1785: .  options - options database, use NULL for default global database

1787:    Output Parameters:
1788: +  N - count of options not used
1789: .  names - names of options not used
1790: -  values - values of options not used

1792:    Level: advanced

1794:    Notes:
1795:    Users should call PetscOptionsLeftRestore() to free the memory allocated in this routine
1796:    Notes: The value returned may be different on each process and depends on which options have been processed
1797:    on the given process

1799: .seealso: PetscOptionsAllUsed(), PetscOptionsLeft()
1800: @*/
1801: PetscErrorCode PetscOptionsLeftGet(PetscOptions options,PetscInt *N,char **names[],char **values[])
1802: {
1803:   PetscInt       i,n;

1808:   options = options ? options : defaultoptions;

1810:   /* The number of unused PETSc options */
1811:   n = 0;
1812:   for (i=0; i<options->N; i++) {
1813:     if (!options->used[i]) n++;
1814:   }
1815:   if (N) { *N = n; }
1816:   if (names)  PetscMalloc1(n,names);
1817:   if (values) PetscMalloc1(n,values);

1819:   n = 0;
1820:   if (names || values) {
1821:     for (i=0; i<options->N; i++) {
1822:       if (!options->used[i]) {
1823:         if (names)  (*names)[n]  = options->names[i];
1824:         if (values) (*values)[n] = options->values[i];
1825:         n++;
1826:       }
1827:     }
1828:   }
1829:   return 0;
1830: }

1832: /*@C
1833:    PetscOptionsLeftRestore - Free memory for the unused PETSc options obtained using PetscOptionsLeftGet.

1835:    Not Collective

1837:    Input Parameters:
1838: +  options - options database, use NULL for default global database
1839: .  names - names of options not used
1840: -  values - values of options not used

1842:    Level: advanced

1844: .seealso: PetscOptionsAllUsed(), PetscOptionsLeft(), PetscOptionsLeftGet()
1845: @*/
1846: PetscErrorCode PetscOptionsLeftRestore(PetscOptions options,PetscInt *N,char **names[],char **values[])
1847: {
1851:   if (N) { *N = 0; }
1852:   if (names)  PetscFree(*names);
1853:   if (values) PetscFree(*values);
1854:   return 0;
1855: }

1857: /*@C
1858:    PetscOptionsMonitorDefault - Print all options set value events using the supplied PetscViewer.

1860:    Logically Collective on ctx

1862:    Input Parameters:
1863: +  name  - option name string
1864: .  value - option value string
1865: -  ctx - an ASCII viewer or NULL

1867:    Level: intermediate

1869:    Notes:
1870:      If ctx=NULL, PetscPrintf() is used.
1871:      The first MPI rank in the PetscViewer viewer actually prints the values, other
1872:      processes may have different values set

1874: .seealso: PetscOptionsMonitorSet()
1875: @*/
1876: PetscErrorCode PetscOptionsMonitorDefault(const char name[],const char value[],void *ctx)
1877: {
1878:   if (ctx) {
1879:     PetscViewer viewer = (PetscViewer)ctx;
1880:     if (!value) {
1881:       PetscViewerASCIIPrintf(viewer,"Removing option: %s\n",name);
1882:     } else if (!value[0]) {
1883:       PetscViewerASCIIPrintf(viewer,"Setting option: %s (no value)\n",name);
1884:     } else {
1885:       PetscViewerASCIIPrintf(viewer,"Setting option: %s = %s\n",name,value);
1886:     }
1887:   } else {
1888:     MPI_Comm comm = PETSC_COMM_WORLD;
1889:     if (!value) {
1890:       PetscPrintf(comm,"Removing option: %s\n",name);
1891:     } else if (!value[0]) {
1892:       PetscPrintf(comm,"Setting option: %s (no value)\n",name);
1893:     } else {
1894:       PetscPrintf(comm,"Setting option: %s = %s\n",name,value);
1895:     }
1896:   }
1897:   return 0;
1898: }

1900: /*@C
1901:    PetscOptionsMonitorSet - Sets an ADDITIONAL function to be called at every method that
1902:    modified the PETSc options database.

1904:    Not Collective

1906:    Input Parameters:
1907: +  monitor - pointer to function (if this is NULL, it turns off monitoring
1908: .  mctx    - [optional] context for private data for the
1909:              monitor routine (use NULL if no context is desired)
1910: -  monitordestroy - [optional] routine that frees monitor context
1911:           (may be NULL)

1913:    Calling Sequence of monitor:
1914: $     monitor (const char name[], const char value[], void *mctx)

1916: +  name - option name string
1917: .  value - option value string
1918: -  mctx  - optional monitoring context, as set by PetscOptionsMonitorSet()

1920:    Options Database Keys:
1921:    See PetscInitialize() for options related to option database monitoring.

1923:    Notes:
1924:    The default is to do nothing.  To print the name and value of options
1925:    being inserted into the database, use PetscOptionsMonitorDefault() as the monitoring routine,
1926:    with a null monitoring context.

1928:    Several different monitoring routines may be set by calling
1929:    PetscOptionsMonitorSet() multiple times; all will be called in the
1930:    order in which they were set.

1932:    Level: intermediate

1934: .seealso: PetscOptionsMonitorDefault(), PetscInitialize()
1935: @*/
1936: PetscErrorCode PetscOptionsMonitorSet(PetscErrorCode (*monitor)(const char name[], const char value[], void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
1937: {
1938:   PetscOptions options = defaultoptions;

1940:   if (options->monitorCancel) return 0;
1942:   options->monitor[options->numbermonitors]          = monitor;
1943:   options->monitordestroy[options->numbermonitors]   = monitordestroy;
1944:   options->monitorcontext[options->numbermonitors++] = (void*)mctx;
1945:   return 0;
1946: }

1948: /*
1949:    PetscOptionsStringToBool - Converts string to PetscBool, handles cases like "yes", "no", "true", "false", "0", "1", "off", "on".
1950:      Empty string is considered as true.
1951: */
1952: PetscErrorCode PetscOptionsStringToBool(const char value[],PetscBool *a)
1953: {
1954:   PetscBool istrue,isfalse;
1955:   size_t    len;

1957:   /* PetscStrlen() returns 0 for NULL or "" */
1958:   PetscStrlen(value,&len);
1959:   if (!len)  {*a = PETSC_TRUE; return 0;}
1960:   PetscStrcasecmp(value,"TRUE",&istrue);
1961:   if (istrue) {*a = PETSC_TRUE; return 0;}
1962:   PetscStrcasecmp(value,"YES",&istrue);
1963:   if (istrue) {*a = PETSC_TRUE; return 0;}
1964:   PetscStrcasecmp(value,"1",&istrue);
1965:   if (istrue) {*a = PETSC_TRUE; return 0;}
1966:   PetscStrcasecmp(value,"on",&istrue);
1967:   if (istrue) {*a = PETSC_TRUE; return 0;}
1968:   PetscStrcasecmp(value,"FALSE",&isfalse);
1969:   if (isfalse) {*a = PETSC_FALSE; return 0;}
1970:   PetscStrcasecmp(value,"NO",&isfalse);
1971:   if (isfalse) {*a = PETSC_FALSE; return 0;}
1972:   PetscStrcasecmp(value,"0",&isfalse);
1973:   if (isfalse) {*a = PETSC_FALSE; return 0;}
1974:   PetscStrcasecmp(value,"off",&isfalse);
1975:   if (isfalse) {*a = PETSC_FALSE; return 0;}
1976:   SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unknown logical value: %s",value);
1977: }

1979: /*
1980:    PetscOptionsStringToInt - Converts a string to an integer value. Handles special cases such as "default" and "decide"
1981: */
1982: PetscErrorCode PetscOptionsStringToInt(const char name[],PetscInt *a)
1983: {
1984:   size_t    len;
1985:   PetscBool decide,tdefault,mouse;

1987:   PetscStrlen(name,&len);

1990:   PetscStrcasecmp(name,"PETSC_DEFAULT",&tdefault);
1991:   if (!tdefault) {
1992:     PetscStrcasecmp(name,"DEFAULT",&tdefault);
1993:   }
1994:   PetscStrcasecmp(name,"PETSC_DECIDE",&decide);
1995:   if (!decide) {
1996:     PetscStrcasecmp(name,"DECIDE",&decide);
1997:   }
1998:   PetscStrcasecmp(name,"mouse",&mouse);

2000:   if (tdefault)    *a = PETSC_DEFAULT;
2001:   else if (decide) *a = PETSC_DECIDE;
2002:   else if (mouse)  *a = -1;
2003:   else {
2004:     char *endptr;
2005:     long strtolval;

2007:     strtolval = strtol(name,&endptr,10);

2010: #if defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE_ATOLL)
2011:     (void) strtolval;
2012:     *a = atoll(name);
2013: #elif defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE___INT64)
2014:     (void) strtolval;
2015:     *a = _atoi64(name);
2016: #else
2017:     *a = (PetscInt)strtolval;
2018: #endif
2019:   }
2020:   return 0;
2021: }

2023: #if defined(PETSC_USE_REAL___FLOAT128)
2024: #include <quadmath.h>
2025: #endif

2027: static PetscErrorCode PetscStrtod(const char name[],PetscReal *a,char **endptr)
2028: {
2029: #if defined(PETSC_USE_REAL___FLOAT128)
2030:   *a = strtoflt128(name,endptr);
2031: #else
2032:   *a = (PetscReal)strtod(name,endptr);
2033: #endif
2034:   return 0;
2035: }

2037: static PetscErrorCode PetscStrtoz(const char name[],PetscScalar *a,char **endptr,PetscBool *isImaginary)
2038: {
2039:   PetscBool      hasi = PETSC_FALSE;
2040:   char           *ptr;
2041:   PetscReal      strtoval;

2043:   PetscStrtod(name,&strtoval,&ptr);
2044:   if (ptr == name) {
2045:     strtoval = 1.;
2046:     hasi = PETSC_TRUE;
2047:     if (name[0] == 'i') {
2048:       ptr++;
2049:     } else if (name[0] == '+' && name[1] == 'i') {
2050:       ptr += 2;
2051:     } else if (name[0] == '-' && name[1] == 'i') {
2052:       strtoval = -1.;
2053:       ptr += 2;
2054:     }
2055:   } else if (*ptr == 'i') {
2056:     hasi = PETSC_TRUE;
2057:     ptr++;
2058:   }
2059:   *endptr = ptr;
2060:   *isImaginary = hasi;
2061:   if (hasi) {
2062: #if !defined(PETSC_USE_COMPLEX)
2063:     SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s contains imaginary but complex not supported ",name);
2064: #else
2065:     *a = PetscCMPLX(0.,strtoval);
2066: #endif
2067:   } else {
2068:     *a = strtoval;
2069:   }
2070:   return 0;
2071: }

2073: /*
2074:    Converts a string to PetscReal value. Handles special cases like "default" and "decide"
2075: */
2076: PetscErrorCode PetscOptionsStringToReal(const char name[],PetscReal *a)
2077: {
2078:   size_t     len;
2079:   PetscBool  match;
2080:   char      *endptr;

2082:   PetscStrlen(name,&len);

2085:   PetscStrcasecmp(name,"PETSC_DEFAULT",&match);
2086:   if (!match) PetscStrcasecmp(name,"DEFAULT",&match);
2087:   if (match) {*a = PETSC_DEFAULT; return 0;}

2089:   PetscStrcasecmp(name,"PETSC_DECIDE",&match);
2090:   if (!match) PetscStrcasecmp(name,"DECIDE",&match);
2091:   if (match) {*a = PETSC_DECIDE; return 0;}

2093:   PetscStrtod(name,a,&endptr);
2095:   return 0;
2096: }

2098: PetscErrorCode PetscOptionsStringToScalar(const char name[],PetscScalar *a)
2099: {
2100:   PetscBool    imag1;
2101:   size_t       len;
2102:   PetscScalar  val = 0.;
2103:   char        *ptr = NULL;

2105:   PetscStrlen(name,&len);
2107:   PetscStrtoz(name,&val,&ptr,&imag1);
2108: #if defined(PETSC_USE_COMPLEX)
2109:   if ((size_t) (ptr - name) < len) {
2110:     PetscBool   imag2;
2111:     PetscScalar val2;

2113:     PetscStrtoz(ptr,&val2,&ptr,&imag2);
2115:     val = PetscCMPLX(PetscRealPart(val),PetscImaginaryPart(val2));
2116:   }
2117: #endif
2119:   *a = val;
2120:   return 0;
2121: }

2123: /*@C
2124:    PetscOptionsGetBool - Gets the Logical (true or false) value for a particular
2125:             option in the database.

2127:    Not Collective

2129:    Input Parameters:
2130: +  options - options database, use NULL for default global database
2131: .  pre - the string to prepend to the name or NULL
2132: -  name - the option one is seeking

2134:    Output Parameters:
2135: +  ivalue - the logical value to return
2136: -  set - PETSC_TRUE  if found, else PETSC_FALSE

2138:    Level: beginner

2140:    Notes:
2141:        TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE
2142:        FALSE, false, NO, no, and 0 all translate to PETSC_FALSE

2144:       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
2145:      is equivalent to -requested_bool true

2147:        If the user does not supply the option at all ivalue is NOT changed. Thus
2148:      you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.

2150: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
2151:           PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsGetInt(), PetscOptionsBool(),
2152:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2153:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2154:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2155:           PetscOptionsFList(), PetscOptionsEList()
2156: @*/
2157: PetscErrorCode PetscOptionsGetBool(PetscOptions options,const char pre[],const char name[],PetscBool *ivalue,PetscBool *set)
2158: {
2159:   const char     *value;
2160:   PetscBool      flag;

2164:   PetscOptionsFindPair(options,pre,name,&value,&flag);
2165:   if (flag) {
2166:     if (set) *set = PETSC_TRUE;
2167:     PetscOptionsStringToBool(value, &flag);
2168:     if (ivalue) *ivalue = flag;
2169:   } else {
2170:     if (set) *set = PETSC_FALSE;
2171:   }
2172:   return 0;
2173: }

2175: /*@C
2176:    PetscOptionsGetEList - Puts a list of option values that a single one may be selected from

2178:    Not Collective

2180:    Input Parameters:
2181: +  options - options database, use NULL for default global database
2182: .  pre - the string to prepend to the name or NULL
2183: .  opt - option name
2184: .  list - the possible choices (one of these must be selected, anything else is invalid)
2185: -  ntext - number of choices

2187:    Output Parameters:
2188: +  value - the index of the value to return (defaults to zero if the option name is given but no choice is listed)
2189: -  set - PETSC_TRUE if found, else PETSC_FALSE

2191:    Level: intermediate

2193:    Notes:
2194:     If the user does not supply the option value is NOT changed. Thus
2195:      you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.

2197:    See PetscOptionsFList() for when the choices are given in a PetscFunctionList()

2199: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
2200:           PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2201:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2202:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2203:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2204:           PetscOptionsFList(), PetscOptionsEList()
2205: @*/
2206: PetscErrorCode PetscOptionsGetEList(PetscOptions options,const char pre[],const char opt[],const char * const *list,PetscInt ntext,PetscInt *value,PetscBool *set)
2207: {
2208:   size_t         alen,len = 0, tlen = 0;
2209:   char           *svalue;
2210:   PetscBool      aset,flg = PETSC_FALSE;
2211:   PetscInt       i;

2214:   for (i=0; i<ntext; i++) {
2215:     PetscStrlen(list[i],&alen);
2216:     if (alen > len) len = alen;
2217:     tlen += len + 1;
2218:   }
2219:   len += 5; /* a little extra space for user mistypes */
2220:   PetscMalloc1(len,&svalue);
2221:   PetscOptionsGetString(options,pre,opt,svalue,len,&aset);
2222:   if (aset) {
2223:     PetscEListFind(ntext,list,svalue,value,&flg);
2224:     if (!flg) {
2225:       char *avail,*pavl;

2227:       PetscMalloc1(tlen,&avail);
2228:       pavl = avail;
2229:       for (i=0; i<ntext; i++) {
2230:         PetscStrlen(list[i],&alen);
2231:         PetscStrcpy(pavl,list[i]);
2232:         pavl += alen;
2233:         PetscStrcpy(pavl," ");
2234:         pavl += 1;
2235:       }
2236:       PetscStrtolower(avail);
2237:       SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Unknown option %s for -%s%s. Available options: %s",svalue,pre ? pre : "",opt+1,avail);
2238:     }
2239:     if (set) *set = PETSC_TRUE;
2240:   } else if (set) *set = PETSC_FALSE;
2241:   PetscFree(svalue);
2242:   return 0;
2243: }

2245: /*@C
2246:    PetscOptionsGetEnum - Gets the enum value for a particular option in the database.

2248:    Not Collective

2250:    Input Parameters:
2251: +  options - options database, use NULL for default global database
2252: .  pre - option prefix or NULL
2253: .  opt - option name
2254: -  list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null

2256:    Output Parameters:
2257: +  value - the  value to return
2258: -  set - PETSC_TRUE if found, else PETSC_FALSE

2260:    Level: beginner

2262:    Notes:
2263:     If the user does not supply the option value is NOT changed. Thus
2264:      you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.

2266:           List is usually something like PCASMTypes or some other predefined list of enum names

2268: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), PetscOptionsGetInt(),
2269:           PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
2270:           PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
2271:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2272:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2273:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2274:           PetscOptionsFList(), PetscOptionsEList(), PetscOptionsGetEList(), PetscOptionsEnum()
2275: @*/
2276: PetscErrorCode PetscOptionsGetEnum(PetscOptions options,const char pre[],const char opt[],const char * const *list,PetscEnum *value,PetscBool *set)
2277: {
2278:   PetscInt       ntext = 0,tval;
2279:   PetscBool      fset;

2282:   while (list[ntext++]) {
2284:   }
2286:   ntext -= 3;
2287:   PetscOptionsGetEList(options,pre,opt,list,ntext,&tval,&fset);
2288:   /* with PETSC_USE_64BIT_INDICES sizeof(PetscInt) != sizeof(PetscEnum) */
2289:   if (fset) *value = (PetscEnum)tval;
2290:   if (set) *set = fset;
2291:   return 0;
2292: }

2294: /*@C
2295:    PetscOptionsGetInt - Gets the integer value for a particular option in the database.

2297:    Not Collective

2299:    Input Parameters:
2300: +  options - options database, use NULL for default global database
2301: .  pre - the string to prepend to the name or NULL
2302: -  name - the option one is seeking

2304:    Output Parameters:
2305: +  ivalue - the integer value to return
2306: -  set - PETSC_TRUE if found, else PETSC_FALSE

2308:    Level: beginner

2310:    Notes:
2311:    If the user does not supply the option ivalue is NOT changed. Thus
2312:    you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.

2314: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
2315:           PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
2316:           PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
2317:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2318:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2319:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2320:           PetscOptionsFList(), PetscOptionsEList()
2321: @*/
2322: PetscErrorCode PetscOptionsGetInt(PetscOptions options,const char pre[],const char name[],PetscInt *ivalue,PetscBool *set)
2323: {
2324:   const char     *value;
2325:   PetscBool      flag;

2329:   PetscOptionsFindPair(options,pre,name,&value,&flag);
2330:   if (flag) {
2331:     if (!value) {
2332:       if (set) *set = PETSC_FALSE;
2333:     } else {
2334:       if (set) *set = PETSC_TRUE;
2335:       PetscOptionsStringToInt(value,ivalue);
2336:     }
2337:   } else {
2338:     if (set) *set = PETSC_FALSE;
2339:   }
2340:   return 0;
2341: }

2343: /*@C
2344:    PetscOptionsGetReal - Gets the double precision value for a particular
2345:    option in the database.

2347:    Not Collective

2349:    Input Parameters:
2350: +  options - options database, use NULL for default global database
2351: .  pre - string to prepend to each name or NULL
2352: -  name - the option one is seeking

2354:    Output Parameters:
2355: +  dvalue - the double value to return
2356: -  set - PETSC_TRUE if found, PETSC_FALSE if not found

2358:    Notes:
2359:     If the user does not supply the option dvalue is NOT changed. Thus
2360:      you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.

2362:    Level: beginner

2364: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2365:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(),PetscOptionsBool(),
2366:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2367:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2368:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2369:           PetscOptionsFList(), PetscOptionsEList()
2370: @*/
2371: PetscErrorCode PetscOptionsGetReal(PetscOptions options,const char pre[],const char name[],PetscReal *dvalue,PetscBool *set)
2372: {
2373:   const char     *value;
2374:   PetscBool      flag;

2378:   PetscOptionsFindPair(options,pre,name,&value,&flag);
2379:   if (flag) {
2380:     if (!value) {
2381:       if (set) *set = PETSC_FALSE;
2382:     } else {
2383:       if (set) *set = PETSC_TRUE;
2384:       PetscOptionsStringToReal(value,dvalue);
2385:     }
2386:   } else {
2387:     if (set) *set = PETSC_FALSE;
2388:   }
2389:   return 0;
2390: }

2392: /*@C
2393:    PetscOptionsGetScalar - Gets the scalar value for a particular
2394:    option in the database.

2396:    Not Collective

2398:    Input Parameters:
2399: +  options - options database, use NULL for default global database
2400: .  pre - string to prepend to each name or NULL
2401: -  name - the option one is seeking

2403:    Output Parameters:
2404: +  dvalue - the double value to return
2405: -  set - PETSC_TRUE if found, else PETSC_FALSE

2407:    Level: beginner

2409:    Usage:
2410:    A complex number 2+3i must be specified with NO spaces

2412:    Notes:
2413:     If the user does not supply the option dvalue is NOT changed. Thus
2414:      you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true.

2416: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2417:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2418:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2419:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2420:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2421:           PetscOptionsFList(), PetscOptionsEList()
2422: @*/
2423: PetscErrorCode PetscOptionsGetScalar(PetscOptions options,const char pre[],const char name[],PetscScalar *dvalue,PetscBool *set)
2424: {
2425:   const char     *value;
2426:   PetscBool      flag;

2430:   PetscOptionsFindPair(options,pre,name,&value,&flag);
2431:   if (flag) {
2432:     if (!value) {
2433:       if (set) *set = PETSC_FALSE;
2434:     } else {
2435: #if !defined(PETSC_USE_COMPLEX)
2436:       PetscOptionsStringToReal(value,dvalue);
2437: #else
2438:       PetscOptionsStringToScalar(value,dvalue);
2439: #endif
2440:       if (set) *set = PETSC_TRUE;
2441:     }
2442:   } else { /* flag */
2443:     if (set) *set = PETSC_FALSE;
2444:   }
2445:   return 0;
2446: }

2448: /*@C
2449:    PetscOptionsGetString - Gets the string value for a particular option in
2450:    the database.

2452:    Not Collective

2454:    Input Parameters:
2455: +  options - options database, use NULL for default global database
2456: .  pre - string to prepend to name or NULL
2457: .  name - the option one is seeking
2458: -  len - maximum length of the string including null termination

2460:    Output Parameters:
2461: +  string - location to copy string
2462: -  set - PETSC_TRUE if found, else PETSC_FALSE

2464:    Level: beginner

2466:    Fortran Note:
2467:    The Fortran interface is slightly different from the C/C++
2468:    interface (len is not used).  Sample usage in Fortran follows
2469: .vb
2470:       character *20    string
2471:       PetscErrorCode   ierr
2472:       PetscBool        set
2473:       call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s',string,set,ierr)
2474: .ve

2476:    Notes:
2477:     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

2479:            If the user does not use the option then the string is not changed. Thus
2480:            you should ALWAYS initialize the string if you access it without first checking if the set flag is true.

2482:     Note:
2483:       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).

2485: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
2486:           PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2487:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2488:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2489:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2490:           PetscOptionsFList(), PetscOptionsEList()
2491: @*/
2492: PetscErrorCode PetscOptionsGetString(PetscOptions options,const char pre[],const char name[],char string[],size_t len,PetscBool *set)
2493: {
2494:   const char     *value;
2495:   PetscBool      flag;

2499:   PetscOptionsFindPair(options,pre,name,&value,&flag);
2500:   if (!flag) {
2501:     if (set) *set = PETSC_FALSE;
2502:   } else {
2503:     if (set) *set = PETSC_TRUE;
2504:     if (value) PetscStrncpy(string,value,len);
2505:     else PetscArrayzero(string,len);
2506:   }
2507:   return 0;
2508: }

2510: char *PetscOptionsGetStringMatlab(PetscOptions options,const char pre[],const char name[])
2511: {
2512:   const char *value;
2513:   PetscBool   flag;

2515:   if (PetscOptionsFindPair(options,pre,name,&value,&flag)) return NULL;
2516:   if (flag) PetscFunctionReturn((char*)value);
2517:   return NULL;
2518: }

2520: /*@C
2521:   PetscOptionsGetBoolArray - Gets an array of Logical (true or false) values for a particular
2522:   option in the database.  The values must be separated with commas with no intervening spaces.

2524:   Not Collective

2526:   Input Parameters:
2527: + options - options database, use NULL for default global database
2528: . pre - string to prepend to each name or NULL
2529: - name - the option one is seeking

2531:   Output Parameters:
2532: + dvalue - the integer values to return
2533: . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved
2534: - set - PETSC_TRUE if found, else PETSC_FALSE

2536:   Level: beginner

2538:   Notes:
2539:   TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE. FALSE, false, NO, no, and 0 all translate to PETSC_FALSE

2541: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2542:           PetscOptionsGetString(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2543:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2544:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2545:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2546:           PetscOptionsFList(), PetscOptionsEList()
2547: @*/
2548: PetscErrorCode PetscOptionsGetBoolArray(PetscOptions options,const char pre[],const char name[],PetscBool dvalue[],PetscInt *nmax,PetscBool *set)
2549: {
2550:   const char     *svalue;
2551:   char           *value;
2552:   PetscInt       n = 0;
2553:   PetscBool      flag;
2554:   PetscToken     token;


2560:   PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2561:   if (!flag || !svalue)  { if (set) *set = PETSC_FALSE; *nmax = 0; return 0;}
2562:   if (set) *set = PETSC_TRUE;
2563:   PetscTokenCreate(svalue,',',&token);
2564:   PetscTokenFind(token,&value);
2565:   while (value && n < *nmax) {
2566:     PetscOptionsStringToBool(value,dvalue);
2567:     PetscTokenFind(token,&value);
2568:     dvalue++;
2569:     n++;
2570:   }
2571:   PetscTokenDestroy(&token);
2572:   *nmax = n;
2573:   return 0;
2574: }

2576: /*@C
2577:   PetscOptionsGetEnumArray - Gets an array of enum values for a particular option in the database.

2579:   Not Collective

2581:   Input Parameters:
2582: + options - options database, use NULL for default global database
2583: . pre - option prefix or NULL
2584: . name - option name
2585: - list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null

2587:   Output Parameters:
2588: + ivalue - the  enum values to return
2589: . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved
2590: - set - PETSC_TRUE if found, else PETSC_FALSE

2592:   Level: beginner

2594:   Notes:
2595:   The array must be passed as a comma separated list.

2597:   There must be no intervening spaces between the values.

2599:   list is usually something like PCASMTypes or some other predefined list of enum names.

2601: .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), PetscOptionsGetInt(),
2602:           PetscOptionsGetEnum(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
2603:           PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(), PetscOptionsName(),
2604:           PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), PetscOptionsStringArray(),PetscOptionsRealArray(),
2605:           PetscOptionsScalar(), PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2606:           PetscOptionsFList(), PetscOptionsEList(), PetscOptionsGetEList(), PetscOptionsEnum()
2607: @*/
2608: PetscErrorCode PetscOptionsGetEnumArray(PetscOptions options,const char pre[],const char name[],const char *const *list,PetscEnum ivalue[],PetscInt *nmax,PetscBool *set)
2609: {
2610:   const char     *svalue;
2611:   char           *value;
2612:   PetscInt       n = 0;
2613:   PetscEnum      evalue;
2614:   PetscBool      flag;
2615:   PetscToken     token;


2622:   PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2623:   if (!flag || !svalue)  { if (set) *set = PETSC_FALSE; *nmax = 0; return 0;}
2624:   if (set) *set = PETSC_TRUE;
2625:   PetscTokenCreate(svalue,',',&token);
2626:   PetscTokenFind(token,&value);
2627:   while (value && n < *nmax) {
2628:     PetscEnumFind(list,value,&evalue,&flag);
2630:     ivalue[n++] = evalue;
2631:     PetscTokenFind(token,&value);
2632:   }
2633:   PetscTokenDestroy(&token);
2634:   *nmax = n;
2635:   return 0;
2636: }

2638: /*@C
2639:   PetscOptionsGetIntArray - Gets an array of integer values for a particular option in the database.

2641:   Not Collective

2643:   Input Parameters:
2644: + options - options database, use NULL for default global database
2645: . pre - string to prepend to each name or NULL
2646: - name - the option one is seeking

2648:   Output Parameters:
2649: + ivalue - the integer values to return
2650: . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved
2651: - set - PETSC_TRUE if found, else PETSC_FALSE

2653:   Level: beginner

2655:   Notes:
2656:   The array can be passed as
2657: .vb
2658:   a comma separated list:                                 0,1,2,3,4,5,6,7
2659:   a range (start-end+1):                                  0-8
2660:   a range with given increment (start-end+1:inc):         0-7:2
2661:   a combination of values and ranges separated by commas: 0,1-8,8-15:2
2662: .ve

2664:   There must be no intervening spaces between the values.

2666: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2667:           PetscOptionsGetString(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2668:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2669:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2670:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2671:           PetscOptionsFList(), PetscOptionsEList()
2672: @*/
2673: PetscErrorCode PetscOptionsGetIntArray(PetscOptions options,const char pre[],const char name[],PetscInt ivalue[],PetscInt *nmax,PetscBool *set)
2674: {
2675:   const char     *svalue;
2676:   char           *value;
2677:   PetscInt       n = 0,i,j,start,end,inc,nvalues;
2678:   size_t         len;
2679:   PetscBool      flag,foundrange;
2680:   PetscToken     token;


2686:   PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2687:   if (!flag || !svalue)  { if (set) *set = PETSC_FALSE; *nmax = 0; return 0;}
2688:   if (set) *set = PETSC_TRUE;
2689:   PetscTokenCreate(svalue,',',&token);
2690:   PetscTokenFind(token,&value);
2691:   while (value && n < *nmax) {
2692:     /* look for form  d-D where d and D are integers */
2693:     foundrange = PETSC_FALSE;
2694:     PetscStrlen(value,&len);
2695:     if (value[0] == '-') i=2;
2696:     else i=1;
2697:     for (;i<(int)len; i++) {
2698:       if (value[i] == '-') {
2700:         value[i] = 0;

2702:         PetscOptionsStringToInt(value,&start);
2703:         inc  = 1;
2704:         j    = i+1;
2705:         for (;j<(int)len; j++) {
2706:           if (value[j] == ':') {
2707:             value[j] = 0;

2709:             PetscOptionsStringToInt(value+j+1,&inc);
2711:             break;
2712:           }
2713:         }
2714:         PetscOptionsStringToInt(value+i+1,&end);
2716:         nvalues = (end-start)/inc + (end-start)%inc;
2718:         for (;start<end; start+=inc) {
2719:           *ivalue = start; ivalue++;n++;
2720:         }
2721:         foundrange = PETSC_TRUE;
2722:         break;
2723:       }
2724:     }
2725:     if (!foundrange) {
2726:       PetscOptionsStringToInt(value,ivalue);
2727:       ivalue++;
2728:       n++;
2729:     }
2730:     PetscTokenFind(token,&value);
2731:   }
2732:   PetscTokenDestroy(&token);
2733:   *nmax = n;
2734:   return 0;
2735: }

2737: /*@C
2738:   PetscOptionsGetRealArray - Gets an array of double precision values for a
2739:   particular option in the database.  The values must be separated with commas with no intervening spaces.

2741:   Not Collective

2743:   Input Parameters:
2744: + options - options database, use NULL for default global database
2745: . pre - string to prepend to each name or NULL
2746: - name - the option one is seeking

2748:   Output Parameters:
2749: + dvalue - the double values to return
2750: . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved
2751: - set - PETSC_TRUE if found, else PETSC_FALSE

2753:   Level: beginner

2755: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2756:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
2757:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2758:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2759:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2760:           PetscOptionsFList(), PetscOptionsEList()
2761: @*/
2762: PetscErrorCode PetscOptionsGetRealArray(PetscOptions options,const char pre[],const char name[],PetscReal dvalue[],PetscInt *nmax,PetscBool *set)
2763: {
2764:   const char     *svalue;
2765:   char           *value;
2766:   PetscInt       n = 0;
2767:   PetscBool      flag;
2768:   PetscToken     token;


2774:   PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2775:   if (!flag || !svalue)  { if (set) *set = PETSC_FALSE; *nmax = 0; return 0;}
2776:   if (set) *set = PETSC_TRUE;
2777:   PetscTokenCreate(svalue,',',&token);
2778:   PetscTokenFind(token,&value);
2779:   while (value && n < *nmax) {
2780:     PetscOptionsStringToReal(value,dvalue++);
2781:     PetscTokenFind(token,&value);
2782:     n++;
2783:   }
2784:   PetscTokenDestroy(&token);
2785:   *nmax = n;
2786:   return 0;
2787: }

2789: /*@C
2790:   PetscOptionsGetScalarArray - Gets an array of scalars for a
2791:   particular option in the database.  The values must be separated with commas with no intervening spaces.

2793:   Not Collective

2795:   Input Parameters:
2796: + options - options database, use NULL for default global database
2797: . pre - string to prepend to each name or NULL
2798: - name - the option one is seeking

2800:   Output Parameters:
2801: + dvalue - the scalar values to return
2802: . nmax - On input maximum number of values to retrieve, on output the actual number of values retrieved
2803: - set - PETSC_TRUE if found, else PETSC_FALSE

2805:   Level: beginner

2807: .seealso: PetscOptionsGetInt(), PetscOptionsHasName(),
2808:           PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(),
2809:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2810:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2811:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2812:           PetscOptionsFList(), PetscOptionsEList()
2813: @*/
2814: PetscErrorCode PetscOptionsGetScalarArray(PetscOptions options,const char pre[],const char name[],PetscScalar dvalue[],PetscInt *nmax,PetscBool *set)
2815: {
2816:   const char     *svalue;
2817:   char           *value;
2818:   PetscInt       n = 0;
2819:   PetscBool      flag;
2820:   PetscToken     token;


2826:   PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2827:   if (!flag || !svalue)  { if (set) *set = PETSC_FALSE; *nmax = 0; return 0;}
2828:   if (set) *set = PETSC_TRUE;
2829:   PetscTokenCreate(svalue,',',&token);
2830:   PetscTokenFind(token,&value);
2831:   while (value && n < *nmax) {
2832:     PetscOptionsStringToScalar(value,dvalue++);
2833:     PetscTokenFind(token,&value);
2834:     n++;
2835:   }
2836:   PetscTokenDestroy(&token);
2837:   *nmax = n;
2838:   return 0;
2839: }

2841: /*@C
2842:   PetscOptionsGetStringArray - Gets an array of string values for a particular
2843:   option in the database. The values must be separated with commas with no intervening spaces.

2845:   Not Collective

2847:   Input Parameters:
2848: + options - options database, use NULL for default global database
2849: . pre - string to prepend to name or NULL
2850: - name - the option one is seeking

2852:   Output Parameters:
2853: + strings - location to copy strings
2854: . nmax - On input maximum number of strings, on output the actual number of strings found
2855: - set - PETSC_TRUE if found, else PETSC_FALSE

2857:   Level: beginner

2859:   Notes:
2860:   The nmax parameter is used for both input and output.

2862:   The user should pass in an array of pointers to char, to hold all the
2863:   strings returned by this function.

2865:   The user is responsible for deallocating the strings that are
2866:   returned. The Fortran interface for this routine is not supported.

2868: .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),
2869:           PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(),
2870:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
2871:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
2872:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
2873:           PetscOptionsFList(), PetscOptionsEList()
2874: @*/
2875: PetscErrorCode PetscOptionsGetStringArray(PetscOptions options,const char pre[],const char name[],char *strings[],PetscInt *nmax,PetscBool *set)
2876: {
2877:   const char     *svalue;
2878:   char           *value;
2879:   PetscInt       n = 0;
2880:   PetscBool      flag;
2881:   PetscToken     token;


2887:   PetscOptionsFindPair(options,pre,name,&svalue,&flag);
2888:   if (!flag || !svalue)  { if (set) *set = PETSC_FALSE; *nmax = 0; return 0;}
2889:   if (set) *set = PETSC_TRUE;
2890:   PetscTokenCreate(svalue,',',&token);
2891:   PetscTokenFind(token,&value);
2892:   while (value && n < *nmax) {
2893:     PetscStrallocpy(value,&strings[n]);
2894:     PetscTokenFind(token,&value);
2895:     n++;
2896:   }
2897:   PetscTokenDestroy(&token);
2898:   *nmax = n;
2899:   return 0;
2900: }

2902: /*@C
2903:    PetscOptionsDeprecated - mark an option as deprecated, optionally replacing it with a new one

2905:    Prints a deprecation warning, unless an option is supplied to suppress.

2907:    Logically Collective

2909:    Input Parameters:
2910: +  pre - string to prepend to name or NULL
2911: .  oldname - the old, deprecated option
2912: .  newname - the new option, or NULL if option is purely removed
2913: .  version - a string describing the version of first deprecation, e.g. "3.9"
2914: -  info - additional information string, or NULL.

2916:    Options Database Keys:
2917: . -options_suppress_deprecated_warnings - do not print deprecation warnings

2919:    Notes:
2920:    Must be called between PetscOptionsBegin() (or PetscObjectOptionsBegin()) and PetscOptionsEnd().
2921:    Only the proces of rank zero that owns the PetscOptionsItems are argument (managed by PetscOptionsBegin() or
2922:    PetscObjectOptionsBegin() prints the information
2923:    If newname is provided, the old option is replaced. Otherwise, it remains
2924:    in the options database.
2925:    If an option is not replaced, the info argument should be used to advise the user
2926:    on how to proceed.
2927:    There is a limit on the length of the warning printed, so very long strings
2928:    provided as info may be truncated.

2930:    Level: developer

2932: .seealso: PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsScalar(), PetscOptionsBool(), PetscOptionsString(), PetscOptionsSetValue()

2934: @*/
2935: PetscErrorCode PetscOptionsDeprecated_Private(PetscOptionItems *PetscOptionsObject,const char oldname[],const char newname[],const char version[],const char info[])
2936: {
2937:   PetscBool          found,quiet;
2938:   const char         *value;
2939:   const char * const quietopt="-options_suppress_deprecated_warnings";
2940:   char               msg[4096];
2941:   char               *prefix = NULL;
2942:   PetscOptions       options = NULL;
2943:   MPI_Comm           comm = PETSC_COMM_SELF;

2947:   if (PetscOptionsObject) {
2948:     prefix  = PetscOptionsObject->prefix;
2949:     options = PetscOptionsObject->options;
2950:     comm    = PetscOptionsObject->comm;
2951:   }
2952:   PetscOptionsFindPair(options,prefix,oldname,&value,&found);
2953:   if (found) {
2954:     if (newname) {
2955:       if (prefix) {
2956:         PetscOptionsPrefixPush(options,prefix);
2957:       }
2958:       PetscOptionsSetValue(options,newname,value);
2959:       if (prefix) {
2960:         PetscOptionsPrefixPop(options);
2961:       }
2962:       PetscOptionsClearValue(options,oldname);
2963:     }
2964:     quiet = PETSC_FALSE;
2965:     PetscOptionsGetBool(options,NULL,quietopt,&quiet,NULL);
2966:     if (!quiet) {
2967:       PetscStrcpy(msg,"** PETSc DEPRECATION WARNING ** : the option ");
2968:       PetscStrcat(msg,oldname);
2969:       PetscStrcat(msg," is deprecated as of version ");
2970:       PetscStrcat(msg,version);
2971:       PetscStrcat(msg," and will be removed in a future release.");
2972:       if (newname) {
2973:         PetscStrcat(msg," Please use the option ");
2974:         PetscStrcat(msg,newname);
2975:         PetscStrcat(msg," instead.");
2976:       }
2977:       if (info) {
2978:         PetscStrcat(msg," ");
2979:         PetscStrcat(msg,info);
2980:       }
2981:       PetscStrcat(msg," (Silence this warning with ");
2982:       PetscStrcat(msg,quietopt);
2983:       PetscStrcat(msg,")\n");
2984:       PetscPrintf(comm,"%s",msg);
2985:     }
2986:   }
2987:   return 0;
2988: }