Actual source code: opts.c

  1: /*
  2:       EPS routines related to options that can be set via the command-line 
  3:       or procedurally.
  4: */
 5:  #include src/eps/epsimpl.h

  7: EXTERN PetscTruth EPSRegisterAllCalled;

 11: /*@
 12:    EPSSetFromOptions - Sets EPS options from the options database.
 13:    This routine must be called before EPSSetUp() if the user is to be 
 14:    allowed to set the solver type. 

 16:    Collective on EPS

 18:    Input Parameters:
 19: .  eps - the eigensolver context

 21:    Notes:  
 22:    To see all options, run your program with the -help option.

 24:    Level: beginner

 26: .seealso: 
 27: @*/
 28: PetscErrorCode EPSSetFromOptions(EPS eps)
 29: {
 30:   PetscErrorCode ierr;
 31:   char           type[256];
 32:   PetscTruth     flg;
 33:   const char     *orth_list[2] = { "mgs" , "cgs" };
 34:   const char     *ref_list[3] = { "never" , "ifneeded", "always" };
 35:   PetscReal      eta;
 36:   EPSOrthogonalizationType orth_type;
 37:   EPSOrthogonalizationRefinementType ref_type;

 41:   if (!EPSRegisterAllCalled) {EPSRegisterAll(PETSC_NULL);}
 42:   PetscOptionsBegin(eps->comm,eps->prefix,"Eigenproblem Solver (EPS) Options","EPS");
 43:     PetscOptionsList("-eps_type","Eigenproblem Solver method","EPSSetType",EPSList,(char*)(eps->type_name?eps->type_name:EPSARNOLDI),type,256,&flg);
 44:     if (flg) {
 45:       EPSSetType(eps,type);
 46:     }
 47:     /*
 48:       Set the type if it was never set.
 49:     */
 50:     if (!eps->type_name) {
 51:       EPSSetType(eps,EPSARNOLDI);
 52:     }

 54:     PetscOptionsLogicalGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
 55:     if (flg) {EPSSetProblemType(eps,EPS_HEP);}
 56:     PetscOptionsLogicalGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
 57:     if (flg) {EPSSetProblemType(eps,EPS_GHEP);}
 58:     PetscOptionsLogicalGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
 59:     if (flg) {EPSSetProblemType(eps,EPS_NHEP);}
 60:     PetscOptionsLogicalGroupEnd("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
 61:     if (flg) {EPSSetProblemType(eps,EPS_GNHEP);}

 63:     orth_type = eps->orthog_type;
 64:     PetscOptionsEList("-eps_orthog_type","Orthogonalization method","EPSSetOrthogonalization",orth_list,2,orth_list[orth_type],(int*)&orth_type,&flg);
 65:     ref_type = eps->orthog_ref;
 66:     PetscOptionsEList("-eps_orthog_refinement","Iterative refinement mode during orthogonalization","EPSSetOrthogonalizationRefinement",ref_list,3,ref_list[ref_type],(int*)&ref_type,&flg);
 67:     eta = eps->orthog_eta;
 68:     PetscOptionsReal("-eps_orthog_eta","Parameter of iterative refinement during orthogonalization","EPSSetOrthogonalizationRefinement",eta,&eta,PETSC_NULL);
 69:     EPSSetOrthogonalization(eps,orth_type,ref_type,eta);

 71:     PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&eps->max_it,PETSC_NULL);
 72:     PetscOptionsReal("-eps_tol","Tolerance","KSPSetTolerances",eps->tol,&eps->tol,PETSC_NULL);
 73:     PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&eps->nev,&flg);
 74:     if( eps->nev<1 ) SETERRQ(1,"Illegal value for option -eps_nev. Must be > 0");
 75:     PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&eps->ncv,&flg);
 76:     if( flg && eps->ncv<1 ) SETERRQ(1,"Illegal value for option -eps_ncv. Must be > 0");

 78:     /* -----------------------------------------------------------------------*/
 79:     /*
 80:       Cancels all monitors hardwired into code before call to EPSSetFromOptions()
 81:     */
 82:     PetscOptionsName("-eps_cancelmonitors","Remove any hardwired monitor routines","EPSClearMonitor",&flg);
 83:     if (flg) {
 84:       EPSClearMonitor(eps);
 85:     }
 86:     /*
 87:       Prints approximate eigenvalues and error estimates at each iteration
 88:     */
 89:     PetscOptionsName("-eps_monitor","Monitor approximate eigenvalues and error estimates","EPSSetMonitor",&flg);
 90:     if (flg) {
 91:       EPSSetMonitor(eps,EPSDefaultMonitor,PETSC_NULL);
 92:     }
 93:   /* -----------------------------------------------------------------------*/
 94:     PetscOptionsLogicalGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
 95:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);}
 96:     PetscOptionsLogicalGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
 97:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE);}
 98:     PetscOptionsLogicalGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
 99:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL);}
