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() */ 3 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 7 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 */ 11 12 #include <petsc/private/petscimpl.h> /*I "petscsys.h" I*/ 13 #include <petscviewer.h> 14 #include <ctype.h> 15 #if defined(PETSC_HAVE_MALLOC_H) 16 #include <malloc.h> 17 #endif 18 #if defined(PETSC_HAVE_STRINGS_H) 19 # include <strings.h> /* strcasecmp */ 20 #endif 21 #if defined(PETSC_HAVE_YAML) 22 #include <yaml.h> 23 #endif 24 25 #if defined(PETSC_HAVE_STRCASECMP) 26 #define PetscOptNameCmp(a,b) strcasecmp(a,b) 27 #elif defined(PETSC_HAVE_STRICMP) 28 #define PetscOptNameCmp(a,b) stricmp(a,b) 29 #else 30 #define PetscOptNameCmp(a,b) Error_strcasecmp_not_found 31 #endif 32 33 #include <petsc/private/hashtable.h> 34 35 /* This assumes ASCII encoding and ignores locale settings */ 36 /* Using tolower() is about 2X slower in microbenchmarks */ 37 PETSC_STATIC_INLINE int PetscToLower(int c) 38 { 39 return ((c >= 'A') & (c <= 'Z')) ? c + 'a' - 'A' : c; 40 } 41 42 /* Bob Jenkins's one at a time hash function (case-insensitive) */ 43 PETSC_STATIC_INLINE unsigned int PetscOptHash(const char key[]) 44 { 45 unsigned int hash = 0; 46 while (*key) { 47 hash += PetscToLower(*key++); 48 hash += hash << 10; 49 hash ^= hash >> 6; 50 } 51 hash += hash << 3; 52 hash ^= hash >> 11; 53 hash += hash << 15; 54 return hash; 55 } 56 57 PETSC_STATIC_INLINE int PetscOptEqual(const char a[],const char b[]) 58 { 59 return !PetscOptNameCmp(a,b); 60 } 61 62 KHASH_INIT(HO, kh_cstr_t, int, 1, PetscOptHash, PetscOptEqual) 63 64 /* 65 This table holds all the options set by the user. For simplicity, we use a static size database 66 */ 67 #define MAXOPTNAME 512 68 #define MAXOPTIONS 512 69 #define MAXALIASES 25 70 #define MAXPREFIXES 25 71 #define MAXOPTIONSMONITORS 5 72 73 struct _n_PetscOptions { 74 PetscOptions previous; 75 int N; /* number of options */ 76 char *names[MAXOPTIONS]; /* option names */ 77 char *values[MAXOPTIONS]; /* option values */ 78 PetscBool used[MAXOPTIONS]; /* flag option use */ 79 80 /* Hash table */ 81 khash_t(HO) *ht; 82 83 /* Prefixes */ 84 int prefixind; 85 int prefixstack[MAXPREFIXES]; 86 char prefix[MAXOPTNAME]; 87 88 /* Aliases */ 89 int Naliases; /* number or aliases */ 90 char *aliases1[MAXALIASES]; /* aliased */ 91 char *aliases2[MAXALIASES]; /* aliasee */ 92 93 /* Help */ 94 PetscBool help; /* flag whether "-help" is in the database */ 95 96 /* Monitors */ 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 }; 102 103 static PetscOptions defaultoptions = NULL; /* the options database routines query this object for options */ 104 105 /* 106 Options events monitor 107 */ 108 static PetscErrorCode PetscOptionsMonitor(PetscOptions options,const char name[],const char value[]) 109 { 110 PetscInt i; 111 PetscErrorCode ierr; 112 113 if (!PetscInitializeCalled) return 0; 114 PetscFunctionBegin; 115 for (i=0; i<options->numbermonitors; i++) { 116 ierr = (*options->monitor[i])(name,value,options->monitorcontext[i]);CHKERRQ(ierr); 117 } 118 PetscFunctionReturn(0); 119 } 120 121 /*@ 122 PetscOptionsCreate - Creates an empty options database. 123 124 Logically collective 125 126 Output Parameter: 127 . options - Options database object 128 129 Level: advanced 130 131 Developer Note: We may want eventually to pass a MPI_Comm to determine the ownership of the object 132 133 .seealso: PetscOptionsDestroy(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsInsert(), PetscOptionsSetValue() 134 @*/ 135 PetscErrorCode PetscOptionsCreate(PetscOptions *options) 136 { 137 if (!options) return PETSC_ERR_ARG_NULL; 138 *options = (PetscOptions)calloc(1,sizeof(**options)); 139 if (!*options) return PETSC_ERR_MEM; 140 return 0; 141 } 142 143 /*@ 144 PetscOptionsDestroy - Destroys an option database. 145 146 Logically collective on whatever communicator was associated with the call to PetscOptionsCreate() 147 148 Input Parameter: 149 . options - the PetscOptions object 150 151 Level: advanced 152 153 .seealso: PetscOptionsInsert(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsInsert(), PetscOptionsSetValue() 154 @*/ 155 PetscErrorCode PetscOptionsDestroy(PetscOptions *options) 156 { 157 PetscErrorCode ierr; 158 159 if (!*options) return 0; 160 if ((*options)->previous) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"You are destroying an option that has been used with PetscOptionsPush() but does not have a corresponding PetscOptionsPop()"); 161 ierr = PetscOptionsClear(*options);if (ierr) return ierr; 162 /* XXX what about monitors ? */ 163 free(*options); 164 *options = NULL; 165 PetscFunctionReturn(0); 166 } 167 168 /* 169 PetscOptionsCreateDefault - Creates the default global options database 170 */ 171 PetscErrorCode PetscOptionsCreateDefault(void) 172 { 173 PetscErrorCode ierr; 174 175 if (!defaultoptions) { 176 ierr = PetscOptionsCreate(&defaultoptions);if (ierr) return ierr; 177 } 178 return 0; 179 } 180 181 /*@ 182 PetscOptionsPush - Push a new PetscOptions object as the default provider of options 183 Allows using different parts of a code to use different options databases 184 185 Logically Collective 186 187 Input Parameter: 188 . opt - the options obtained with PetscOptionsCreate() 189 190 Notes: 191 Use PetscOptionsPop() to return to the previous default options database 192 193 The collectivity of this routine is complex; only the MPI processes that call this routine will 194 have the affect of these options. If some processes that create objects call this routine and others do 195 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 196 on different ranks. 197 198 Level: advanced 199 200 .seealso: PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsInsert(), PetscOptionsSetValue(), PetscOptionsLeft() 201 202 @*/ 203 PetscErrorCode PetscOptionsPush(PetscOptions opt) 204 { 205 PetscErrorCode ierr; 206 207 PetscFunctionBegin; 208 if (!defaultoptions) { 209 ierr = PetscOptionsCreateDefault();CHKERRQ(ierr); 210 } 211 opt->previous = defaultoptions; 212 defaultoptions = opt; 213 PetscFunctionReturn(0); 214 } 215 216 /*@ 217 PetscOptionsPop - Pop the most recent PetscOptionsPush() to return to the previous default options 218 219 Logically collective on whatever communicator was associated with the call to PetscOptionsCreate() 220 221 Notes: 222 Use PetscOptionsPop() to return to the previous default options database 223 Allows using different parts of a code to use different options databases 224 225 Level: advanced 226 227 .seealso: PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsInsert(), PetscOptionsSetValue(), PetscOptionsLeft() 228 229 @*/ 230 PetscErrorCode PetscOptionsPop(void) 231 { 232 PetscOptions current = defaultoptions; 233 234 PetscFunctionBegin; 235 if (!defaultoptions) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Missing default options"); 236 if (!defaultoptions->previous) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"PetscOptionsPop() called too many times"); 237 defaultoptions = defaultoptions->previous; 238 current->previous = NULL; 239 PetscFunctionReturn(0); 240 } 241 242 /* 243 PetscOptionsDestroyDefault - Destroys the default global options database 244 */ 245 PetscErrorCode PetscOptionsDestroyDefault(void) 246 { 247 PetscErrorCode ierr; 248 PetscOptions tmp; 249 250 /* Destroy any options that the user forgot to pop */ 251 while (defaultoptions->previous) { 252 tmp = defaultoptions; 253 ierr = PetscOptionsPop();CHKERRQ(ierr); 254 ierr = PetscOptionsDestroy(&tmp);CHKERRQ(ierr); 255 } 256 ierr = PetscOptionsDestroy(&defaultoptions);if (ierr) return ierr; 257 return 0; 258 } 259 260 /*@C 261 PetscOptionsValidKey - PETSc Options database keys must begin with one or two dashes (-) followed by a letter. 262 263 Not collective 264 265 Input Parameter: 266 . key - string to check if valid 267 268 Output Parameter: 269 . valid - PETSC_TRUE if a valid key 270 271 Level: intermediate 272 @*/ 273 PetscErrorCode PetscOptionsValidKey(const char key[],PetscBool *valid) 274 { 275 char *ptr; 276 277 PetscFunctionBegin; 278 if (key) PetscValidCharPointer(key,1); 279 PetscValidPointer(valid,2); 280 *valid = PETSC_FALSE; 281 if (!key) PetscFunctionReturn(0); 282 if (key[0] != '-') PetscFunctionReturn(0); 283 if (key[1] == '-') key++; 284 if (!isalpha((int)key[1])) PetscFunctionReturn(0); 285 (void) strtod(key,&ptr); 286 if (ptr != key && !(*ptr == '_' || isalnum((int)*ptr))) PetscFunctionReturn(0); 287 *valid = PETSC_TRUE; 288 PetscFunctionReturn(0); 289 } 290 291 /*@C 292 PetscOptionsInsertString - Inserts options into the database from a string 293 294 Logically Collective 295 296 Input Parameter: 297 . in_str - string that contains options separated by blanks 298 299 Level: intermediate 300 301 The collectivity of this routine is complex; only the MPI processes that call this routine will 302 have the affect of these options. If some processes that create objects call this routine and others do 303 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 304 on different ranks. 305 306 Contributed by Boyana Norris 307 308 .seealso: PetscOptionsSetValue(), PetscOptionsView(), PetscOptionsHasName(), PetscOptionsGetInt(), 309 PetscOptionsGetReal(), PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(), 310 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 311 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 312 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 313 PetscOptionsFList(), PetscOptionsEList(), PetscOptionsInsertFile() 314 @*/ 315 PetscErrorCode PetscOptionsInsertString(PetscOptions options,const char in_str[]) 316 { 317 char *first,*second; 318 PetscErrorCode ierr; 319 PetscToken token; 320 PetscBool key,ispush,ispop,isopts; 321 322 PetscFunctionBegin; 323 ierr = PetscTokenCreate(in_str,' ',&token);CHKERRQ(ierr); 324 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 325 while (first) { 326 ierr = PetscStrcasecmp(first,"-prefix_push",&ispush);CHKERRQ(ierr); 327 ierr = PetscStrcasecmp(first,"-prefix_pop",&ispop);CHKERRQ(ierr); 328 ierr = PetscStrcasecmp(first,"-options_file",&isopts);CHKERRQ(ierr); 329 ierr = PetscOptionsValidKey(first,&key);CHKERRQ(ierr); 330 if (ispush) { 331 ierr = PetscTokenFind(token,&second);CHKERRQ(ierr); 332 ierr = PetscOptionsPrefixPush(options,second);CHKERRQ(ierr); 333 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 334 } else if (ispop) { 335 ierr = PetscOptionsPrefixPop(options);CHKERRQ(ierr); 336 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 337 } else if (isopts) { 338 ierr = PetscTokenFind(token,&second);CHKERRQ(ierr); 339 ierr = PetscOptionsInsertFile(PETSC_COMM_SELF,options,second,PETSC_TRUE);CHKERRQ(ierr); 340 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 341 } else if (key) { 342 ierr = PetscTokenFind(token,&second);CHKERRQ(ierr); 343 ierr = PetscOptionsValidKey(second,&key);CHKERRQ(ierr); 344 if (!key) { 345 ierr = PetscOptionsSetValue(options,first,second);CHKERRQ(ierr); 346 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 347 } else { 348 ierr = PetscOptionsSetValue(options,first,NULL);CHKERRQ(ierr); 349 first = second; 350 } 351 } else { 352 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 353 } 354 } 355 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 356 PetscFunctionReturn(0); 357 } 358 359 /* 360 Returns a line (ended by a \n, \r or null character of any length. Result should be freed with free() 361 */ 362 static char *Petscgetline(FILE * f) 363 { 364 size_t size = 0; 365 size_t len = 0; 366 size_t last = 0; 367 char *buf = NULL; 368 369 if (feof(f)) return NULL; 370 do { 371 size += 1024; /* BUFSIZ is defined as "the optimal read size for this platform" */ 372 buf = (char*)realloc((void*)buf,size); /* realloc(NULL,n) is the same as malloc(n) */ 373 /* Actually do the read. Note that fgets puts a terminal '\0' on the 374 end of the string, so we make sure we overwrite this */ 375 if (!fgets(buf+len,1024,f)) buf[len]=0; 376 PetscStrlen(buf,&len); 377 last = len - 1; 378 } while (!feof(f) && buf[last] != '\n' && buf[last] != '\r'); 379 if (len) return buf; 380 free(buf); 381 return NULL; 382 } 383 384 /*@C 385 PetscOptionsInsertFile - Inserts options into the database from a file. 386 387 Collective 388 389 Input Parameter: 390 + comm - the processes that will share the options (usually PETSC_COMM_WORLD) 391 . options - options database, use NULL for default global database 392 . file - name of file 393 - require - if PETSC_TRUE will generate an error if the file does not exist 394 395 396 Notes: 397 Use # for lines that are comments and which should be ignored. 398 Usually, instead of using this command, one should list the file name in the call to PetscInitialize(), this insures that certain options 399 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 400 calls to XXXSetFromOptions() it should not be used for options listed under PetscInitialize(). 401 The collectivity of this routine is complex; only the MPI processes in comm will 402 have the affect of these options. If some processes that create objects call this routine and others do 403 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 404 on different ranks. 405 406 Level: developer 407 408 .seealso: PetscOptionsSetValue(), PetscOptionsView(), PetscOptionsHasName(), PetscOptionsGetInt(), 409 PetscOptionsGetReal(), PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(), 410 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 411 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 412 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 413 PetscOptionsFList(), PetscOptionsEList() 414 415 @*/ 416 PetscErrorCode PetscOptionsInsertFile(MPI_Comm comm,PetscOptions options,const char file[],PetscBool require) 417 { 418 char *string,fname[PETSC_MAX_PATH_LEN],*first,*second,*third,*vstring = NULL,*astring = NULL,*packed = NULL; 419 PetscErrorCode ierr; 420 size_t i,len,bytes; 421 FILE *fd; 422 PetscToken token; 423 int err; 424 char cmt[1]={'#'},*cmatch; 425 PetscMPIInt rank,cnt=0,acnt=0,counts[2]; 426 PetscBool isdir; 427 428 PetscFunctionBegin; 429 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 430 if (!rank) { 431 cnt = 0; 432 acnt = 0; 433 434 ierr = PetscFixFilename(file,fname);CHKERRQ(ierr); 435 fd = fopen(fname,"r"); 436 ierr = PetscTestDirectory(fname,'r',&isdir);CHKERRQ(ierr); 437 if (isdir && require) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Specified options file %s is a directory",fname); 438 if (fd && !isdir) { 439 PetscSegBuffer vseg,aseg; 440 ierr = PetscSegBufferCreate(1,4000,&vseg);CHKERRQ(ierr); 441 ierr = PetscSegBufferCreate(1,2000,&aseg);CHKERRQ(ierr); 442 443 /* the following line will not work when opening initial files (like .petscrc) since info is not yet set */ 444 ierr = PetscInfo1(NULL,"Opened options file %s\n",file);CHKERRQ(ierr); 445 446 while ((string = Petscgetline(fd))) { 447 /* eliminate comments from each line */ 448 for (i=0; i<1; i++) { 449 ierr = PetscStrchr(string,cmt[i],&cmatch);CHKERRQ(ierr); 450 if (cmatch) *cmatch = 0; 451 } 452 ierr = PetscStrlen(string,&len);CHKERRQ(ierr); 453 /* replace tabs, ^M, \n with " " */ 454 for (i=0; i<len; i++) { 455 if (string[i] == '\t' || string[i] == '\r' || string[i] == '\n') { 456 string[i] = ' '; 457 } 458 } 459 ierr = PetscTokenCreate(string,' ',&token);CHKERRQ(ierr); 460 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 461 if (!first) { 462 goto destroy; 463 } else if (!first[0]) { /* if first token is empty spaces, redo first token */ 464 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 465 } 466 ierr = PetscTokenFind(token,&second);CHKERRQ(ierr); 467 if (!first) { 468 goto destroy; 469 } else if (first[0] == '-') { 470 ierr = PetscStrlen(first,&len);CHKERRQ(ierr); 471 ierr = PetscSegBufferGet(vseg,len+1,&vstring);CHKERRQ(ierr); 472 ierr = PetscArraycpy(vstring,first,len);CHKERRQ(ierr); 473 vstring[len] = ' '; 474 if (second) { 475 ierr = PetscStrlen(second,&len);CHKERRQ(ierr); 476 ierr = PetscSegBufferGet(vseg,len+3,&vstring);CHKERRQ(ierr); 477 vstring[0] = '"'; 478 ierr = PetscArraycpy(vstring+1,second,len);CHKERRQ(ierr); 479 vstring[len+1] = '"'; 480 vstring[len+2] = ' '; 481 } 482 } else { 483 PetscBool match; 484 485 ierr = PetscStrcasecmp(first,"alias",&match);CHKERRQ(ierr); 486 if (match) { 487 ierr = PetscTokenFind(token,&third);CHKERRQ(ierr); 488 if (!third) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Error in options file:alias missing (%s)",second); 489 ierr = PetscStrlen(second,&len);CHKERRQ(ierr); 490 ierr = PetscSegBufferGet(aseg,len+1,&astring);CHKERRQ(ierr); 491 ierr = PetscArraycpy(astring,second,len);CHKERRQ(ierr); 492 astring[len] = ' '; 493 494 ierr = PetscStrlen(third,&len);CHKERRQ(ierr); 495 ierr = PetscSegBufferGet(aseg,len+1,&astring);CHKERRQ(ierr); 496 ierr = PetscArraycpy(astring,third,len);CHKERRQ(ierr); 497 astring[len] = ' '; 498 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unknown statement in options file: (%s)",string); 499 } 500 destroy: 501 free(string); 502 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 503 } 504 err = fclose(fd); 505 if (err) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fclose() failed on file"); 506 ierr = PetscSegBufferGetSize(aseg,&bytes);CHKERRQ(ierr); /* size without null termination */ 507 ierr = PetscMPIIntCast(bytes,&acnt);CHKERRQ(ierr); 508 ierr = PetscSegBufferGet(aseg,1,&astring);CHKERRQ(ierr); 509 astring[0] = 0; 510 ierr = PetscSegBufferGetSize(vseg,&bytes);CHKERRQ(ierr); /* size without null termination */ 511 ierr = PetscMPIIntCast(bytes,&cnt);CHKERRQ(ierr); 512 ierr = PetscSegBufferGet(vseg,1,&vstring);CHKERRQ(ierr); 513 vstring[0] = 0; 514 ierr = PetscMalloc1(2+acnt+cnt,&packed);CHKERRQ(ierr); 515 ierr = PetscSegBufferExtractTo(aseg,packed);CHKERRQ(ierr); 516 ierr = PetscSegBufferExtractTo(vseg,packed+acnt+1);CHKERRQ(ierr); 517 ierr = PetscSegBufferDestroy(&aseg);CHKERRQ(ierr); 518 ierr = PetscSegBufferDestroy(&vseg);CHKERRQ(ierr); 519 } else if (require) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Unable to open Options File %s",fname); 520 } 521 522 counts[0] = acnt; 523 counts[1] = cnt; 524 ierr = MPI_Bcast(counts,2,MPI_INT,0,comm);CHKERRQ(ierr); 525 acnt = counts[0]; 526 cnt = counts[1]; 527 if (rank) { 528 ierr = PetscMalloc1(2+acnt+cnt,&packed);CHKERRQ(ierr); 529 } 530 if (acnt || cnt) { 531 ierr = MPI_Bcast(packed,2+acnt+cnt,MPI_CHAR,0,comm);CHKERRQ(ierr); 532 astring = packed; 533 vstring = packed + acnt + 1; 534 } 535 536 if (acnt) { 537 PetscToken token; 538 char *first,*second; 539 540 ierr = PetscTokenCreate(astring,' ',&token);CHKERRQ(ierr); 541 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 542 while (first) { 543 ierr = PetscTokenFind(token,&second);CHKERRQ(ierr); 544 ierr = PetscOptionsSetAlias(options,first,second);CHKERRQ(ierr); 545 ierr = PetscTokenFind(token,&first);CHKERRQ(ierr); 546 } 547 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 548 } 549 550 if (cnt) { 551 ierr = PetscOptionsInsertString(options,vstring);CHKERRQ(ierr); 552 } 553 ierr = PetscFree(packed);CHKERRQ(ierr); 554 PetscFunctionReturn(0); 555 } 556 557 static PetscErrorCode PetscOptionsInsertArgs(PetscOptions options,int argc,char *args[]) 558 { 559 PetscErrorCode ierr; 560 int left = argc - 1; 561 char **eargs = args + 1; 562 563 PetscFunctionBegin; 564 while (left) { 565 PetscBool isoptions_file,isprefixpush,isprefixpop,isp4,tisp4,isp4yourname,isp4rmrank,key; 566 ierr = PetscStrcasecmp(eargs[0],"-options_file",&isoptions_file);CHKERRQ(ierr); 567 ierr = PetscStrcasecmp(eargs[0],"-prefix_push",&isprefixpush);CHKERRQ(ierr); 568 ierr = PetscStrcasecmp(eargs[0],"-prefix_pop",&isprefixpop);CHKERRQ(ierr); 569 ierr = PetscStrcasecmp(eargs[0],"-p4pg",&isp4);CHKERRQ(ierr); 570 ierr = PetscStrcasecmp(eargs[0],"-p4yourname",&isp4yourname);CHKERRQ(ierr); 571 ierr = PetscStrcasecmp(eargs[0],"-p4rmrank",&isp4rmrank);CHKERRQ(ierr); 572 ierr = PetscStrcasecmp(eargs[0],"-p4wd",&tisp4);CHKERRQ(ierr); 573 isp4 = (PetscBool) (isp4 || tisp4); 574 ierr = PetscStrcasecmp(eargs[0],"-np",&tisp4);CHKERRQ(ierr); 575 isp4 = (PetscBool) (isp4 || tisp4); 576 ierr = PetscStrcasecmp(eargs[0],"-p4amslave",&tisp4);CHKERRQ(ierr); 577 ierr = PetscOptionsValidKey(eargs[0],&key);CHKERRQ(ierr); 578 579 if (!key) { 580 eargs++; left--; 581 } else if (isoptions_file) { 582 if (left <= 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing filename for -options_file filename option"); 583 if (eargs[1][0] == '-') SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing filename for -options_file filename option"); 584 ierr = PetscOptionsInsertFile(PETSC_COMM_WORLD,options,eargs[1],PETSC_TRUE);CHKERRQ(ierr); 585 eargs += 2; left -= 2; 586 } else if (isprefixpush) { 587 if (left <= 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing prefix for -prefix_push option"); 588 if (eargs[1][0] == '-') SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Missing prefix for -prefix_push option (prefixes cannot start with '-')"); 589 ierr = PetscOptionsPrefixPush(options,eargs[1]);CHKERRQ(ierr); 590 eargs += 2; left -= 2; 591 } else if (isprefixpop) { 592 ierr = PetscOptionsPrefixPop(options);CHKERRQ(ierr); 593 eargs++; left--; 594 595 /* 596 These are "bad" options that MPICH, etc put on the command line 597 we strip them out here. 598 */ 599 } else if (tisp4 || isp4rmrank) { 600 eargs += 1; left -= 1; 601 } else if (isp4 || isp4yourname) { 602 eargs += 2; left -= 2; 603 } else { 604 PetscBool nextiskey = PETSC_FALSE; 605 if (left >= 2) {ierr = PetscOptionsValidKey(eargs[1],&nextiskey);CHKERRQ(ierr);} 606 if (left < 2 || nextiskey) { 607 ierr = PetscOptionsSetValue(options,eargs[0],NULL);CHKERRQ(ierr); 608 eargs++; left--; 609 } else { 610 ierr = PetscOptionsSetValue(options,eargs[0],eargs[1]);CHKERRQ(ierr); 611 eargs += 2; left -= 2; 612 } 613 } 614 } 615 PetscFunctionReturn(0); 616 } 617 618 619 /*@C 620 PetscOptionsInsert - Inserts into the options database from the command line, 621 the environmental variable and a file. 622 623 Collective on PETSC_COMM_WORLD 624 625 Input Parameters: 626 + options - options database or NULL for the default global database 627 . argc - count of number of command line arguments 628 . args - the command line arguments 629 - file - optional filename, defaults to ~username/.petscrc 630 631 Note: 632 Since PetscOptionsInsert() is automatically called by PetscInitialize(), 633 the user does not typically need to call this routine. PetscOptionsInsert() 634 can be called several times, adding additional entries into the database. 635 636 Options Database Keys: 637 + -options_monitor <optional filename> - print options names and values as they are set 638 - -options_file <filename> - read options from a file 639 640 Level: advanced 641 642 .seealso: PetscOptionsDestroy(), PetscOptionsView(), PetscOptionsInsertString(), PetscOptionsInsertFile(), 643 PetscInitialize() 644 @*/ 645 PetscErrorCode PetscOptionsInsert(PetscOptions options,int *argc,char ***args,const char file[]) 646 { 647 PetscErrorCode ierr; 648 PetscMPIInt rank; 649 char filename[PETSC_MAX_PATH_LEN]; 650 PetscBool flag = PETSC_FALSE; 651 652 653 PetscFunctionBegin; 654 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 655 656 if (file && file[0]) { 657 ierr = PetscStrreplace(PETSC_COMM_WORLD,file,filename,PETSC_MAX_PATH_LEN);CHKERRQ(ierr); 658 ierr = PetscOptionsInsertFile(PETSC_COMM_WORLD,options,filename,PETSC_TRUE);CHKERRQ(ierr); 659 } 660 /* 661 We want to be able to give -skip_petscrc on the command line, but need to parse it first. Since the command line 662 should take precedence, we insert it twice. It would be sufficient to just scan for -skip_petscrc. 663 */ 664 if (argc && args && *argc) {ierr = PetscOptionsInsertArgs(options,*argc,*args);CHKERRQ(ierr);} 665 ierr = PetscOptionsGetBool(NULL,NULL,"-skip_petscrc",&flag,NULL);CHKERRQ(ierr); 666 if (!flag) { 667 ierr = PetscGetHomeDirectory(filename,PETSC_MAX_PATH_LEN-16);CHKERRQ(ierr); 668 /* PetscOptionsInsertFile() does a fopen() on rank0 only - so only rank0 HomeDir value is relavent */ 669 if (filename[0]) { ierr = PetscStrcat(filename,"/.petscrc");CHKERRQ(ierr); } 670 ierr = PetscOptionsInsertFile(PETSC_COMM_WORLD,options,filename,PETSC_FALSE);CHKERRQ(ierr); 671 ierr = PetscOptionsInsertFile(PETSC_COMM_WORLD,options,".petscrc",PETSC_FALSE);CHKERRQ(ierr); 672 ierr = PetscOptionsInsertFile(PETSC_COMM_WORLD,options,"petscrc",PETSC_FALSE);CHKERRQ(ierr); 673 } 674 675 /* insert environment options */ 676 { 677 char *eoptions = NULL; 678 size_t len = 0; 679 if (!rank) { 680 eoptions = (char*)getenv("PETSC_OPTIONS"); 681 ierr = PetscStrlen(eoptions,&len);CHKERRQ(ierr); 682 ierr = MPI_Bcast(&len,1,MPIU_SIZE_T,0,PETSC_COMM_WORLD);CHKERRQ(ierr); 683 } else { 684 ierr = MPI_Bcast(&len,1,MPIU_SIZE_T,0,PETSC_COMM_WORLD);CHKERRQ(ierr); 685 if (len) { 686 ierr = PetscMalloc1(len+1,&eoptions);CHKERRQ(ierr); 687 } 688 } 689 if (len) { 690 ierr = MPI_Bcast(eoptions,len,MPI_CHAR,0,PETSC_COMM_WORLD);CHKERRQ(ierr); 691 if (rank) eoptions[len] = 0; 692 ierr = PetscOptionsInsertString(options,eoptions);CHKERRQ(ierr); 693 if (rank) {ierr = PetscFree(eoptions);CHKERRQ(ierr);} 694 } 695 } 696 697 #if defined(PETSC_HAVE_YAML) 698 { 699 char yaml_file[PETSC_MAX_PATH_LEN]; 700 PetscBool yaml_flg; 701 ierr = PetscOptionsGetString(NULL,NULL,"-options_file_yaml",yaml_file,PETSC_MAX_PATH_LEN,&yaml_flg);CHKERRQ(ierr); 702 if (yaml_flg) { 703 ierr = PetscOptionsInsertFileYAML(PETSC_COMM_WORLD,yaml_file,PETSC_TRUE);CHKERRQ(ierr); 704 } 705 } 706 #endif 707 708 /* insert command line options again because they take precedence over arguments in petscrc/environment */ 709 if (argc && args && *argc) {ierr = PetscOptionsInsertArgs(options,*argc,*args);CHKERRQ(ierr);} 710 PetscFunctionReturn(0); 711 } 712 713 /*@C 714 PetscOptionsView - Prints the options that have been loaded. This is 715 useful for debugging purposes. 716 717 Logically Collective on PetscViewer 718 719 Input Parameter: 720 + options - options database, use NULL for default global database 721 - viewer - must be an PETSCVIEWERASCII viewer 722 723 Options Database Key: 724 . -options_view - Activates PetscOptionsView() within PetscFinalize() 725 726 Notes: 727 Only the rank zero process of MPI_Comm used to create view prints the option values. Other processes 728 may have different values but they are not printed. 729 730 Level: advanced 731 732 .seealso: PetscOptionsAllUsed() 733 @*/ 734 PetscErrorCode PetscOptionsView(PetscOptions options,PetscViewer viewer) 735 { 736 PetscErrorCode ierr; 737 PetscInt i; 738 PetscBool isascii; 739 740 PetscFunctionBegin; 741 if (viewer) PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); 742 options = options ? options : defaultoptions; 743 if (!viewer) viewer = PETSC_VIEWER_STDOUT_WORLD; 744 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr); 745 if (!isascii) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Only supports ASCII viewer"); 746 747 if (!options->N) { 748 ierr = PetscViewerASCIIPrintf(viewer,"#No PETSc Option Table entries\n");CHKERRQ(ierr); 749 PetscFunctionReturn(0); 750 } 751 752 ierr = PetscViewerASCIIPrintf(viewer,"#PETSc Option Table entries:\n");CHKERRQ(ierr); 753 for (i=0; i<options->N; i++) { 754 if (options->values[i]) { 755 ierr = PetscViewerASCIIPrintf(viewer,"-%s %s\n",options->names[i],options->values[i]);CHKERRQ(ierr); 756 } else { 757 ierr = PetscViewerASCIIPrintf(viewer,"-%s\n",options->names[i]);CHKERRQ(ierr); 758 } 759 } 760 ierr = PetscViewerASCIIPrintf(viewer,"#End of PETSc Option Table entries\n");CHKERRQ(ierr); 761 PetscFunctionReturn(0); 762 } 763 764 /* 765 Called by error handlers to print options used in run 766 */ 767 PETSC_EXTERN PetscErrorCode PetscOptionsViewError(void) 768 { 769 PetscInt i; 770 PetscOptions options = defaultoptions; 771 772 PetscFunctionBegin; 773 if (options->N) { 774 (*PetscErrorPrintf)("PETSc Option Table entries:\n"); 775 } else { 776 (*PetscErrorPrintf)("No PETSc Option Table entries\n"); 777 } 778 for (i=0; i<options->N; i++) { 779 if (options->values[i]) { 780 (*PetscErrorPrintf)("-%s %s\n",options->names[i],options->values[i]); 781 } else { 782 (*PetscErrorPrintf)("-%s\n",options->names[i]); 783 } 784 } 785 PetscFunctionReturn(0); 786 } 787 788 /*@C 789 PetscOptionsPrefixPush - Designate a prefix to be used by all options insertions to follow. 790 791 Logically Collective 792 793 Input Parameter: 794 + options - options database, or NULL for the default global database 795 - prefix - The string to append to the existing prefix 796 797 Options Database Keys: 798 + -prefix_push <some_prefix_> - push the given prefix 799 - -prefix_pop - pop the last prefix 800 801 Notes: 802 It is common to use this in conjunction with -options_file as in 803 804 $ -prefix_push system1_ -options_file system1rc -prefix_pop -prefix_push system2_ -options_file system2rc -prefix_pop 805 806 where the files no longer require all options to be prefixed with -system2_. 807 808 The collectivity of this routine is complex; only the MPI processes that call this routine will 809 have the affect of these options. If some processes that create objects call this routine and others do 810 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 811 on different ranks. 812 813 Level: advanced 814 815 .seealso: PetscOptionsPrefixPop(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsSetValue() 816 @*/ 817 PetscErrorCode PetscOptionsPrefixPush(PetscOptions options,const char prefix[]) 818 { 819 PetscErrorCode ierr; 820 size_t n; 821 PetscInt start; 822 char key[MAXOPTNAME+1]; 823 PetscBool valid; 824 825 PetscFunctionBegin; 826 PetscValidCharPointer(prefix,1); 827 options = options ? options : defaultoptions; 828 if (options->prefixind >= MAXPREFIXES) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Maximum depth of prefix stack %d exceeded, recompile \n src/sys/objects/options.c with larger value for MAXPREFIXES",MAXPREFIXES); 829 key[0] = '-'; /* keys must start with '-' */ 830 ierr = PetscStrncpy(key+1,prefix,sizeof(key)-1);CHKERRQ(ierr); 831 ierr = PetscOptionsValidKey(key,&valid);CHKERRQ(ierr); 832 if (!valid) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Given prefix \"%s\" not valid (the first character must be a letter, do not include leading '-')",prefix); 833 start = options->prefixind ? options->prefixstack[options->prefixind-1] : 0; 834 ierr = PetscStrlen(prefix,&n);CHKERRQ(ierr); 835 if (n+1 > sizeof(options->prefix)-start) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Maximum prefix length %d exceeded",sizeof(options->prefix)); 836 ierr = PetscArraycpy(options->prefix+start,prefix,n+1);CHKERRQ(ierr); 837 options->prefixstack[options->prefixind++] = start+n; 838 PetscFunctionReturn(0); 839 } 840 841 /*@C 842 PetscOptionsPrefixPop - Remove the latest options prefix, see PetscOptionsPrefixPush() for details 843 844 Logically Collective on the MPI_Comm that called PetscOptionsPrefixPush() 845 846 Input Parameters: 847 . options - options database, or NULL for the default global database 848 849 Level: advanced 850 851 .seealso: PetscOptionsPrefixPush(), PetscOptionsPush(), PetscOptionsPop(), PetscOptionsCreate(), PetscOptionsSetValue() 852 @*/ 853 PetscErrorCode PetscOptionsPrefixPop(PetscOptions options) 854 { 855 PetscInt offset; 856 857 PetscFunctionBegin; 858 options = options ? options : defaultoptions; 859 if (options->prefixind < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"More prefixes popped than pushed"); 860 options->prefixind--; 861 offset = options->prefixind ? options->prefixstack[options->prefixind-1] : 0; 862 options->prefix[offset] = 0; 863 PetscFunctionReturn(0); 864 } 865 866 /*@C 867 PetscOptionsClear - Removes all options form the database leaving it empty. 868 869 Logically Collective 870 871 Input Parameters: 872 . options - options database, use NULL for the default global database 873 874 The collectivity of this routine is complex; only the MPI processes that call this routine will 875 have the affect of these options. If some processes that create objects call this routine and others do 876 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 877 on different ranks. 878 879 Level: developer 880 881 .seealso: PetscOptionsInsert() 882 @*/ 883 PetscErrorCode PetscOptionsClear(PetscOptions options) 884 { 885 PetscInt i; 886 887 options = options ? options : defaultoptions; 888 if (!options) return 0; 889 890 for (i=0; i<options->N; i++) { 891 if (options->names[i]) free(options->names[i]); 892 if (options->values[i]) free(options->values[i]); 893 } 894 options->N = 0; 895 896 for (i=0; i<options->Naliases; i++) { 897 free(options->aliases1[i]); 898 free(options->aliases2[i]); 899 } 900 options->Naliases = 0; 901 902 /* destroy hash table */ 903 kh_destroy(HO,options->ht); 904 options->ht = NULL; 905 906 options->prefixind = 0; 907 options->prefix[0] = 0; 908 options->help = PETSC_FALSE; 909 return 0; 910 } 911 912 /*@C 913 PetscOptionsSetAlias - Makes a key and alias for another key 914 915 Logically Collective 916 917 Input Parameters: 918 + options - options database, or NULL for default global database 919 . newname - the alias 920 - oldname - the name that alias will refer to 921 922 Level: advanced 923 924 The collectivity of this routine is complex; only the MPI processes that call this routine will 925 have the affect of these options. If some processes that create objects call this routine and others do 926 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 927 on different ranks. 928 929 .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),OptionsHasName(), 930 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(),PetscOptionsBool(), 931 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 932 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 933 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 934 PetscOptionsFList(), PetscOptionsEList() 935 @*/ 936 PetscErrorCode PetscOptionsSetAlias(PetscOptions options,const char newname[],const char oldname[]) 937 { 938 PetscInt n; 939 size_t len; 940 PetscErrorCode ierr; 941 942 PetscFunctionBegin; 943 PetscValidCharPointer(newname,2); 944 PetscValidCharPointer(oldname,3); 945 options = options ? options : defaultoptions; 946 if (newname[0] != '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"aliased must start with '-': Instead %s",newname); 947 if (oldname[0] != '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"aliasee must start with '-': Instead %s",oldname); 948 949 n = options->Naliases; 950 if (n >= MAXALIASES) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_MEM,"You have defined to many PETSc options aliases, limit %d recompile \n src/sys/objects/options.c with larger value for MAXALIASES",MAXALIASES); 951 952 newname++; oldname++; 953 ierr = PetscStrlen(newname,&len);CHKERRQ(ierr); 954 options->aliases1[n] = (char*)malloc((len+1)*sizeof(char)); 955 ierr = PetscStrcpy(options->aliases1[n],newname);CHKERRQ(ierr); 956 ierr = PetscStrlen(oldname,&len);CHKERRQ(ierr); 957 options->aliases2[n] = (char*)malloc((len+1)*sizeof(char)); 958 ierr = PetscStrcpy(options->aliases2[n],oldname);CHKERRQ(ierr); 959 options->Naliases++; 960 PetscFunctionReturn(0); 961 } 962 963 /*@C 964 PetscOptionsSetValue - Sets an option name-value pair in the options 965 database, overriding whatever is already present. 966 967 Logically Collective 968 969 Input Parameters: 970 + options - options database, use NULL for the default global database 971 . name - name of option, this SHOULD have the - prepended 972 - value - the option value (not used for all options, so can be NULL) 973 974 Level: intermediate 975 976 Note: 977 This function can be called BEFORE PetscInitialize() 978 979 The collectivity of this routine is complex; only the MPI processes that call this routine will 980 have the affect of these options. If some processes that create objects call this routine and others do 981 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 982 on different ranks. 983 984 Developers Note: Uses malloc() directly because PETSc may not be initialized yet. 985 986 .seealso: PetscOptionsInsert(), PetscOptionsClearValue() 987 @*/ 988 PetscErrorCode PetscOptionsSetValue(PetscOptions options,const char name[],const char value[]) 989 { 990 size_t len; 991 int N,n,i; 992 char **names; 993 char fullname[MAXOPTNAME] = ""; 994 PetscErrorCode ierr; 995 996 if (!options && !defaultoptions) { 997 ierr = PetscOptionsCreateDefault();if (ierr) return ierr; 998 } 999 options = options ? options : defaultoptions; 1000 1001 if (name[0] != '-') return PETSC_ERR_ARG_OUTOFRANGE; 1002 1003 /* this is so that -h and -help are equivalent (p4 does not like -help)*/ 1004 if (!strcmp(name,"-h")) name = "-help"; 1005 if (!PetscOptNameCmp(name,"-help")) options->help = PETSC_TRUE; 1006 1007 name++; /* skip starting dash */ 1008 1009 if (options->prefixind > 0) { 1010 strncpy(fullname,options->prefix,sizeof(fullname)); 1011 fullname[sizeof(fullname)-1] = 0; 1012 strncat(fullname,name,sizeof(fullname)-strlen(fullname)-1); 1013 fullname[sizeof(fullname)-1] = 0; 1014 name = fullname; 1015 } 1016 1017 /* check against aliases */ 1018 N = options->Naliases; 1019 for (i=0; i<N; i++) { 1020 int result = PetscOptNameCmp(options->aliases1[i],name); 1021 if (!result) { name = options->aliases2[i]; break; } 1022 } 1023 1024 /* slow search */ 1025 N = n = options->N; 1026 names = options->names; 1027 for (i=0; i<N; i++) { 1028 int result = PetscOptNameCmp(names[i],name); 1029 if (!result) { 1030 n = i; goto setvalue; 1031 } else if (result > 0) { 1032 n = i; break; 1033 } 1034 } 1035 if (N >= MAXOPTIONS) return PETSC_ERR_MEM; 1036 /* shift remaining values up 1 */ 1037 for (i=N; i>n; i--) { 1038 options->names[i] = options->names[i-1]; 1039 options->values[i] = options->values[i-1]; 1040 options->used[i] = options->used[i-1]; 1041 } 1042 options->names[n] = NULL; 1043 options->values[n] = NULL; 1044 options->used[n] = PETSC_FALSE; 1045 options->N++; 1046 1047 /* destroy hash table */ 1048 kh_destroy(HO,options->ht); 1049 options->ht = NULL; 1050 1051 /* set new name */ 1052 len = strlen(name); 1053 options->names[n] = (char*)malloc((len+1)*sizeof(char)); 1054 if (!options->names[n]) return PETSC_ERR_MEM; 1055 strcpy(options->names[n],name); 1056 1057 setvalue: 1058 /* set new value */ 1059 if (options->values[n]) free(options->values[n]); 1060 len = value ? strlen(value) : 0; 1061 if (len) { 1062 options->values[n] = (char*)malloc((len+1)*sizeof(char)); 1063 if (!options->values[n]) return PETSC_ERR_MEM; 1064 strcpy(options->values[n],value); 1065 } else { 1066 options->values[n] = NULL; 1067 } 1068 1069 ierr = PetscOptionsMonitor(options,name,value?value:"");if (ierr) return ierr; 1070 return 0; 1071 } 1072 1073 /*@C 1074 PetscOptionsClearValue - Clears an option name-value pair in the options 1075 database, overriding whatever is already present. 1076 1077 Logically Collective 1078 1079 Input Parameter: 1080 + options - options database, use NULL for the default global database 1081 - name - name of option, this SHOULD have the - prepended 1082 1083 Level: intermediate 1084 1085 The collectivity of this routine is complex; only the MPI processes that call this routine will 1086 have the affect of these options. If some processes that create objects call this routine and others do 1087 not the code may fail in complicated ways because the same parallel solvers may incorrectly use different options 1088 on different ranks. 1089 1090 .seealso: PetscOptionsInsert() 1091 @*/ 1092 PetscErrorCode PetscOptionsClearValue(PetscOptions options,const char name[]) 1093 { 1094 int N,n,i; 1095 char **names; 1096 PetscErrorCode ierr; 1097 1098 PetscFunctionBegin; 1099 options = options ? options : defaultoptions; 1100 if (name[0] != '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Name must begin with '-': Instead %s",name); 1101 1102 /* this is so that -h and -help are equivalent (p4 does not like -help)*/ 1103 if (!strcmp(name,"-h")) name = "-help"; 1104 if (!PetscOptNameCmp(name,"-help")) options->help = PETSC_FALSE; 1105 1106 name++; /* skip starting dash */ 1107 1108 /* slow search */ 1109 N = n = options->N; 1110 names = options->names; 1111 for (i=0; i<N; i++) { 1112 int result = PetscOptNameCmp(names[i],name); 1113 if (!result) { 1114 n = i; break; 1115 } else if (result > 0) { 1116 n = N; break; 1117 } 1118 } 1119 if (n == N) PetscFunctionReturn(0); /* it was not present */ 1120 1121 /* remove name and value */ 1122 if (options->names[n]) free(options->names[n]); 1123 if (options->values[n]) free(options->values[n]); 1124 /* shift remaining values down 1 */ 1125 for (i=n; i<N-1; i++) { 1126 options->names[i] = options->names[i+1]; 1127 options->values[i] = options->values[i+1]; 1128 options->used[i] = options->used[i+1]; 1129 } 1130 options->N--; 1131 1132 /* destroy hash table */ 1133 kh_destroy(HO,options->ht); 1134 options->ht = NULL; 1135 1136 ierr = PetscOptionsMonitor(options,name,NULL);CHKERRQ(ierr); 1137 PetscFunctionReturn(0); 1138 } 1139 1140 /*@C 1141 PetscOptionsFindPair - Gets an option name-value pair from the options database. 1142 1143 Not Collective 1144 1145 Input Parameters: 1146 + options - options database, use NULL for the default global database 1147 . pre - the string to prepend to the name or NULL, this SHOULD NOT have the "-" prepended 1148 - name - name of option, this SHOULD have the "-" prepended 1149 1150 Output Parameters: 1151 + value - the option value (optional, not used for all options) 1152 - set - whether the option is set (optional) 1153 1154 Notes: 1155 Each process may find different values or no value depending on how options were inserted into the database 1156 1157 Level: developer 1158 1159 .seealso: PetscOptionsSetValue(), PetscOptionsClearValue() 1160 @*/ 1161 PetscErrorCode PetscOptionsFindPair(PetscOptions options,const char pre[],const char name[],const char *value[],PetscBool *set) 1162 { 1163 char buf[MAXOPTNAME]; 1164 PetscBool usehashtable = PETSC_TRUE; 1165 PetscBool matchnumbers = PETSC_TRUE; 1166 PetscErrorCode ierr; 1167 1168 PetscFunctionBegin; 1169 options = options ? options : defaultoptions; 1170 if (pre && PetscUnlikely(pre[0] == '-')) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Prefix cannot begin with '-': Instead %s",pre); 1171 if (PetscUnlikely(name[0] != '-')) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Name must begin with '-': Instead %s",name); 1172 1173 name++; /* skip starting dash */ 1174 1175 /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */ 1176 if (pre && pre[0]) { 1177 char *ptr = buf; 1178 if (name[0] == '-') { *ptr++ = '-'; name++; } 1179 ierr = PetscStrncpy(ptr,pre,buf+sizeof(buf)-ptr);CHKERRQ(ierr); 1180 ierr = PetscStrlcat(buf,name,sizeof(buf));CHKERRQ(ierr); 1181 name = buf; 1182 } 1183 1184 if (PetscDefined(USE_DEBUG)) { 1185 PetscBool valid; 1186 char key[MAXOPTNAME+1] = "-"; 1187 ierr = PetscStrncpy(key+1,name,sizeof(key)-1);CHKERRQ(ierr); 1188 ierr = PetscOptionsValidKey(key,&valid);CHKERRQ(ierr); 1189 if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid option '%s' obtained from pre='%s' and name='%s'",key,pre?pre:"",name); 1190 } 1191 1192 if (!options->ht && usehashtable) { 1193 int i,ret; 1194 khiter_t it; 1195 khash_t(HO) *ht; 1196 ht = kh_init(HO); 1197 if (PetscUnlikely(!ht)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_MEM,"Hash table allocation failed"); 1198 ret = kh_resize(HO,ht,options->N*2); /* twice the required size to reduce risk of collisions */ 1199 if (PetscUnlikely(ret)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_MEM,"Hash table allocation failed"); 1200 for (i=0; i<options->N; i++) { 1201 it = kh_put(HO,ht,options->names[i],&ret); 1202 if (PetscUnlikely(ret != 1)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_MEM,"Hash table allocation failed"); 1203 kh_val(ht,it) = i; 1204 } 1205 options->ht = ht; 1206 } 1207 1208 if (usehashtable) 1209 { /* fast search */ 1210 khash_t(HO) *ht = options->ht; 1211 khiter_t it = kh_get(HO,ht,name); 1212 if (it != kh_end(ht)) { 1213 int i = kh_val(ht,it); 1214 options->used[i] = PETSC_TRUE; 1215 if (value) *value = options->values[i]; 1216 if (set) *set = PETSC_TRUE; 1217 PetscFunctionReturn(0); 1218 } 1219 } else 1220 { /* slow search */ 1221 int i, N = options->N; 1222 for (i=0; i<N; i++) { 1223 int result = PetscOptNameCmp(options->names[i],name); 1224 if (!result) { 1225 options->used[i] = PETSC_TRUE; 1226 if (value) *value = options->values[i]; 1227 if (set) *set = PETSC_TRUE; 1228 PetscFunctionReturn(0); 1229 } else if (result > 0) { 1230 break; 1231 } 1232 } 1233 } 1234 1235 /* 1236 The following block slows down all lookups in the most frequent path (most lookups are unsuccessful). 1237 Maybe this special lookup mode should be enabled on request with a push/pop API. 1238 The feature of matching _%d_ used sparingly in the codebase. 1239 */ 1240 if (matchnumbers) { 1241 int i,j,cnt = 0,locs[16],loce[16]; 1242 /* determine the location and number of all _%d_ in the key */ 1243 for (i=0; name[i]; i++) { 1244 if (name[i] == '_') { 1245 for (j=i+1; name[j]; j++) { 1246 if (name[j] >= '0' && name[j] <= '9') continue; 1247 if (name[j] == '_' && j > i+1) { /* found a number */ 1248 locs[cnt] = i+1; 1249 loce[cnt++] = j+1; 1250 } 1251 i = j-1; 1252 break; 1253 } 1254 } 1255 } 1256 for (i=0; i<cnt; i++) { 1257 PetscBool found; 1258 char opt[MAXOPTNAME+1] = "-", tmp[MAXOPTNAME]; 1259 ierr = PetscStrncpy(tmp,name,PetscMin((size_t)(locs[i]+1),sizeof(tmp)));CHKERRQ(ierr); 1260 ierr = PetscStrlcat(opt,tmp,sizeof(opt));CHKERRQ(ierr); 1261 ierr = PetscStrlcat(opt,name+loce[i],sizeof(opt));CHKERRQ(ierr); 1262 ierr = PetscOptionsFindPair(options,NULL,opt,value,&found);CHKERRQ(ierr); 1263 if (found) {if (set) *set = PETSC_TRUE; PetscFunctionReturn(0);} 1264 } 1265 } 1266 1267 if (set) *set = PETSC_FALSE; 1268 PetscFunctionReturn(0); 1269 } 1270 1271 /* Check whether any option begins with pre+name */ 1272 PETSC_EXTERN PetscErrorCode PetscOptionsFindPairPrefix_Private(PetscOptions options,const char pre[], const char name[],const char *value[],PetscBool *set) 1273 { 1274 char buf[MAXOPTNAME]; 1275 int numCnt = 0, locs[16],loce[16]; 1276 PetscErrorCode ierr; 1277 1278 PetscFunctionBegin; 1279 options = options ? options : defaultoptions; 1280 if (pre && pre[0] == '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Prefix cannot begin with '-': Instead %s",pre); 1281 if (name[0] != '-') SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Name must begin with '-': Instead %s",name); 1282 1283 name++; /* skip starting dash */ 1284 1285 /* append prefix to name, if prefix="foo_" and option='--bar", prefixed option is --foo_bar */ 1286 if (pre && pre[0]) { 1287 char *ptr = buf; 1288 if (name[0] == '-') { *ptr++ = '-'; name++; } 1289 ierr = PetscStrncpy(ptr,pre,sizeof(buf)+(size_t)(ptr-buf));CHKERRQ(ierr); 1290 ierr = PetscStrlcat(buf,name,sizeof(buf));CHKERRQ(ierr); 1291 name = buf; 1292 } 1293 1294 if (PetscDefined(USE_DEBUG)) { 1295 PetscBool valid; 1296 char key[MAXOPTNAME+1] = "-"; 1297 ierr = PetscStrncpy(key+1,name,sizeof(key)-1);CHKERRQ(ierr); 1298 ierr = PetscOptionsValidKey(key,&valid);CHKERRQ(ierr); 1299 if (!valid) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid option '%s' obtained from pre='%s' and name='%s'",key,pre?pre:"",name); 1300 } 1301 1302 /* determine the location and number of all _%d_ in the key */ 1303 { 1304 int i,j; 1305 for (i=0; name[i]; i++) { 1306 if (name[i] == '_') { 1307 for (j=i+1; name[j]; j++) { 1308 if (name[j] >= '0' && name[j] <= '9') continue; 1309 if (name[j] == '_' && j > i+1) { /* found a number */ 1310 locs[numCnt] = i+1; 1311 loce[numCnt++] = j+1; 1312 } 1313 i = j-1; 1314 break; 1315 } 1316 } 1317 } 1318 } 1319 1320 { /* slow search */ 1321 int c, i; 1322 size_t len; 1323 PetscBool match; 1324 1325 for (c = -1; c < numCnt; ++c) { 1326 char opt[MAXOPTNAME+1] = "", tmp[MAXOPTNAME]; 1327 1328 if (c < 0) { 1329 ierr = PetscStrcpy(opt,name);CHKERRQ(ierr); 1330 } else { 1331 ierr = PetscStrncpy(tmp,name,PetscMin((size_t)(locs[c]+1),sizeof(tmp)));CHKERRQ(ierr); 1332 ierr = PetscStrlcat(opt,tmp,sizeof(opt));CHKERRQ(ierr); 1333 ierr = PetscStrlcat(opt,name+loce[c],sizeof(opt));CHKERRQ(ierr); 1334 } 1335 ierr = PetscStrlen(opt,&len);CHKERRQ(ierr); 1336 for (i=0; i<options->N; i++) { 1337 ierr = PetscStrncmp(options->names[i],opt,len,&match);CHKERRQ(ierr); 1338 if (match) { 1339 options->used[i] = PETSC_TRUE; 1340 if (value) *value = options->values[i]; 1341 if (set) *set = PETSC_TRUE; 1342 PetscFunctionReturn(0); 1343 } 1344 } 1345 } 1346 } 1347 1348 if (set) *set = PETSC_FALSE; 1349 PetscFunctionReturn(0); 1350 } 1351 1352 /*@C 1353 PetscOptionsReject - Generates an error if a certain option is given. 1354 1355 Not Collective 1356 1357 Input Parameters: 1358 + options - options database, use NULL for default global database 1359 . pre - the option prefix (may be NULL) 1360 . name - the option name one is seeking 1361 - mess - error message (may be NULL) 1362 1363 Level: advanced 1364 1365 .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(),OptionsHasName(), 1366 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(), 1367 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 1368 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 1369 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 1370 PetscOptionsFList(), PetscOptionsEList() 1371 @*/ 1372 PetscErrorCode PetscOptionsReject(PetscOptions options,const char pre[],const char name[],const char mess[]) 1373 { 1374 PetscErrorCode ierr; 1375 PetscBool flag = PETSC_FALSE; 1376 1377 PetscFunctionBegin; 1378 ierr = PetscOptionsHasName(options,pre,name,&flag);CHKERRQ(ierr); 1379 if (flag) { 1380 if (mess && mess[0]) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Program has disabled option: -%s%s with %s",pre?pre:"",name+1,mess); 1381 else SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Program has disabled option: -%s%s",pre?pre:"",name+1); 1382 } 1383 PetscFunctionReturn(0); 1384 } 1385 1386 /*@C 1387 PetscOptionsHasHelp - Determines whether the "-help" option is in the database. 1388 1389 Not Collective 1390 1391 Input Parameters: 1392 . options - options database, use NULL for default global database 1393 1394 Output Parameters: 1395 . set - PETSC_TRUE if found else PETSC_FALSE. 1396 1397 Level: advanced 1398 1399 .seealso: PetscOptionsHasName() 1400 @*/ 1401 PetscErrorCode PetscOptionsHasHelp(PetscOptions options,PetscBool *set) 1402 { 1403 PetscFunctionBegin; 1404 PetscValidPointer(set,2); 1405 options = options ? options : defaultoptions; 1406 *set = options->help; 1407 PetscFunctionReturn(0); 1408 } 1409 1410 /*@C 1411 PetscOptionsHasName - Determines whether a certain option is given in the database. This returns true whether the option is a number, string or boolean, even 1412 its value is set to false. 1413 1414 Not Collective 1415 1416 Input Parameters: 1417 + options - options database, use NULL for default global database 1418 . pre - string to prepend to the name or NULL 1419 - name - the option one is seeking 1420 1421 Output Parameters: 1422 . set - PETSC_TRUE if found else PETSC_FALSE. 1423 1424 Level: beginner 1425 1426 Notes: 1427 Name cannot be simply "-h". 1428 1429 In many cases you probably want to use PetscOptionsGetBool() instead of calling this, to allowing toggling values. 1430 1431 .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(), 1432 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(), 1433 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 1434 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 1435 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 1436 PetscOptionsFList(), PetscOptionsEList() 1437 @*/ 1438 PetscErrorCode PetscOptionsHasName(PetscOptions options,const char pre[],const char name[],PetscBool *set) 1439 { 1440 const char *value; 1441 PetscErrorCode ierr; 1442 PetscBool flag; 1443 1444 PetscFunctionBegin; 1445 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);CHKERRQ(ierr); 1446 if (set) *set = flag; 1447 PetscFunctionReturn(0); 1448 } 1449 1450 /*@C 1451 PetscOptionsGetAll - Lists all the options the program was run with in a single string. 1452 1453 Not Collective 1454 1455 Input Parameter: 1456 . options - the options database, use NULL for the default global database 1457 1458 Output Parameter: 1459 . copts - pointer where string pointer is stored 1460 1461 Notes: 1462 The array and each entry in the array should be freed with PetscFree() 1463 Each process may have different values depending on how the options were inserted into the database 1464 1465 Level: advanced 1466 1467 .seealso: PetscOptionsAllUsed(), PetscOptionsView(), PetscOptionsPush(), PetscOptionsPop() 1468 @*/ 1469 PetscErrorCode PetscOptionsGetAll(PetscOptions options,char *copts[]) 1470 { 1471 PetscErrorCode ierr; 1472 PetscInt i; 1473 size_t len = 1,lent = 0; 1474 char *coptions = NULL; 1475 1476 PetscFunctionBegin; 1477 PetscValidPointer(copts,2); 1478 options = options ? options : defaultoptions; 1479 /* count the length of the required string */ 1480 for (i=0; i<options->N; i++) { 1481 ierr = PetscStrlen(options->names[i],&lent);CHKERRQ(ierr); 1482 len += 2 + lent; 1483 if (options->values[i]) { 1484 ierr = PetscStrlen(options->values[i],&lent);CHKERRQ(ierr); 1485 len += 1 + lent; 1486 } 1487 } 1488 ierr = PetscMalloc1(len,&coptions);CHKERRQ(ierr); 1489 coptions[0] = 0; 1490 for (i=0; i<options->N; i++) { 1491 ierr = PetscStrcat(coptions,"-");CHKERRQ(ierr); 1492 ierr = PetscStrcat(coptions,options->names[i]);CHKERRQ(ierr); 1493 ierr = PetscStrcat(coptions," ");CHKERRQ(ierr); 1494 if (options->values[i]) { 1495 ierr = PetscStrcat(coptions,options->values[i]);CHKERRQ(ierr); 1496 ierr = PetscStrcat(coptions," ");CHKERRQ(ierr); 1497 } 1498 } 1499 *copts = coptions; 1500 PetscFunctionReturn(0); 1501 } 1502 1503 /*@C 1504 PetscOptionsUsed - Indicates if PETSc has used a particular option set in the database 1505 1506 Not Collective 1507 1508 Input Parameter: 1509 + options - options database, use NULL for default global database 1510 - name - string name of option 1511 1512 Output Parameter: 1513 . used - PETSC_TRUE if the option was used, otherwise false, including if option was not found in options database 1514 1515 Level: advanced 1516 1517 Notes: 1518 The value returned may be different on each process and depends on which options have been processed 1519 on the given process 1520 1521 .seealso: PetscOptionsView(), PetscOptionsLeft(), PetscOptionsAllUsed() 1522 @*/ 1523 PetscErrorCode PetscOptionsUsed(PetscOptions options,const char *name,PetscBool *used) 1524 { 1525 PetscInt i; 1526 PetscErrorCode ierr; 1527 1528 PetscFunctionBegin; 1529 PetscValidCharPointer(name,2); 1530 PetscValidPointer(used,3); 1531 options = options ? options : defaultoptions; 1532 *used = PETSC_FALSE; 1533 for (i=0; i<options->N; i++) { 1534 ierr = PetscStrcmp(options->names[i],name,used);CHKERRQ(ierr); 1535 if (*used) { 1536 *used = options->used[i]; 1537 break; 1538 } 1539 } 1540 PetscFunctionReturn(0); 1541 } 1542 1543 /*@ 1544 PetscOptionsAllUsed - Returns a count of the number of options in the 1545 database that have never been selected. 1546 1547 Not Collective 1548 1549 Input Parameter: 1550 . options - options database, use NULL for default global database 1551 1552 Output Parameter: 1553 . N - count of options not used 1554 1555 Level: advanced 1556 1557 Notes: 1558 The value returned may be different on each process and depends on which options have been processed 1559 on the given process 1560 1561 .seealso: PetscOptionsView() 1562 @*/ 1563 PetscErrorCode PetscOptionsAllUsed(PetscOptions options,PetscInt *N) 1564 { 1565 PetscInt i,n = 0; 1566 1567 PetscFunctionBegin; 1568 PetscValidIntPointer(N,2); 1569 options = options ? options : defaultoptions; 1570 for (i=0; i<options->N; i++) { 1571 if (!options->used[i]) n++; 1572 } 1573 *N = n; 1574 PetscFunctionReturn(0); 1575 } 1576 1577 /*@ 1578 PetscOptionsLeft - Prints to screen any options that were set and never used. 1579 1580 Not Collective 1581 1582 Input Parameter: 1583 . options - options database; use NULL for default global database 1584 1585 Options Database Key: 1586 . -options_left - activates PetscOptionsAllUsed() within PetscFinalize() 1587 1588 Notes: 1589 This is rarely used directly, it is called by PetscFinalize() in debug more or if -options_left 1590 is passed otherwise to help users determine possible mistakes in their usage of options. This 1591 only prints values on process zero of PETSC_COMM_WORLD. Other processes depending the objects 1592 used may have different options that are left unused. 1593 1594 Level: advanced 1595 1596 .seealso: PetscOptionsAllUsed() 1597 @*/ 1598 PetscErrorCode PetscOptionsLeft(PetscOptions options) 1599 { 1600 PetscErrorCode ierr; 1601 PetscInt i; 1602 PetscInt cnt = 0; 1603 PetscOptions toptions; 1604 1605 PetscFunctionBegin; 1606 toptions = options ? options : defaultoptions; 1607 for (i=0; i<toptions->N; i++) { 1608 if (!toptions->used[i]) { 1609 if (toptions->values[i]) { 1610 ierr = PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s value: %s\n",toptions->names[i],toptions->values[i]);CHKERRQ(ierr); 1611 } else { 1612 ierr = PetscPrintf(PETSC_COMM_WORLD,"Option left: name:-%s (no value)\n",toptions->names[i]);CHKERRQ(ierr); 1613 } 1614 } 1615 } 1616 if (!options) { 1617 toptions = defaultoptions; 1618 while (toptions->previous) { 1619 cnt++; 1620 toptions = toptions->previous; 1621 } 1622 if (cnt) { 1623 ierr = PetscPrintf(PETSC_COMM_WORLD,"Option left: You may have forgotten some calls to PetscOptionsPop(),\n PetscOptionsPop() has been called %D less times than PetscOptionsPush()\n",cnt);CHKERRQ(ierr); 1624 } 1625 } 1626 PetscFunctionReturn(0); 1627 } 1628 1629 /*@C 1630 PetscOptionsLeftGet - Returns all options that were set and never used. 1631 1632 Not Collective 1633 1634 Input Parameter: 1635 . options - options database, use NULL for default global database 1636 1637 Output Parameter: 1638 + N - count of options not used 1639 . names - names of options not used 1640 - values - values of options not used 1641 1642 Level: advanced 1643 1644 Notes: 1645 Users should call PetscOptionsLeftRestore() to free the memory allocated in this routine 1646 Notes: The value returned may be different on each process and depends on which options have been processed 1647 on the given process 1648 1649 .seealso: PetscOptionsAllUsed(), PetscOptionsLeft() 1650 @*/ 1651 PetscErrorCode PetscOptionsLeftGet(PetscOptions options,PetscInt *N,char **names[],char **values[]) 1652 { 1653 PetscErrorCode ierr; 1654 PetscInt i,n; 1655 1656 PetscFunctionBegin; 1657 if (N) PetscValidIntPointer(N,2); 1658 if (names) PetscValidPointer(names,3); 1659 if (values) PetscValidPointer(values,4); 1660 options = options ? options : defaultoptions; 1661 1662 /* The number of unused PETSc options */ 1663 n = 0; 1664 for (i=0; i<options->N; i++) { 1665 if (!options->used[i]) n++; 1666 } 1667 if (N) { *N = n; } 1668 if (names) { ierr = PetscMalloc1(n,names);CHKERRQ(ierr); } 1669 if (values) { ierr = PetscMalloc1(n,values);CHKERRQ(ierr); } 1670 1671 n = 0; 1672 if (names || values) { 1673 for (i=0; i<options->N; i++) { 1674 if (!options->used[i]) { 1675 if (names) (*names)[n] = options->names[i]; 1676 if (values) (*values)[n] = options->values[i]; 1677 n++; 1678 } 1679 } 1680 } 1681 PetscFunctionReturn(0); 1682 } 1683 1684 /*@C 1685 PetscOptionsLeftRestore - Free memory for the unused PETSc options obtained using PetscOptionsLeftGet. 1686 1687 Not Collective 1688 1689 Input Parameter: 1690 + options - options database, use NULL for default global database 1691 . names - names of options not used 1692 - values - values of options not used 1693 1694 Level: advanced 1695 1696 .seealso: PetscOptionsAllUsed(), PetscOptionsLeft(), PetscOptionsLeftGet() 1697 @*/ 1698 PetscErrorCode PetscOptionsLeftRestore(PetscOptions options,PetscInt *N,char **names[],char **values[]) 1699 { 1700 PetscErrorCode ierr; 1701 1702 PetscFunctionBegin; 1703 if (N) PetscValidIntPointer(N,2); 1704 if (names) PetscValidPointer(names,3); 1705 if (values) PetscValidPointer(values,4); 1706 if (N) { *N = 0; } 1707 if (names) { ierr = PetscFree(*names);CHKERRQ(ierr); } 1708 if (values) { ierr = PetscFree(*values);CHKERRQ(ierr); } 1709 PetscFunctionReturn(0); 1710 } 1711 1712 /*@C 1713 PetscOptionsSetFromOptions - Sets options related to the handling of options in PETSc 1714 1715 Collective on PETSC_COMM_WORLD 1716 1717 Input Parameter: 1718 . options - options database, use NULL for default global database 1719 1720 Options Database Keys: 1721 + -options_monitor <optional filename> - prints the names and values of all runtime options as they are set. The monitor functionality is not 1722 available for options set through a file, environment variable, or on 1723 the command line. Only options set after PetscInitialize() completes will 1724 be monitored. 1725 - -options_monitor_cancel - cancel all options database monitors 1726 1727 Notes: 1728 To see all options, run your program with the -help option 1729 1730 Level: intermediate 1731 1732 @*/ 1733 PetscErrorCode PetscOptionsSetFromOptions(PetscOptions options) 1734 { 1735 PetscBool flgc = PETSC_FALSE,flgm; 1736 PetscErrorCode ierr; 1737 char monfilename[PETSC_MAX_PATH_LEN]; 1738 PetscViewer monviewer; 1739 1740 PetscFunctionBegin; 1741 /* 1742 The options argument is currently ignored since we currently maintain only a single options database 1743 1744 options = options ? options : defaultoptions; 1745 */ 1746 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Options for handling options","PetscOptions");CHKERRQ(ierr); 1747 ierr = PetscOptionsString("-options_monitor","Monitor options database","PetscOptionsMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flgm);CHKERRQ(ierr); 1748 ierr = PetscOptionsBool("-options_monitor_cancel","Cancel all options database monitors","PetscOptionsMonitorCancel",flgc,&flgc,NULL);CHKERRQ(ierr); 1749 ierr = PetscOptionsEnd();CHKERRQ(ierr); 1750 if (flgm) { 1751 ierr = PetscViewerASCIIOpen(PETSC_COMM_WORLD,monfilename,&monviewer);CHKERRQ(ierr); 1752 ierr = PetscOptionsMonitorSet(PetscOptionsMonitorDefault,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);CHKERRQ(ierr); 1753 } 1754 if (flgc) { ierr = PetscOptionsMonitorCancel();CHKERRQ(ierr); } 1755 PetscFunctionReturn(0); 1756 } 1757 1758 /*@C 1759 PetscOptionsMonitorDefault - Print all options set value events. 1760 1761 Logically Collective on ctx 1762 1763 Input Parameters: 1764 + name - option name string 1765 . value - option value string 1766 - ctx - an ASCII viewer 1767 1768 Level: intermediate 1769 1770 Notes: 1771 The first MPI rank in the PetscViewer viewer actually prints the values, other 1772 processes may have different values set 1773 1774 .seealso: PetscOptionsMonitorSet() 1775 @*/ 1776 PetscErrorCode PetscOptionsMonitorDefault(const char name[],const char value[],void *ctx) 1777 { 1778 PetscErrorCode ierr; 1779 PetscViewer viewer = (PetscViewer)ctx; 1780 1781 PetscFunctionBegin; 1782 if (!value) { 1783 ierr = PetscViewerASCIIPrintf(viewer,"Removing option: %s\n",name,value);CHKERRQ(ierr); 1784 } else if (!value[0]) { 1785 ierr = PetscViewerASCIIPrintf(viewer,"Setting option: %s (no value)\n",name);CHKERRQ(ierr); 1786 } else { 1787 ierr = PetscViewerASCIIPrintf(viewer,"Setting option: %s = %s\n",name,value);CHKERRQ(ierr); 1788 } 1789 PetscFunctionReturn(0); 1790 } 1791 1792 /*@C 1793 PetscOptionsMonitorSet - Sets an ADDITIONAL function to be called at every method that 1794 modified the PETSc options database. 1795 1796 Not Collective 1797 1798 Input Parameters: 1799 + monitor - pointer to function (if this is NULL, it turns off monitoring 1800 . mctx - [optional] context for private data for the 1801 monitor routine (use NULL if no context is desired) 1802 - monitordestroy - [optional] routine that frees monitor context 1803 (may be NULL) 1804 1805 Calling Sequence of monitor: 1806 $ monitor (const char name[], const char value[], void *mctx) 1807 1808 + name - option name string 1809 . value - option value string 1810 - mctx - optional monitoring context, as set by PetscOptionsMonitorSet() 1811 1812 Options Database Keys: 1813 + -options_monitor - sets PetscOptionsMonitorDefault() 1814 - -options_monitor_cancel - cancels all monitors that have 1815 been hardwired into a code by 1816 calls to PetscOptionsMonitorSet(), but 1817 does not cancel those set via 1818 the options database. 1819 1820 Notes: 1821 The default is to do nothing. To print the name and value of options 1822 being inserted into the database, use PetscOptionsMonitorDefault() as the monitoring routine, 1823 with a null monitoring context. 1824 1825 Several different monitoring routines may be set by calling 1826 PetscOptionsMonitorSet() multiple times; all will be called in the 1827 order in which they were set. 1828 1829 Level: beginner 1830 1831 .seealso: PetscOptionsMonitorDefault(), PetscOptionsMonitorCancel() 1832 @*/ 1833 PetscErrorCode PetscOptionsMonitorSet(PetscErrorCode (*monitor)(const char name[], const char value[], void*),void *mctx,PetscErrorCode (*monitordestroy)(void**)) 1834 { 1835 PetscOptions options = defaultoptions; 1836 1837 PetscFunctionBegin; 1838 if (options->numbermonitors >= MAXOPTIONSMONITORS) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many PetscOptions monitors set"); 1839 options->monitor[options->numbermonitors] = monitor; 1840 options->monitordestroy[options->numbermonitors] = monitordestroy; 1841 options->monitorcontext[options->numbermonitors++] = (void*)mctx; 1842 PetscFunctionReturn(0); 1843 } 1844 1845 /*@ 1846 PetscOptionsMonitorCancel - Clears all monitors for a PetscOptions object. 1847 1848 Not Collective 1849 1850 Options Database Key: 1851 . -options_monitor_cancel - Cancels all monitors that have 1852 been hardwired into a code by calls to PetscOptionsMonitorSet(), 1853 but does not cancel those set via the options database. 1854 1855 Level: intermediate 1856 1857 .seealso: PetscOptionsMonitorDefault(), PetscOptionsMonitorSet() 1858 @*/ 1859 PetscErrorCode PetscOptionsMonitorCancel(void) 1860 { 1861 PetscErrorCode ierr; 1862 PetscInt i; 1863 PetscOptions options = defaultoptions; 1864 1865 PetscFunctionBegin; 1866 for (i=0; i<options->numbermonitors; i++) { 1867 if (options->monitordestroy[i]) { 1868 ierr = (*options->monitordestroy[i])(&options->monitorcontext[i]);CHKERRQ(ierr); 1869 } 1870 } 1871 options->numbermonitors = 0; 1872 PetscFunctionReturn(0); 1873 } 1874 1875 /* 1876 PetscOptionsStringToBool - Converts string to PetscBool , handles cases like "yes", "no", "true", "false", "0", "1", "off", "on". 1877 */ 1878 PetscErrorCode PetscOptionsStringToBool(const char value[],PetscBool *a) 1879 { 1880 PetscBool istrue,isfalse; 1881 size_t len; 1882 PetscErrorCode ierr; 1883 1884 PetscFunctionBegin; 1885 ierr = PetscStrlen(value,&len);CHKERRQ(ierr); 1886 if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Character string of length zero has no logical value"); 1887 ierr = PetscStrcasecmp(value,"TRUE",&istrue);CHKERRQ(ierr); 1888 if (istrue) {*a = PETSC_TRUE; PetscFunctionReturn(0);} 1889 ierr = PetscStrcasecmp(value,"YES",&istrue);CHKERRQ(ierr); 1890 if (istrue) {*a = PETSC_TRUE; PetscFunctionReturn(0);} 1891 ierr = PetscStrcasecmp(value,"1",&istrue);CHKERRQ(ierr); 1892 if (istrue) {*a = PETSC_TRUE; PetscFunctionReturn(0);} 1893 ierr = PetscStrcasecmp(value,"on",&istrue);CHKERRQ(ierr); 1894 if (istrue) {*a = PETSC_TRUE; PetscFunctionReturn(0);} 1895 ierr = PetscStrcasecmp(value,"FALSE",&isfalse);CHKERRQ(ierr); 1896 if (isfalse) {*a = PETSC_FALSE; PetscFunctionReturn(0);} 1897 ierr = PetscStrcasecmp(value,"NO",&isfalse);CHKERRQ(ierr); 1898 if (isfalse) {*a = PETSC_FALSE; PetscFunctionReturn(0);} 1899 ierr = PetscStrcasecmp(value,"0",&isfalse);CHKERRQ(ierr); 1900 if (isfalse) {*a = PETSC_FALSE; PetscFunctionReturn(0);} 1901 ierr = PetscStrcasecmp(value,"off",&isfalse);CHKERRQ(ierr); 1902 if (isfalse) {*a = PETSC_FALSE; PetscFunctionReturn(0);} 1903 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unknown logical value: %s",value); 1904 } 1905 1906 /* 1907 PetscOptionsStringToInt - Converts a string to an integer value. Handles special cases such as "default" and "decide" 1908 */ 1909 PetscErrorCode PetscOptionsStringToInt(const char name[],PetscInt *a) 1910 { 1911 PetscErrorCode ierr; 1912 size_t len; 1913 PetscBool decide,tdefault,mouse; 1914 1915 PetscFunctionBegin; 1916 ierr = PetscStrlen(name,&len);CHKERRQ(ierr); 1917 if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"character string of length zero has no numerical value"); 1918 1919 ierr = PetscStrcasecmp(name,"PETSC_DEFAULT",&tdefault);CHKERRQ(ierr); 1920 if (!tdefault) { 1921 ierr = PetscStrcasecmp(name,"DEFAULT",&tdefault);CHKERRQ(ierr); 1922 } 1923 ierr = PetscStrcasecmp(name,"PETSC_DECIDE",&decide);CHKERRQ(ierr); 1924 if (!decide) { 1925 ierr = PetscStrcasecmp(name,"DECIDE",&decide);CHKERRQ(ierr); 1926 } 1927 ierr = PetscStrcasecmp(name,"mouse",&mouse);CHKERRQ(ierr); 1928 1929 if (tdefault) *a = PETSC_DEFAULT; 1930 else if (decide) *a = PETSC_DECIDE; 1931 else if (mouse) *a = -1; 1932 else { 1933 char *endptr; 1934 long strtolval; 1935 1936 strtolval = strtol(name,&endptr,10); 1937 if ((size_t) (endptr - name) != len) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s has no integer value (do not include . in it)",name); 1938 1939 #if defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE_ATOLL) 1940 (void) strtolval; 1941 *a = atoll(name); 1942 #elif defined(PETSC_USE_64BIT_INDICES) && defined(PETSC_HAVE___INT64) 1943 (void) strtolval; 1944 *a = _atoi64(name); 1945 #else 1946 *a = (PetscInt)strtolval; 1947 #endif 1948 } 1949 PetscFunctionReturn(0); 1950 } 1951 1952 #if defined(PETSC_USE_REAL___FLOAT128) 1953 #include <quadmath.h> 1954 #endif 1955 1956 static PetscErrorCode PetscStrtod(const char name[],PetscReal *a,char **endptr) 1957 { 1958 PetscFunctionBegin; 1959 #if defined(PETSC_USE_REAL___FLOAT128) 1960 *a = strtoflt128(name,endptr); 1961 #else 1962 *a = (PetscReal)strtod(name,endptr); 1963 #endif 1964 PetscFunctionReturn(0); 1965 } 1966 1967 static PetscErrorCode PetscStrtoz(const char name[],PetscScalar *a,char **endptr,PetscBool *isImaginary) 1968 { 1969 PetscBool hasi = PETSC_FALSE; 1970 char *ptr; 1971 PetscReal strtoval; 1972 PetscErrorCode ierr; 1973 1974 PetscFunctionBegin; 1975 ierr = PetscStrtod(name,&strtoval,&ptr);CHKERRQ(ierr); 1976 if (ptr == name) { 1977 strtoval = 1.; 1978 hasi = PETSC_TRUE; 1979 if (name[0] == 'i') { 1980 ptr++; 1981 } else if (name[0] == '+' && name[1] == 'i') { 1982 ptr += 2; 1983 } else if (name[0] == '-' && name[1] == 'i') { 1984 strtoval = -1.; 1985 ptr += 2; 1986 } 1987 } else if (*ptr == 'i') { 1988 hasi = PETSC_TRUE; 1989 ptr++; 1990 } 1991 *endptr = ptr; 1992 *isImaginary = hasi; 1993 if (hasi) { 1994 #if !defined(PETSC_USE_COMPLEX) 1995 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s contains imaginary but complex not supported ",name); 1996 #else 1997 *a = PetscCMPLX(0.,strtoval); 1998 #endif 1999 } else { 2000 *a = strtoval; 2001 } 2002 PetscFunctionReturn(0); 2003 } 2004 2005 /* 2006 Converts a string to PetscReal value. Handles special cases like "default" and "decide" 2007 */ 2008 PetscErrorCode PetscOptionsStringToReal(const char name[],PetscReal *a) 2009 { 2010 size_t len; 2011 PetscBool match; 2012 char *endptr; 2013 PetscErrorCode ierr; 2014 2015 PetscFunctionBegin; 2016 ierr = PetscStrlen(name,&len);CHKERRQ(ierr); 2017 if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"String of length zero has no numerical value"); 2018 2019 ierr = PetscStrcasecmp(name,"PETSC_DEFAULT",&match);CHKERRQ(ierr); 2020 if (!match) { 2021 ierr = PetscStrcasecmp(name,"DEFAULT",&match);CHKERRQ(ierr); 2022 } 2023 if (match) {*a = PETSC_DEFAULT; PetscFunctionReturn(0);} 2024 2025 ierr = PetscStrcasecmp(name,"PETSC_DECIDE",&match);CHKERRQ(ierr); 2026 if (!match) { 2027 ierr = PetscStrcasecmp(name,"DECIDE",&match);CHKERRQ(ierr); 2028 } 2029 if (match) {*a = PETSC_DECIDE; PetscFunctionReturn(0);} 2030 2031 ierr = PetscStrtod(name,a,&endptr);CHKERRQ(ierr); 2032 if ((size_t) (endptr - name) != len) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s has no numeric value",name); 2033 PetscFunctionReturn(0); 2034 } 2035 2036 PetscErrorCode PetscOptionsStringToScalar(const char name[],PetscScalar *a) 2037 { 2038 PetscBool imag1; 2039 size_t len; 2040 PetscScalar val = 0.; 2041 char *ptr = NULL; 2042 PetscErrorCode ierr; 2043 2044 PetscFunctionBegin; 2045 ierr = PetscStrlen(name,&len);CHKERRQ(ierr); 2046 if (!len) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"character string of length zero has no numerical value"); 2047 ierr = PetscStrtoz(name,&val,&ptr,&imag1);CHKERRQ(ierr); 2048 #if defined(PETSC_USE_COMPLEX) 2049 if ((size_t) (ptr - name) < len) { 2050 PetscBool imag2; 2051 PetscScalar val2; 2052 2053 ierr = PetscStrtoz(ptr,&val2,&ptr,&imag2);CHKERRQ(ierr); 2054 if (imag1 || !imag2) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s: must specify imaginary component second",name); 2055 val = PetscCMPLX(PetscRealPart(val),PetscImaginaryPart(val2)); 2056 } 2057 #endif 2058 if ((size_t) (ptr - name) != len) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Input string %s has no numeric value ",name); 2059 *a = val; 2060 PetscFunctionReturn(0); 2061 } 2062 2063 /*@C 2064 PetscOptionsGetBool - Gets the Logical (true or false) value for a particular 2065 option in the database. 2066 2067 Not Collective 2068 2069 Input Parameters: 2070 + options - options database, use NULL for default global database 2071 . pre - the string to prepend to the name or NULL 2072 - name - the option one is seeking 2073 2074 Output Parameter: 2075 + ivalue - the logical value to return 2076 - set - PETSC_TRUE if found, else PETSC_FALSE 2077 2078 Level: beginner 2079 2080 Notes: 2081 TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE 2082 FALSE, false, NO, no, and 0 all translate to PETSC_FALSE 2083 2084 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 2085 is equivalent to -requested_bool true 2086 2087 If the user does not supply the option at all ivalue is NOT changed. Thus 2088 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2089 2090 .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), 2091 PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsGetInt(), PetscOptionsBool(), 2092 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2093 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2094 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2095 PetscOptionsFList(), PetscOptionsEList() 2096 @*/ 2097 PetscErrorCode PetscOptionsGetBool(PetscOptions options,const char pre[],const char name[],PetscBool *ivalue,PetscBool *set) 2098 { 2099 const char *value; 2100 PetscBool flag; 2101 PetscErrorCode ierr; 2102 2103 PetscFunctionBegin; 2104 PetscValidCharPointer(name,3); 2105 if (ivalue) PetscValidIntPointer(ivalue,4); 2106 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);CHKERRQ(ierr); 2107 if (flag) { 2108 if (set) *set = PETSC_TRUE; 2109 if (!value) { 2110 if (ivalue) *ivalue = PETSC_TRUE; 2111 } else { 2112 ierr = PetscOptionsStringToBool(value, &flag);CHKERRQ(ierr); 2113 if (ivalue) *ivalue = flag; 2114 } 2115 } else { 2116 if (set) *set = PETSC_FALSE; 2117 } 2118 PetscFunctionReturn(0); 2119 } 2120 2121 /*@C 2122 PetscOptionsGetEList - Puts a list of option values that a single one may be selected from 2123 2124 Not Collective 2125 2126 Input Parameters: 2127 + options - options database, use NULL for default global database 2128 . pre - the string to prepend to the name or NULL 2129 . opt - option name 2130 . list - the possible choices (one of these must be selected, anything else is invalid) 2131 - ntext - number of choices 2132 2133 Output Parameter: 2134 + value - the index of the value to return (defaults to zero if the option name is given but no choice is listed) 2135 - set - PETSC_TRUE if found, else PETSC_FALSE 2136 2137 Level: intermediate 2138 2139 Notes: 2140 If the user does not supply the option value is NOT changed. Thus 2141 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2142 2143 See PetscOptionsFList() for when the choices are given in a PetscFunctionList() 2144 2145 .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(), 2146 PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(), 2147 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2148 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2149 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2150 PetscOptionsFList(), PetscOptionsEList() 2151 @*/ 2152 PetscErrorCode PetscOptionsGetEList(PetscOptions options,const char pre[],const char opt[],const char * const *list,PetscInt ntext,PetscInt *value,PetscBool *set) 2153 { 2154 PetscErrorCode ierr; 2155 size_t alen,len = 0, tlen = 0; 2156 char *svalue; 2157 PetscBool aset,flg = PETSC_FALSE; 2158 PetscInt i; 2159 2160 PetscFunctionBegin; 2161 PetscValidCharPointer(opt,3); 2162 for (i=0; i<ntext; i++) { 2163 ierr = PetscStrlen(list[i],&alen);CHKERRQ(ierr); 2164 if (alen > len) len = alen; 2165 tlen += len + 1; 2166 } 2167 len += 5; /* a little extra space for user mistypes */ 2168 ierr = PetscMalloc1(len,&svalue);CHKERRQ(ierr); 2169 ierr = PetscOptionsGetString(options,pre,opt,svalue,len,&aset);CHKERRQ(ierr); 2170 if (aset) { 2171 ierr = PetscEListFind(ntext,list,svalue,value,&flg);CHKERRQ(ierr); 2172 if (!flg) { 2173 char *avail,*pavl; 2174 2175 ierr = PetscMalloc1(tlen,&avail);CHKERRQ(ierr); 2176 pavl = avail; 2177 for (i=0; i<ntext; i++) { 2178 ierr = PetscStrlen(list[i],&alen);CHKERRQ(ierr); 2179 ierr = PetscStrcpy(pavl,list[i]);CHKERRQ(ierr); 2180 pavl += alen; 2181 ierr = PetscStrcpy(pavl," ");CHKERRQ(ierr); 2182 pavl += 1; 2183 } 2184 ierr = PetscStrtolower(avail);CHKERRQ(ierr); 2185 SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_USER,"Unknown option %s for -%s%s. Available options: %s",svalue,pre ? pre : "",opt+1,avail); 2186 } 2187 if (set) *set = PETSC_TRUE; 2188 } else if (set) *set = PETSC_FALSE; 2189 ierr = PetscFree(svalue);CHKERRQ(ierr); 2190 PetscFunctionReturn(0); 2191 } 2192 2193 /*@C 2194 PetscOptionsGetEnum - Gets the enum value for a particular option in the database. 2195 2196 Not Collective 2197 2198 Input Parameters: 2199 + options - options database, use NULL for default global database 2200 . pre - option prefix or NULL 2201 . opt - option name 2202 . list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null 2203 - defaultv - the default (current) value 2204 2205 Output Parameter: 2206 + value - the value to return 2207 - set - PETSC_TRUE if found, else PETSC_FALSE 2208 2209 Level: beginner 2210 2211 Notes: 2212 If the user does not supply the option value is NOT changed. Thus 2213 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2214 2215 List is usually something like PCASMTypes or some other predefined list of enum names 2216 2217 .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), PetscOptionsGetInt(), 2218 PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool() 2219 PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(), 2220 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2221 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2222 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2223 PetscOptionsFList(), PetscOptionsEList(), PetscOptionsGetEList(), PetscOptionsEnum() 2224 @*/ 2225 PetscErrorCode PetscOptionsGetEnum(PetscOptions options,const char pre[],const char opt[],const char * const *list,PetscEnum *value,PetscBool *set) 2226 { 2227 PetscErrorCode ierr; 2228 PetscInt ntext = 0,tval; 2229 PetscBool fset; 2230 2231 PetscFunctionBegin; 2232 PetscValidCharPointer(opt,3); 2233 while (list[ntext++]) { 2234 if (ntext > 50) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"List argument appears to be wrong or have more than 50 entries"); 2235 } 2236 if (ntext < 3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"List argument must have at least two entries: typename and type prefix"); 2237 ntext -= 3; 2238 ierr = PetscOptionsGetEList(options,pre,opt,list,ntext,&tval,&fset);CHKERRQ(ierr); 2239 /* with PETSC_USE_64BIT_INDICES sizeof(PetscInt) != sizeof(PetscEnum) */ 2240 if (fset) *value = (PetscEnum)tval; 2241 if (set) *set = fset; 2242 PetscFunctionReturn(0); 2243 } 2244 2245 /*@C 2246 PetscOptionsGetInt - Gets the integer value for a particular option in the database. 2247 2248 Not Collective 2249 2250 Input Parameters: 2251 + options - options database, use NULL for default global database 2252 . pre - the string to prepend to the name or NULL 2253 - name - the option one is seeking 2254 2255 Output Parameter: 2256 + ivalue - the integer value to return 2257 - set - PETSC_TRUE if found, else PETSC_FALSE 2258 2259 Level: beginner 2260 2261 Notes: 2262 If the user does not supply the option ivalue is NOT changed. Thus 2263 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2264 2265 .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), 2266 PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool() 2267 PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(), 2268 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2269 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2270 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2271 PetscOptionsFList(), PetscOptionsEList() 2272 @*/ 2273 PetscErrorCode PetscOptionsGetInt(PetscOptions options,const char pre[],const char name[],PetscInt *ivalue,PetscBool *set) 2274 { 2275 const char *value; 2276 PetscErrorCode ierr; 2277 PetscBool flag; 2278 2279 PetscFunctionBegin; 2280 PetscValidCharPointer(name,3); 2281 PetscValidIntPointer(ivalue,4); 2282 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);CHKERRQ(ierr); 2283 if (flag) { 2284 if (!value) { 2285 if (set) *set = PETSC_FALSE; 2286 } else { 2287 if (set) *set = PETSC_TRUE; 2288 ierr = PetscOptionsStringToInt(value,ivalue);CHKERRQ(ierr); 2289 } 2290 } else { 2291 if (set) *set = PETSC_FALSE; 2292 } 2293 PetscFunctionReturn(0); 2294 } 2295 2296 /*@C 2297 PetscOptionsGetReal - Gets the double precision value for a particular 2298 option in the database. 2299 2300 Not Collective 2301 2302 Input Parameters: 2303 + options - options database, use NULL for default global database 2304 . pre - string to prepend to each name or NULL 2305 - name - the option one is seeking 2306 2307 Output Parameter: 2308 + dvalue - the double value to return 2309 - set - PETSC_TRUE if found, PETSC_FALSE if not found 2310 2311 Notes: 2312 If the user does not supply the option dvalue is NOT changed. Thus 2313 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2314 2315 Level: beginner 2316 2317 .seealso: PetscOptionsGetInt(), PetscOptionsHasName(), 2318 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(),PetscOptionsBool(), 2319 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2320 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2321 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2322 PetscOptionsFList(), PetscOptionsEList() 2323 @*/ 2324 PetscErrorCode PetscOptionsGetReal(PetscOptions options,const char pre[],const char name[],PetscReal *dvalue,PetscBool *set) 2325 { 2326 const char *value; 2327 PetscBool flag; 2328 PetscErrorCode ierr; 2329 2330 PetscFunctionBegin; 2331 PetscValidCharPointer(name,3); 2332 PetscValidRealPointer(dvalue,4); 2333 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);CHKERRQ(ierr); 2334 if (flag) { 2335 if (!value) { 2336 if (set) *set = PETSC_FALSE; 2337 } else { 2338 if (set) *set = PETSC_TRUE; 2339 ierr = PetscOptionsStringToReal(value,dvalue);CHKERRQ(ierr); 2340 } 2341 } else { 2342 if (set) *set = PETSC_FALSE; 2343 } 2344 PetscFunctionReturn(0); 2345 } 2346 2347 /*@C 2348 PetscOptionsGetScalar - Gets the scalar value for a particular 2349 option in the database. 2350 2351 Not Collective 2352 2353 Input Parameters: 2354 + options - options database, use NULL for default global database 2355 . pre - string to prepend to each name or NULL 2356 - name - the option one is seeking 2357 2358 Output Parameter: 2359 + dvalue - the double value to return 2360 - set - PETSC_TRUE if found, else PETSC_FALSE 2361 2362 Level: beginner 2363 2364 Usage: 2365 A complex number 2+3i must be specified with NO spaces 2366 2367 Notes: 2368 If the user does not supply the option dvalue is NOT changed. Thus 2369 you should ALWAYS initialize the ivalue if you access it without first checking if the set flag is true. 2370 2371 .seealso: PetscOptionsGetInt(), PetscOptionsHasName(), 2372 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(), 2373 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2374 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2375 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2376 PetscOptionsFList(), PetscOptionsEList() 2377 @*/ 2378 PetscErrorCode PetscOptionsGetScalar(PetscOptions options,const char pre[],const char name[],PetscScalar *dvalue,PetscBool *set) 2379 { 2380 const char *value; 2381 PetscBool flag; 2382 PetscErrorCode ierr; 2383 2384 PetscFunctionBegin; 2385 PetscValidCharPointer(name,3); 2386 PetscValidScalarPointer(dvalue,4); 2387 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);CHKERRQ(ierr); 2388 if (flag) { 2389 if (!value) { 2390 if (set) *set = PETSC_FALSE; 2391 } else { 2392 #if !defined(PETSC_USE_COMPLEX) 2393 ierr = PetscOptionsStringToReal(value,dvalue);CHKERRQ(ierr); 2394 #else 2395 ierr = PetscOptionsStringToScalar(value,dvalue);CHKERRQ(ierr); 2396 #endif 2397 if (set) *set = PETSC_TRUE; 2398 } 2399 } else { /* flag */ 2400 if (set) *set = PETSC_FALSE; 2401 } 2402 PetscFunctionReturn(0); 2403 } 2404 2405 /*@C 2406 PetscOptionsGetString - Gets the string value for a particular option in 2407 the database. 2408 2409 Not Collective 2410 2411 Input Parameters: 2412 + options - options database, use NULL for default global database 2413 . pre - string to prepend to name or NULL 2414 . name - the option one is seeking 2415 - len - maximum length of the string including null termination 2416 2417 Output Parameters: 2418 + string - location to copy string 2419 - set - PETSC_TRUE if found, else PETSC_FALSE 2420 2421 Level: beginner 2422 2423 Fortran Note: 2424 The Fortran interface is slightly different from the C/C++ 2425 interface (len is not used). Sample usage in Fortran follows 2426 .vb 2427 character *20 string 2428 PetscErrorCode ierr 2429 PetscBool set 2430 call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-s',string,set,ierr) 2431 .ve 2432 2433 Notes: 2434 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 2435 2436 If the user does not use the option then the string is not changed. Thus 2437 you should ALWAYS initialize the string if you access it without first checking if the set flag is true. 2438 2439 Note: 2440 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). 2441 2442 .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(), 2443 PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(), 2444 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2445 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2446 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2447 PetscOptionsFList(), PetscOptionsEList() 2448 @*/ 2449 PetscErrorCode PetscOptionsGetString(PetscOptions options,const char pre[],const char name[],char string[],size_t len,PetscBool *set) 2450 { 2451 const char *value; 2452 PetscBool flag; 2453 PetscErrorCode ierr; 2454 2455 PetscFunctionBegin; 2456 PetscValidCharPointer(name,3); 2457 PetscValidCharPointer(string,4); 2458 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);CHKERRQ(ierr); 2459 if (!flag) { 2460 if (set) *set = PETSC_FALSE; 2461 } else { 2462 if (set) *set = PETSC_TRUE; 2463 if (value) { 2464 ierr = PetscStrncpy(string,value,len);CHKERRQ(ierr); 2465 } else { 2466 ierr = PetscArrayzero(string,len);CHKERRQ(ierr); 2467 } 2468 } 2469 PetscFunctionReturn(0); 2470 } 2471 2472 char *PetscOptionsGetStringMatlab(PetscOptions options,const char pre[],const char name[]) 2473 { 2474 const char *value; 2475 PetscBool flag; 2476 PetscErrorCode ierr; 2477 2478 PetscFunctionBegin; 2479 ierr = PetscOptionsFindPair(options,pre,name,&value,&flag);if (ierr) PetscFunctionReturn(NULL); 2480 if (flag) PetscFunctionReturn((char*)value); 2481 else PetscFunctionReturn(NULL); 2482 } 2483 2484 /*@C 2485 PetscOptionsGetBoolArray - Gets an array of Logical (true or false) values for a particular 2486 option in the database. The values must be separated with commas with 2487 no intervening spaces. 2488 2489 Not Collective 2490 2491 Input Parameters: 2492 + options - options database, use NULL for default global database 2493 . pre - string to prepend to each name or NULL 2494 . name - the option one is seeking 2495 - nmax - maximum number of values to retrieve 2496 2497 Output Parameter: 2498 + dvalue - the integer values to return 2499 . nmax - actual number of values retreived 2500 - set - PETSC_TRUE if found, else PETSC_FALSE 2501 2502 Level: beginner 2503 2504 Notes: 2505 TRUE, true, YES, yes, nostring, and 1 all translate to PETSC_TRUE 2506 FALSE, false, NO, no, and 0 all translate to PETSC_FALSE 2507 2508 .seealso: PetscOptionsGetInt(), PetscOptionsHasName(), 2509 PetscOptionsGetString(), PetscOptionsGetRealArray(), PetscOptionsBool(), 2510 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2511 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2512 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2513 PetscOptionsFList(), PetscOptionsEList() 2514 @*/ 2515 PetscErrorCode PetscOptionsGetBoolArray(PetscOptions options,const char pre[],const char name[],PetscBool dvalue[],PetscInt *nmax,PetscBool *set) 2516 { 2517 const char *svalue; 2518 char *value; 2519 PetscErrorCode ierr; 2520 PetscInt n = 0; 2521 PetscBool flag; 2522 PetscToken token; 2523 2524 PetscFunctionBegin; 2525 PetscValidCharPointer(name,3); 2526 PetscValidIntPointer(dvalue,4); 2527 PetscValidIntPointer(nmax,5); 2528 2529 ierr = PetscOptionsFindPair(options,pre,name,&svalue,&flag);CHKERRQ(ierr); 2530 if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; PetscFunctionReturn(0);} 2531 if (set) *set = PETSC_TRUE; 2532 ierr = PetscTokenCreate(svalue,',',&token);CHKERRQ(ierr); 2533 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2534 while (value && n < *nmax) { 2535 ierr = PetscOptionsStringToBool(value,dvalue);CHKERRQ(ierr); 2536 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2537 dvalue++; 2538 n++; 2539 } 2540 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 2541 *nmax = n; 2542 PetscFunctionReturn(0); 2543 } 2544 2545 /*@C 2546 PetscOptionsGetEnumArray - Gets an array of enum values for a particular option in the database. 2547 2548 Not Collective 2549 2550 Input Parameters: 2551 + options - options database, use NULL for default global database 2552 . pre - option prefix or NULL 2553 . name - option name 2554 . list - array containing the list of choices, followed by the enum name, followed by the enum prefix, followed by a null 2555 - nmax - maximum number of values to retrieve 2556 2557 Output Parameters: 2558 + ivalue - the enum values to return 2559 . nmax - actual number of values retreived 2560 - set - PETSC_TRUE if found, else PETSC_FALSE 2561 2562 Level: beginner 2563 2564 Notes: 2565 The array must be passed as a comma separated list. 2566 2567 There must be no intervening spaces between the values. 2568 2569 list is usually something like PCASMTypes or some other predefined list of enum names. 2570 2571 .seealso: PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(), PetscOptionsGetInt(), 2572 PetscOptionsGetEnum(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool() 2573 PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(), PetscOptionsName(), 2574 PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), PetscOptionsStringArray(),PetscOptionsRealArray(), 2575 PetscOptionsScalar(), PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2576 PetscOptionsFList(), PetscOptionsEList(), PetscOptionsGetEList(), PetscOptionsEnum() 2577 @*/ 2578 PetscErrorCode PetscOptionsGetEnumArray(PetscOptions options,const char pre[],const char name[],const char *const *list,PetscEnum ivalue[],PetscInt *nmax,PetscBool *set) 2579 { 2580 const char *svalue; 2581 char *value; 2582 PetscInt n = 0; 2583 PetscEnum evalue; 2584 PetscBool flag; 2585 PetscToken token; 2586 PetscErrorCode ierr; 2587 2588 PetscFunctionBegin; 2589 PetscValidCharPointer(name,3); 2590 PetscValidPointer(list,4); 2591 PetscValidPointer(ivalue,5); 2592 PetscValidIntPointer(nmax,6); 2593 2594 ierr = PetscOptionsFindPair(options,pre,name,&svalue,&flag);CHKERRQ(ierr); 2595 if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; PetscFunctionReturn(0);} 2596 if (set) *set = PETSC_TRUE; 2597 ierr = PetscTokenCreate(svalue,',',&token);CHKERRQ(ierr); 2598 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2599 while (value && n < *nmax) { 2600 ierr = PetscEnumFind(list,value,&evalue,&flag);CHKERRQ(ierr); 2601 if (!flag) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_USER,"Unknown enum value '%s' for -%s%s",svalue,pre ? pre : "",name+1); 2602 ivalue[n++] = evalue; 2603 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2604 } 2605 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 2606 *nmax = n; 2607 PetscFunctionReturn(0); 2608 } 2609 2610 /*@C 2611 PetscOptionsGetIntArray - Gets an array of integer values for a particular 2612 option in the database. 2613 2614 Not Collective 2615 2616 Input Parameters: 2617 + options - options database, use NULL for default global database 2618 . pre - string to prepend to each name or NULL 2619 . name - the option one is seeking 2620 - nmax - maximum number of values to retrieve 2621 2622 Output Parameter: 2623 + ivalue - the integer values to return 2624 . nmax - actual number of values retreived 2625 - set - PETSC_TRUE if found, else PETSC_FALSE 2626 2627 Level: beginner 2628 2629 Notes: 2630 The array can be passed as 2631 a comma separated list: 0,1,2,3,4,5,6,7 2632 a range (start-end+1): 0-8 2633 a range with given increment (start-end+1:inc): 0-7:2 2634 a combination of values and ranges separated by commas: 0,1-8,8-15:2 2635 2636 There must be no intervening spaces between the values. 2637 2638 .seealso: PetscOptionsGetInt(), PetscOptionsHasName(), 2639 PetscOptionsGetString(), PetscOptionsGetRealArray(), PetscOptionsBool(), 2640 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2641 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2642 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2643 PetscOptionsFList(), PetscOptionsEList() 2644 @*/ 2645 PetscErrorCode PetscOptionsGetIntArray(PetscOptions options,const char pre[],const char name[],PetscInt ivalue[],PetscInt *nmax,PetscBool *set) 2646 { 2647 const char *svalue; 2648 char *value; 2649 PetscErrorCode ierr; 2650 PetscInt n = 0,i,j,start,end,inc,nvalues; 2651 size_t len; 2652 PetscBool flag,foundrange; 2653 PetscToken token; 2654 2655 PetscFunctionBegin; 2656 PetscValidCharPointer(name,3); 2657 PetscValidIntPointer(ivalue,4); 2658 PetscValidIntPointer(nmax,5); 2659 2660 ierr = PetscOptionsFindPair(options,pre,name,&svalue,&flag);CHKERRQ(ierr); 2661 if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; PetscFunctionReturn(0);} 2662 if (set) *set = PETSC_TRUE; 2663 ierr = PetscTokenCreate(svalue,',',&token);CHKERRQ(ierr); 2664 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2665 while (value && n < *nmax) { 2666 /* look for form d-D where d and D are integers */ 2667 foundrange = PETSC_FALSE; 2668 ierr = PetscStrlen(value,&len);CHKERRQ(ierr); 2669 if (value[0] == '-') i=2; 2670 else i=1; 2671 for (;i<(int)len; i++) { 2672 if (value[i] == '-') { 2673 if (i == (int)len-1) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry %s\n",n,value); 2674 value[i] = 0; 2675 2676 ierr = PetscOptionsStringToInt(value,&start);CHKERRQ(ierr); 2677 inc = 1; 2678 j = i+1; 2679 for (;j<(int)len; j++) { 2680 if (value[j] == ':') { 2681 value[j] = 0; 2682 2683 ierr = PetscOptionsStringToInt(value+j+1,&inc);CHKERRQ(ierr); 2684 if (inc <= 0) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry,%s cannot have negative increment",n,value+j+1); 2685 break; 2686 } 2687 } 2688 ierr = PetscOptionsStringToInt(value+i+1,&end);CHKERRQ(ierr); 2689 if (end <= start) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry, %s-%s cannot have decreasing list",n,value,value+i+1); 2690 nvalues = (end-start)/inc + (end-start)%inc; 2691 if (n + nvalues > *nmax) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_USER,"Error in %D-th array entry, not enough space left in array (%D) to contain entire range from %D to %D",n,*nmax-n,start,end); 2692 for (;start<end; start+=inc) { 2693 *ivalue = start; ivalue++;n++; 2694 } 2695 foundrange = PETSC_TRUE; 2696 break; 2697 } 2698 } 2699 if (!foundrange) { 2700 ierr = PetscOptionsStringToInt(value,ivalue);CHKERRQ(ierr); 2701 ivalue++; 2702 n++; 2703 } 2704 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2705 } 2706 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 2707 *nmax = n; 2708 PetscFunctionReturn(0); 2709 } 2710 2711 /*@C 2712 PetscOptionsGetRealArray - Gets an array of double precision values for a 2713 particular option in the database. The values must be separated with 2714 commas with no intervening spaces. 2715 2716 Not Collective 2717 2718 Input Parameters: 2719 + options - options database, use NULL for default global database 2720 . pre - string to prepend to each name or NULL 2721 . name - the option one is seeking 2722 - nmax - maximum number of values to retrieve 2723 2724 Output Parameters: 2725 + dvalue - the double values to return 2726 . nmax - actual number of values retreived 2727 - set - PETSC_TRUE if found, else PETSC_FALSE 2728 2729 Level: beginner 2730 2731 .seealso: PetscOptionsGetInt(), PetscOptionsHasName(), 2732 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(), 2733 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2734 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2735 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2736 PetscOptionsFList(), PetscOptionsEList() 2737 @*/ 2738 PetscErrorCode PetscOptionsGetRealArray(PetscOptions options,const char pre[],const char name[],PetscReal dvalue[],PetscInt *nmax,PetscBool *set) 2739 { 2740 const char *svalue; 2741 char *value; 2742 PetscErrorCode ierr; 2743 PetscInt n = 0; 2744 PetscBool flag; 2745 PetscToken token; 2746 2747 PetscFunctionBegin; 2748 PetscValidCharPointer(name,3); 2749 PetscValidRealPointer(dvalue,4); 2750 PetscValidIntPointer(nmax,5); 2751 2752 ierr = PetscOptionsFindPair(options,pre,name,&svalue,&flag);CHKERRQ(ierr); 2753 if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; PetscFunctionReturn(0);} 2754 if (set) *set = PETSC_TRUE; 2755 ierr = PetscTokenCreate(svalue,',',&token);CHKERRQ(ierr); 2756 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2757 while (value && n < *nmax) { 2758 ierr = PetscOptionsStringToReal(value,dvalue++);CHKERRQ(ierr); 2759 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2760 n++; 2761 } 2762 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 2763 *nmax = n; 2764 PetscFunctionReturn(0); 2765 } 2766 2767 /*@C 2768 PetscOptionsGetScalarArray - Gets an array of scalars for a 2769 particular option in the database. The values must be separated with 2770 commas with no intervening spaces. 2771 2772 Not Collective 2773 2774 Input Parameters: 2775 + options - options database, use NULL for default global database 2776 . pre - string to prepend to each name or NULL 2777 . name - the option one is seeking 2778 - nmax - maximum number of values to retrieve 2779 2780 Output Parameters: 2781 + dvalue - the scalar values to return 2782 . nmax - actual number of values retreived 2783 - set - PETSC_TRUE if found, else PETSC_FALSE 2784 2785 Level: beginner 2786 2787 .seealso: PetscOptionsGetInt(), PetscOptionsHasName(), 2788 PetscOptionsGetString(), PetscOptionsGetIntArray(), PetscOptionsBool(), 2789 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2790 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2791 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2792 PetscOptionsFList(), PetscOptionsEList() 2793 @*/ 2794 PetscErrorCode PetscOptionsGetScalarArray(PetscOptions options,const char pre[],const char name[],PetscScalar dvalue[],PetscInt *nmax,PetscBool *set) 2795 { 2796 const char *svalue; 2797 char *value; 2798 PetscErrorCode ierr; 2799 PetscInt n = 0; 2800 PetscBool flag; 2801 PetscToken token; 2802 2803 PetscFunctionBegin; 2804 PetscValidCharPointer(name,3); 2805 PetscValidRealPointer(dvalue,4); 2806 PetscValidIntPointer(nmax,5); 2807 2808 ierr = PetscOptionsFindPair(options,pre,name,&svalue,&flag);CHKERRQ(ierr); 2809 if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; PetscFunctionReturn(0);} 2810 if (set) *set = PETSC_TRUE; 2811 ierr = PetscTokenCreate(svalue,',',&token);CHKERRQ(ierr); 2812 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2813 while (value && n < *nmax) { 2814 ierr = PetscOptionsStringToScalar(value,dvalue++);CHKERRQ(ierr); 2815 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2816 n++; 2817 } 2818 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 2819 *nmax = n; 2820 PetscFunctionReturn(0); 2821 } 2822 2823 /*@C 2824 PetscOptionsGetStringArray - Gets an array of string values for a particular 2825 option in the database. The values must be separated with commas with 2826 no intervening spaces. 2827 2828 Not Collective 2829 2830 Input Parameters: 2831 + options - options database, use NULL for default global database 2832 . pre - string to prepend to name or NULL 2833 . name - the option one is seeking 2834 - nmax - maximum number of strings 2835 2836 Output Parameter: 2837 + strings - location to copy strings 2838 - set - PETSC_TRUE if found, else PETSC_FALSE 2839 2840 Level: beginner 2841 2842 Notes: 2843 The user should pass in an array of pointers to char, to hold all the 2844 strings returned by this function. 2845 2846 The user is responsible for deallocating the strings that are 2847 returned. The Fortran interface for this routine is not supported. 2848 2849 .seealso: PetscOptionsGetInt(), PetscOptionsGetReal(), 2850 PetscOptionsHasName(), PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool(), 2851 PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(), 2852 PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(), 2853 PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(), 2854 PetscOptionsFList(), PetscOptionsEList() 2855 @*/ 2856 PetscErrorCode PetscOptionsGetStringArray(PetscOptions options,const char pre[],const char name[],char *strings[],PetscInt *nmax,PetscBool *set) 2857 { 2858 const char *svalue; 2859 char *value; 2860 PetscErrorCode ierr; 2861 PetscInt n = 0; 2862 PetscBool flag; 2863 PetscToken token; 2864 2865 PetscFunctionBegin; 2866 PetscValidCharPointer(name,3); 2867 PetscValidPointer(strings,4); 2868 PetscValidIntPointer(nmax,5); 2869 2870 ierr = PetscOptionsFindPair(options,pre,name,&svalue,&flag);CHKERRQ(ierr); 2871 if (!flag || !svalue) { if (set) *set = PETSC_FALSE; *nmax = 0; PetscFunctionReturn(0);} 2872 if (set) *set = PETSC_TRUE; 2873 ierr = PetscTokenCreate(svalue,',',&token);CHKERRQ(ierr); 2874 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2875 while (value && n < *nmax) { 2876 ierr = PetscStrallocpy(value,&strings[n]);CHKERRQ(ierr); 2877 ierr = PetscTokenFind(token,&value);CHKERRQ(ierr); 2878 n++; 2879 } 2880 ierr = PetscTokenDestroy(&token);CHKERRQ(ierr); 2881 *nmax = n; 2882 PetscFunctionReturn(0); 2883 } 2884 2885 /*@C 2886 PetscOptionsDeprecated - mark an option as deprecated, optionally replacing it with a new one 2887 2888 Prints a deprecation warning, unless an option is supplied to suppress. 2889 2890 Logically Collective 2891 2892 Input Parameters: 2893 + pre - string to prepend to name or NULL 2894 . oldname - the old, deprecated option 2895 . newname - the new option, or NULL if option is purely removed 2896 . version - a string describing the version of first deprecation, e.g. "3.9" 2897 - info - additional information string, or NULL. 2898 2899 Options Database Keys: 2900 . -options_suppress_deprecated_warnings - do not print deprecation warnings 2901 2902 Notes: 2903 Must be called between PetscOptionsBegin() (or PetscObjectOptionsBegin()) and PetscOptionsEnd(). 2904 Only the proces of rank zero that owns the PetscOptionsItems are argument (managed by PetscOptionsBegin() or 2905 PetscObjectOptionsBegin() prints the information 2906 If newname is provided, the old option is replaced. Otherwise, it remains 2907 in the options database. 2908 If an option is not replaced, the info argument should be used to advise the user 2909 on how to proceed. 2910 There is a limit on the length of the warning printed, so very long strings 2911 provided as info may be truncated. 2912 2913 Level: developer 2914 2915 .seealso: PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsScalar(), PetscOptionsBool(), PetscOptionsString(), PetscOptionsSetValue() 2916 2917 @*/ 2918 PetscErrorCode PetscOptionsDeprecated_Private(PetscOptionItems *PetscOptionsObject,const char oldname[],const char newname[],const char version[],const char info[]) 2919 { 2920 PetscErrorCode ierr; 2921 PetscBool found,quiet; 2922 const char *value; 2923 const char * const quietopt="-options_suppress_deprecated_warnings"; 2924 char msg[4096]; 2925 char *prefix = NULL; 2926 PetscOptions options = NULL; 2927 MPI_Comm comm = PETSC_COMM_SELF; 2928 2929 PetscFunctionBegin; 2930 PetscValidCharPointer(oldname,2); 2931 PetscValidCharPointer(version,4); 2932 if (PetscOptionsObject) { 2933 prefix = PetscOptionsObject->prefix; 2934 options = PetscOptionsObject->options; 2935 comm = PetscOptionsObject->comm; 2936 } 2937 ierr = PetscOptionsFindPair(options,prefix,oldname,&value,&found);CHKERRQ(ierr); 2938 if (found) { 2939 if (newname) { 2940 if (prefix) { 2941 ierr = PetscOptionsPrefixPush(options,prefix);CHKERRQ(ierr); 2942 } 2943 ierr = PetscOptionsSetValue(options,newname,value);CHKERRQ(ierr); 2944 if (prefix) { 2945 ierr = PetscOptionsPrefixPop(options);CHKERRQ(ierr); 2946 } 2947 ierr = PetscOptionsClearValue(options,oldname);CHKERRQ(ierr); 2948 } 2949 quiet = PETSC_FALSE; 2950 ierr = PetscOptionsGetBool(options,NULL,quietopt,&quiet,NULL);CHKERRQ(ierr); 2951 if (!quiet) { 2952 ierr = PetscStrcpy(msg,"** PETSc DEPRECATION WARNING ** : the option ");CHKERRQ(ierr); 2953 ierr = PetscStrcat(msg,oldname);CHKERRQ(ierr); 2954 ierr = PetscStrcat(msg," is deprecated as of version ");CHKERRQ(ierr); 2955 ierr = PetscStrcat(msg,version);CHKERRQ(ierr); 2956 ierr = PetscStrcat(msg," and will be removed in a future release.");CHKERRQ(ierr); 2957 if (newname) { 2958 ierr = PetscStrcat(msg," Please use the option ");CHKERRQ(ierr); 2959 ierr = PetscStrcat(msg,newname);CHKERRQ(ierr); 2960 ierr = PetscStrcat(msg," instead.");CHKERRQ(ierr); 2961 } 2962 if (info) { 2963 ierr = PetscStrcat(msg," ");CHKERRQ(ierr); 2964 ierr = PetscStrcat(msg,info);CHKERRQ(ierr); 2965 } 2966 ierr = PetscStrcat(msg," (Silence this warning with ");CHKERRQ(ierr); 2967 ierr = PetscStrcat(msg,quietopt);CHKERRQ(ierr); 2968 ierr = PetscStrcat(msg,")\n");CHKERRQ(ierr); 2969 ierr = PetscPrintf(comm,msg);CHKERRQ(ierr); 2970 } 2971 } 2972 PetscFunctionReturn(0); 2973 } 2974