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