100:     PetscOptionsLogicalGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
101:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);}
102:     PetscOptionsLogicalGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
103:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY);}
104:     PetscOptionsLogicalGroupEnd("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
105:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY);}

107:     PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
108:     PetscOptionsName("-eps_view_binary","Save the matrices associated to the eigenproblem","EPSSetFromOptions",0);
109:     PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);

111:     if (eps->ops->setfromoptions) {
112:       (*eps->ops->setfromoptions)(eps);
113:     }
114:   PetscOptionsEnd();

116:   STSetFromOptions(eps->OP);
117:   return(0);
118: }

122: /*@
123:    EPSGetTolerances - Gets the tolerance and maximum
124:    iteration count used by the default EPS convergence tests. 

126:    Not Collective

128:    Input Parameter:
129: .  eps - the eigensolver context
130:   
131:    Output Parameters:
132: +  tol - the convergence tolerance
133: -  maxits - maximum number of iterations

135:    Notes:
136:    The user can specify PETSC_NULL for any parameter that is not needed.

138:    Level: intermediate

140: .seealso: EPSSetTolerances()
141: @*/
142: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,int *maxits)
143: {
146:   if (tol)    *tol    = eps->tol;
147:   if (maxits) *maxits = eps->max_it;
148:   return(0);
149: }

153: /*@
154:    EPSSetTolerances - Sets the tolerance and maximum
155:    iteration count used by the default EPS convergence testers. 

157:    Collective on EPS

159:    Input Parameters:
160: +  eps - the eigensolver context
161: .  tol - the convergence tolerance
162: -  maxits - maximum number of iterations to use

164:    Options Database Keys:
165: +  -eps_tol <tol> - Sets the convergence tolerance
166: -  -eps_max_it <maxits> - Sets the maximum number of iterations allowed

168:    Notes:
169:    Use PETSC_DEFAULT to retain the default value of any of the tolerances.

171:    Level: intermediate

173: .seealso: EPSGetTolerances()
174: @*/
175: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,int maxits)
176: {
179:   if (tol != PETSC_DEFAULT)    eps->tol    = tol;
180:   if (maxits != PETSC_DEFAULT) eps->max_it = maxits;
181:   return(0);
182: }

186: /*@
187:    EPSGetDimensions - Gets the number of eigenvalues to compute
188:    and the dimension of the subspace.

190:    Not Collective

192:    Input Parameter:
193: .  eps - the eigensolver context
194:   
195:    Output Parameters:
196: +  nev - number of eigenvalues to compute
197: -  ncv - the maximum dimension of the subspace to be used by the solver

199:    Notes:
200:    The user can specify PETSC_NULL for any parameter that is not needed.

202:    Level: intermediate

204: .seealso: EPSSetDimensions()
205: @*/
206: PetscErrorCode EPSGetDimensions(EPS eps,int *nev,int *ncv)
207: {
210:   if( nev )   *nev = eps->nev;
211:   if( ncv )   *ncv = eps->ncv;
212:   return(0);
213: }

