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