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 lenght
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:   PetscBool         *set;

701:   PetscCalloc2(n,&val,n,&set);

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

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

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

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

734:   /* Store precedent options in database and mark them as used */
735:   for (o=0; o<n; o++) {
736:     if (set[o]) {
737:       PetscOptionsSetValue_Private(options,opt[o],val[o],&a);
738:       options->used[a] = PETSC_TRUE;
739:     }
740:   }

742:   PetscFree2(val,set);
743:   options->precedentProcessed = PETSC_TRUE;
744:   return 0;
745: }

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

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

767:    Collective on PETSC_COMM_WORLD

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

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

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

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

789:    Level: advanced

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

802:   MPI_Comm_rank(comm,&rank);

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

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

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

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

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

872:    Logically Collective on PetscViewer

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

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

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

885:    Level: advanced

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

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

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

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

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

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

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

943:    Logically Collective

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

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

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

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

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

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

965: Level: advanced

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

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

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

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

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

1000:    Level: advanced

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

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

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

1019:     Logically Collective

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

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

1029:    Level: developer

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

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

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

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

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

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

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

1065:    Logically Collective

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

1072:    Level: advanced

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

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

1094:   options = options ? options : defaultoptions;
1095:   PetscOptionsValidKey(newname,&valid);
1097:   PetscOptionsValidKey(oldname,&valid);

1100:   n = options->Naliases;

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

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

1118:    Logically Collective

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

1125:    Level: intermediate

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

1232:   PetscOptionsMonitor(options,name,value);
1233:   if (pos) *pos = n;
1234:   return 0;
1235: }

1237: /*@C
1238:    PetscOptionsClearValue - Clears an option name-value pair in the options
1239:    database, overriding whatever is already present.

1241:    Logically Collective

1243:    Input Parameters:
1244: +  options - options database, use NULL for the default global database
1245: -  name - name of option, this SHOULD have the - prepended

1247:    Level: intermediate

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

1254: .seealso: PetscOptionsInsert()
1255: @*/
1256: PetscErrorCode PetscOptionsClearValue(PetscOptions options,const char name[])
1257: {
1258:   int            N,n,i;
1259:   char           **names;

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

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

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

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

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

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

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

1302:    Not Collective

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

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

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

1316:    Level: developer

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

1326:   options = options ? options : defaultoptions;

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

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

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

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

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

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

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

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

1434:   options = options ? options : defaultoptions;

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

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

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

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

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

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

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

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

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

1510:    Not Collective

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

1518:    Level: advanced

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

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

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

1542:    Not Collective

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

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

1550:    Level: advanced

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

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

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

1574:    Not Collective

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

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

1584:    Level: beginner

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

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

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

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

1609:    Not Collective

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

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

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

1621:    Level: advanced

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

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

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

1660:    Not Collective

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

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

1669:    Level: advanced

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

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

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

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

1699:    Not Collective

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

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

1707:    Level: advanced

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

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

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

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

1731:    Not Collective

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

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

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

1745:    Level: advanced

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

1755:   toptions = options ? options : defaultoptions;
1756:   for (i=0; i<toptions->N; i++) {
1757:     if (!toptions->used[i]) {
1758:       if (toptions->values[i]) {
1759:         PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s value: %s\n",toptions->names[i],toptions->values[i]);
1760:       } else {
1761:         PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s (no value)\n",toptions->names[i]);
1762:       }
1763:     }
1764:   }
1765:   if (!options) {
1766:     toptions = defaultoptions;
1767:     while (toptions->previous) {
1768:       cnt++;
1769:       toptions = toptions->previous;
1770:     }
1771:     if (cnt) {
1772:       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);
1773:     }
1774:   }
1775:   return 0;
1776: }

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

1781:    Not Collective

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

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

1791:    Level: advanced

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

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

1807:   options = options ? options : defaultoptions;

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

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

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

1834:    Not Collective

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

1841:    Level: advanced

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

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

1859:    Logically Collective on ctx

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

1866:    Level: intermediate

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

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

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

1903:    Not Collective

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

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

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

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

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

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

1931:    Level: intermediate

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

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

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

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

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

1986:   PetscStrlen(name,&len);

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

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

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

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

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

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

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

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

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

2081:   PetscStrlen(name,&len);

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

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

2092:   PetscStrtod(name,a,&endptr);
2094:   return 0;
2095: }

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

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

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

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

2126:    Not Collective

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

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

2137:    Level: beginner

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

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

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

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

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

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

2177:    Not Collective

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

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

2190:    Level: intermediate

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

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

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

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

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

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

2247:    Not Collective

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

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

2259:    Level: beginner

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

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

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

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

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

2296:    Not Collective

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

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

2307:    Level: beginner

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

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

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

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

2346:    Not Collective

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

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

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

2361:    Level: beginner

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

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

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

2395:    Not Collective

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

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

2406:    Level: beginner

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

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

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

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

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

2451:    Not Collective

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

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

2463:    Level: beginner

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

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

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

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

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

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

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

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

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

2523:   Not Collective

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

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

2535:   Level: beginner

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

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


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

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

2578:   Not Collective

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

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

2591:   Level: beginner

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

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

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

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


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

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

2640:   Not Collective

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

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

2652:   Level: beginner

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

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

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


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

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

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

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

2740:   Not Collective

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

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

2752:   Level: beginner

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


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

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

2792:   Not Collective

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

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

2804:   Level: beginner

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


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

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

2844:   Not Collective

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

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

2856:   Level: beginner

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

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

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

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


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

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

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

2906:    Logically Collective

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

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

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

2929:    Level: developer

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

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

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