217: /*@
218:    EPSSetDimensions - Sets the number of eigenvalues to compute
219:    and the dimension of the subspace.

221:    Collective on EPS

223:    Input Parameters:
224: +  eps - the eigensolver context
225: .  nev - number of eigenvalues to compute
226: -  ncv - the maximum dimension of the subspace to be used by the solver

228:    Options Database Keys:
229: +  -eps_nev <nev> - Sets the number of eigenvalues
230: -  -eps_ncv <ncv> - Sets the dimension of the subspace

232:    Notes:
233:    Use PETSC_DEFAULT to retain the previous value of any parameter.

235:    Use PETSC_DECIDE for ncv to assign a reasonably good value, which is 
236:    dependent on the solution method.

238:    Level: intermediate

240: .seealso: EPSGetDimensions()
241: @*/
242: PetscErrorCode EPSSetDimensions(EPS eps,int nev,int ncv)
243: {

247:   if( nev != PETSC_DEFAULT ) {
248:     if (nev<1) SETERRQ(1,"Illegal value of nev. Must be > 0");
249:     eps->nev = nev;
250:     eps->setupcalled = 0;
251:   }
252:   if( ncv != PETSC_DEFAULT ) {
253:     if( ncv == PETSC_DECIDE ) eps->ncv = 0;
254:     else {
255:       if (ncv<1) SETERRQ(1,"Illegal value of ncv. Must be > 0");
256:       eps->ncv = ncv;
257:     }
258:     eps->setupcalled = 0;
259:   }
260:   return(0);
261: }

265: /*@
266:     EPSSetWhichEigenpairs - Specifies which portion of the spectrum is 
267:     to be sought.

269:     Collective on EPS

271:     Input Parameter:
272: .   eps - eigensolver context obtained from EPSCreate()

274:     Output Parameter:
275: .   which - the portion of the spectrum to be sought

277:     Possible values:
278:     The parameter 'which' can have one of these values:
279:     
280: +     EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
281: .     EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
282: .     EPS_LARGEST_REAL - largest real parts
283: .     EPS_SMALLEST_REAL - smallest real parts
284: .     EPS_LARGEST_IMAGINARY - largest imaginary parts
285: -     EPS_SMALLEST_IMAGINARY - smallest imaginary parts

287:     Options Database Keys:
288: +   -eps_largest_magnitude - Sets largest eigenvalues in magnitude
289: .   -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
290: .   -eps_largest_real - Sets largest real parts
291: .   -eps_smallest_real - Sets smallest real parts
292: .   -eps_largest_imaginary - Sets largest imaginary parts in magnitude
293: -   -eps_smallest_imaginary - Sets smallest imaginary parts in magnitude

295:     Notes:
296:     Not all eigensolvers implemented in EPS account for all the possible values
297:     stated above. Also, some values make sense only for certain types of 
298:     problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
299:     and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part 
300:     for eigenvalue selection.     
301:     
302:     Level: intermediate

304: .seealso: EPSGetWhichEigenpairs(), EPSSortEigenvalues()
305: @*/
306: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
307: {
310:   eps->which = which;
311:   return(0);
312: }

316: /*@
317:     EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be 
318:     sought.

320:     Not Collective

322:     Input Parameter:
323: .   eps - eigensolver context obtained from EPSCreate()

325:     Output Parameter:
326: .   which - the portion of the spectrum to be sought

328:     Notes:
329:     See EPSSetWhichEigenpairs() for possible values of which

331:     Level: intermediate

333: .seealso: EPSSetWhichEigenpairs()
334: @*/
335: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
336: {
339:   *which = eps->which;
340:   return(0);
341: }

345: /*@
346:    EPSSetProblemType - Specifies the type of the eigenvalue problem.

348:    Collective on EPS

350:    Input Parameters:
351: +  eps      - the eigensolver context
352: -  type     - a known type of eigenvalue problem 

354:    Options Database Keys:
355: +  -eps_hermitian - Hermitian eigenvalue problem
356: .  -eps_gen_hermitian - generalized Hermitian eigenvalue problem
357: .  -eps_non_hermitian - non-Hermitian eigenvalue problem
358: -  -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem 
359:     
360:    Note:  
361:    Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
362:    (EPS_NHEP), generalized Hermitian (EPS_GHEP) and generalized non-Hermitian 
363:    (EPS_GNHEP).

365:    This function must be used to instruct SLEPc to exploit symmetry. If no
366:    problem type is specified, by default a non-Hermitian problem is assumed
367:    (either standard or generalized). If the user knows that the problem is
368:    Hermitian (i.e. A=A^H) of generalized Hermitian (i.e. A=A^H, B=B^H, and 
369:    B positive definite) then it is recommended to set the problem type so
370:    that eigensolver can exploit these properties. 

372:    Level: beginner

374: .seealso: EPSSetOperators(), EPSSetType(), EPSProblemType
375: @*/
376: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
377: {
378:   PetscErrorCode ierr;
379:   Mat            A,B;


384:   STGetOperators(eps->OP,&A,&B);
385:   if (!A) { SETERRQ(1,"Must call EPSSetOperators() first"); }

387:   switch (type) {
388:     case EPS_HEP:
389:       eps->isgeneralized = PETSC_FALSE;
390:       eps->ishermitian = PETSC_TRUE;
391:       STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
392:       break;
393:     case EPS_NHEP:
394:       eps->isgeneralized = PETSC_FALSE;
395:       eps->ishermitian = PETSC_FALSE;
396:       STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
397:       break;
398:     case EPS_GHEP:
399:       eps->isgeneralized = PETSC_TRUE;
400:       eps->ishermitian = PETSC_TRUE;
401:       STSetBilinearForm(eps->OP,STINNER_B_HERMITIAN);
402:       break;
403:     case EPS_GNHEP:
404:       eps->isgeneralized = PETSC_TRUE;
405:       eps->ishermitian = PETSC_FALSE;
406:       STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
407:       break;
408: /*
409:     case EPS_CSEP: 
410:       eps->isgeneralized = PETSC_FALSE;
411:       eps->ishermitian = PETSC_FALSE;
412:       STSetBilinearForm(eps->OP,STINNER_SYMMETRIC);
413:       break;
414:     case EPS_GCSEP:
415:       eps->isgeneralized = PETSC_TRUE;
416:       eps->ishermitian = PETSC_FALSE;
417:       STSetBilinearForm(eps->OP,STINNER_B_SYMMETRIC);
418:       break;
419: */
420:     default:
421:       SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
422:   }
423:   eps->problem_type = type;

425:   return(0);
426: }

430: /*@
431:    EPSGetProblemType - Gets the problem type from the EPS object.

433:    Not Collective

435:    Input Parameter:
436: .  eps - the eigensolver context 

438:    Output Parameter:
439: .  type - name of EPS problem type 

441:    Level: intermediate

443: .seealso: EPSSetProblemType()
444: @*/
445: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
446: {
449:   *type = eps->problem_type;
450:   return(0);
451: }

455: /*@
456:    EPSSetOrthogonalization - Specifies the type of orthogonalization technique
457:    to be used inside the eigensolver.

459:    Collective on EPS

461:    Input Parameters:
462: +  eps        - the eigensolver context 
463: .  type       - a known type of orthogonalization
464: .  refinement - type of refinement
465: -  eta        - parameter for dynamic refinement

467:    Options Database Keys:
468: +  -eps_orthog_type <type> -  Where <type> is cgs for Classical Gram-Schmidt
469:                               orthogonalization (default)
470:                               or mgs for Modified Gram-Schmidt orthogonalization
471: .  -eps_orthog_refinement <type> -  Where <type> is one of never, ifneeded
472:                               (default) or always 
473: -  -eps_orthog_eta <eta> -  For setting the value of eta (or PETSC_DEFAULT)
474:     
475:    Notes:  
476:    The value of eta is used only when refinement type is "ifneeded". 

478:    The default orthogonalization technique 
479:    works well for most problems. MGS is numerically more robust than CGS,
480:    but CGS may give better scalability.

482:    Level: intermediate

484: .seealso: EPSOrthogonalize(), EPSGetOrthogonalization()
485: @*/
486: PetscErrorCode EPSSetOrthogonalization(EPS eps,EPSOrthogonalizationType type, EPSOrthogonalizationRefinementType refinement, PetscReal eta)
487: {
490:   switch (type) {
491:     case EPS_CGS_ORTH:
492:     case EPS_MGS_ORTH:
493:       eps->orthog_type = type;
494:       break;
495:     default:
496:       SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown orthogonalization type");
497:   }
498:   switch (refinement) {
499:     case EPS_ORTH_REFINE_NEVER:
500:     case EPS_ORTH_REFINE_IFNEEDED:
501:     case EPS_ORTH_REFINE_ALWAYS:
502:       eps->orthog_ref = refinement;
503:       break;
504:     default:
505:       SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown refinement type");
506:   }
507:   if (eta == PETSC_DEFAULT) eps->orthog_eta = 0.7071;
508:   else eps->orthog_eta = eta;
509:   return(0);
510: }

514: /*@
515:    EPSGetOrthogonalization - Gets the orthogonalization type from the 
516:    EPS object.

518:    Not Collective

520:    Input Parameter:
521: .  eps - Eigensolver context 

523:    Output Parameter:
524: +  type       - type of orthogonalization technique
525: .  refinement - type of refinement
526: -  eta        - parameter for dynamic refinement

528:    Level: intermediate

530: .seealso: EPSSetOrthogonalization()
531: @*/
532: PetscErrorCode EPSGetOrthogonalization(EPS eps,EPSOrthogonalizationType *type,EPSOrthogonalizationRefinementType *refinement, PetscReal *eta)
533: {
536:   if (type) *type = eps->orthog_type;
537:   if (refinement) *refinement = eps->orthog_ref;
538:   if (eta) *eta = eps->orthog_eta;
539:   return(0);
540: }

544: /*@C
545:    EPSSetOptionsPrefix - Sets the prefix used for searching for all 
546:    EPS options in the database.

548:    Collective on EPS

550:    Input Parameters:
551: +  eps - the eigensolver context
552: -  prefix - the prefix string to prepend to all EPS option requests

554:    Notes:
555:    A hyphen (-) must NOT be given at the beginning of the prefix name.
556:    The first character of all runtime options is AUTOMATICALLY the
557:    hyphen.

559:    For example, to distinguish between the runtime options for two
560:    different EPS contexts, one could call
561: .vb
562:       EPSSetOptionsPrefix(eps1,"eig1_")
563:       EPSSetOptionsPrefix(eps2,"eig2_")
564: .ve

566:    Level: advanced

568: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
569: @*/
570: PetscErrorCode EPSSetOptionsPrefix(EPS eps,char *prefix)
571: {
572:   PetscErrorCode ierr;
575:   PetscObjectSetOptionsPrefix((PetscObject)eps, prefix);
576:   STSetOptionsPrefix(eps->OP,prefix);
577:   return(0);
578: }
579: 
582: /*@C
583:    EPSAppendOptionsPrefix - Appends to the prefix used for searching for all 
584:    EPS options in the database.

586:    Collective on EPS

588:    Input Parameters:
589: +  eps - the eigensolver context
590: -  prefix - the prefix string to prepend to all EPS option requests

592:    Notes:
593:    A hyphen (-) must NOT be given at the beginning of the prefix name.
594:    The first character of all runtime options is AUTOMATICALLY the hyphen.

596:    Level: advanced

598: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
599: @*/
600: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,char *prefix)
601: {
602:   PetscErrorCode ierr;
605:   PetscObjectAppendOptionsPrefix((PetscObject)eps, prefix);
606:   STAppendOptionsPrefix(eps->OP,prefix);
607:   return(0);
608: }

612: /*@C
613:    EPSGetOptionsPrefix - Gets the prefix used for searching for all 
614:    EPS options in the database.

616:    Not Collective

618:    Input Parameters:
619: .  eps - the eigensolver context

621:    Output Parameters:
622: .  prefix - pointer to the prefix string used is returned

624:    Notes: On the fortran side, the user should pass in a string 'prefix' of
625:    sufficient length to hold the prefix.

627:    Level: advanced

629: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
630: @*/
631: PetscErrorCode EPSGetOptionsPrefix(EPS eps,char **prefix)
632: {
633:   PetscErrorCode ierr;
636:   PetscObjectGetOptionsPrefix((PetscObject)eps, prefix);
637:   return(0);
638: }