Actual source code: pepopts.c

slepc-3.5.2 2014-10-10
Report Typos and Errors
  1: /*
  2:       PEP routines related to options that can be set via the command-line
  3:       or procedurally.

  5:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  6:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  7:    Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain

  9:    This file is part of SLEPc.

 11:    SLEPc is free software: you can redistribute it and/or modify it under  the
 12:    terms of version 3 of the GNU Lesser General Public License as published by
 13:    the Free Software Foundation.

 15:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 16:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 17:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 18:    more details.

 20:    You  should have received a copy of the GNU Lesser General  Public  License
 21:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 22:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 23: */

 25: #include <slepc-private/pepimpl.h>       /*I "slepcpep.h" I*/

 29: /*@
 30:    PEPSetFromOptions - Sets PEP options from the options database.
 31:    This routine must be called before PEPSetUp() if the user is to be
 32:    allowed to set the solver type.

 34:    Collective on PEP

 36:    Input Parameters:
 37: .  pep - the polynomial eigensolver context

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

 42:    Level: beginner
 43: @*/
 44: PetscErrorCode PEPSetFromOptions(PEP pep)
 45: {
 46:   PetscErrorCode   ierr;
 47:   char             type[256],monfilename[PETSC_MAX_PATH_LEN];
 48:   PetscBool        flg,flg1,flg2,flg3,flg4;
 49:   PetscReal        r,t;
 50:   PetscScalar      s;
 51:   PetscInt         i,j,k;
 52:   PetscViewer      monviewer;
 53:   SlepcConvMonitor ctx;

 57:   if (!PEPRegisterAllCalled) { PEPRegisterAll(); }
 58:   PetscObjectOptionsBegin((PetscObject)pep);
 59:     PetscOptionsFList("-pep_type","Polynomial Eigenvalue Problem method","PEPSetType",PEPList,(char*)(((PetscObject)pep)->type_name?((PetscObject)pep)->type_name:PEPTOAR),type,256,&flg);
 60:     if (flg) {
 61:       PEPSetType(pep,type);
 62:     } else if (!((PetscObject)pep)->type_name) {
 63:       PEPSetType(pep,PEPTOAR);
 64:     }

 66:     PetscOptionsBoolGroupBegin("-pep_general","general polynomial eigenvalue problem","PEPSetProblemType",&flg);
 67:     if (flg) { PEPSetProblemType(pep,PEP_GENERAL); }
 68:     PetscOptionsBoolGroup("-pep_hermitian","hermitian polynomial eigenvalue problem","PEPSetProblemType",&flg);
 69:     if (flg) { PEPSetProblemType(pep,PEP_HERMITIAN); }
 70:     PetscOptionsBoolGroupEnd("-pep_gyroscopic","gyroscopic polynomial eigenvalue problem","PEPSetProblemType",&flg);
 71:     if (flg) { PEPSetProblemType(pep,PEP_GYROSCOPIC); }

 73:     PetscOptionsEnum("-pep_scale","Scaling strategy","PEPSetScale",PEPScaleTypes,(PetscEnum)pep->scale,(PetscEnum*)&pep->scale,NULL);

 75:     r = pep->sfactor;
 76:     PetscOptionsReal("-pep_scale_factor","Scale factor","PEPSetScale",pep->sfactor,&r,&flg1);
 77:     j = pep->sits;
 78:     PetscOptionsInt("-pep_scale_its","Number of iterations in diagonal scaling","PEPSetScale",pep->sits,&j,&flg2);
 79:     t = pep->slambda;
 80:     PetscOptionsReal("-pep_scale_lambda","Estimate of eigenvalue (modulus) for diagonal scaling","PEPSetScale",pep->slambda,&t,&flg3);
 81:     if (flg1 || flg2 || flg3) {
 82:       PEPSetScale(pep,pep->scale,r,j,t);
 83:     }

 85:     PetscOptionsEnum("-pep_extract","Extraction method","PEPSetExtract",PEPExtractTypes,(PetscEnum)pep->extract,(PetscEnum*)&pep->extract,NULL);

 87:     PetscOptionsEnum("-pep_refine","Iterative refinement method","PEPSetRefine",PEPRefineTypes,(PetscEnum)pep->refine,(PetscEnum*)&pep->refine,NULL);

 89:     i = pep->npart;
 90:     PetscOptionsInt("-pep_refine_partitions","Number of partitions of the communicator for iterative refinement","PEPSetRefine",pep->npart,&i,&flg1);
 91:     r = pep->rtol;
 92:     PetscOptionsReal("-pep_refine_tol","Tolerance for iterative refinement","PEPSetRefine",pep->rtol,&r,&flg2);
 93:     j = pep->rits;
 94:     PetscOptionsInt("-pep_refine_its","Maximum number of iterations for iterative refinement","PEPSetRefine",pep->rits,&j,&flg3);
 95:     flg = pep->schur;
 96:     PetscOptionsBool("-pep_refine_schur","Use Schur complement for iterative refinement","PEPSetRefine",pep->schur,&flg,&flg4);
 97:     if (flg1 || flg2 || flg3 || flg4) {
 98:       PEPSetRefine(pep,pep->refine,i,r,j,flg);
 99:     }

101:     i = pep->max_it? pep->max_it: PETSC_DEFAULT;
102:     PetscOptionsInt("-pep_max_it","Maximum number of iterations","PEPSetTolerances",pep->max_it,&i,&flg1);
103:     r = pep->tol;
104:     PetscOptionsReal("-pep_tol","Tolerance","PEPSetTolerances",pep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:pep->tol,&r,&flg2);
105:     if (flg1 || flg2) {
106:       PEPSetTolerances(pep,r,i);
107:     }

109:     PetscOptionsBoolGroupBegin("-pep_conv_eig","Relative error convergence test","PEPSetConvergenceTest",&flg);
110:     if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_EIG); }
111:     PetscOptionsBoolGroup("-pep_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","PEPSetConvergenceTest",&flg);
112:     if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_NORM); }
113:     PetscOptionsBoolGroup("-pep_conv_abs","Absolute error convergence test","PEPSetConvergenceTest",&flg);
114:     if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_ABS); }
115:     PetscOptionsBoolGroupEnd("-pep_conv_user","User-defined convergence test","PEPSetConvergenceTest",&flg);
116:     if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_USER); }

118:     i = pep->nev;
119:     PetscOptionsInt("-pep_nev","Number of eigenvalues to compute","PEPSetDimensions",pep->nev,&i,&flg1);
120:     j = pep->ncv? pep->ncv: PETSC_DEFAULT;
121:     PetscOptionsInt("-pep_ncv","Number of basis vectors","PEPSetDimensions",pep->ncv,&j,&flg2);
122:     k = pep->mpd? pep->mpd: PETSC_DEFAULT;
123:     PetscOptionsInt("-pep_mpd","Maximum dimension of projected problem","PEPSetDimensions",pep->mpd,&k,&flg3);
124:     if (flg1 || flg2 || flg3) {
125:       PEPSetDimensions(pep,i,j,k);
126:     }

128:     PetscOptionsScalar("-pep_target","Value of the target","PEPSetTarget",pep->target,&s,&flg);
129:     if (flg) {
130:       PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE);
131:       PEPSetTarget(pep,s);
132:     }

134:     PetscOptionsEnum("-pep_basis","Polynomial basis","PEPSetBasis",PEPBasisTypes,(PetscEnum)pep->basis,(PetscEnum*)&pep->basis,NULL);

136:     /* -----------------------------------------------------------------------*/
137:     /*
138:       Cancels all monitors hardwired into code before call to PEPSetFromOptions()
139:     */
140:     flg  = PETSC_FALSE;
141:     PetscOptionsBool("-pep_monitor_cancel","Remove any hardwired monitor routines","PEPMonitorCancel",flg,&flg,NULL);
142:     if (flg) {
143:       PEPMonitorCancel(pep);
144:     }
145:     /*
146:       Prints approximate eigenvalues and error estimates at each iteration
147:     */
148:     PetscOptionsString("-pep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","PEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
149:     if (flg) {
150:       PetscViewerASCIIOpen(PetscObjectComm((PetscObject)pep),monfilename,&monviewer);
151:       PEPMonitorSet(pep,PEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
152:     }
153:     PetscOptionsString("-pep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","PEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
154:     if (flg) {
155:       PetscNew(&ctx);
156:       PetscViewerASCIIOpen(PetscObjectComm((PetscObject)pep),monfilename,&ctx->viewer);
157:       PEPMonitorSet(pep,PEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
158:     }
159:     PetscOptionsString("-pep_monitor_all","Monitor approximate eigenvalues and error estimates","PEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
160:     if (flg) {
161:       PetscViewerASCIIOpen(PetscObjectComm((PetscObject)pep),monfilename,&monviewer);
162:       PEPMonitorSet(pep,PEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
163:       PEPSetTrackAll(pep,PETSC_TRUE);
164:     }
165:     flg = PETSC_FALSE;
166:     PetscOptionsBool("-pep_monitor_lg","Monitor first unconverged approximate error estimate graphically","PEPMonitorSet",flg,&flg,NULL);
167:     if (flg) {
168:       PEPMonitorSet(pep,PEPMonitorLG,NULL,NULL);
169:     }
170:     flg = PETSC_FALSE;
171:     PetscOptionsBool("-pep_monitor_lg_all","Monitor error estimates graphically","PEPMonitorSet",flg,&flg,NULL);
172:     if (flg) {
173:       PEPMonitorSet(pep,PEPMonitorLGAll,NULL,NULL);
174:       PEPSetTrackAll(pep,PETSC_TRUE);
175:     }
176:   /* -----------------------------------------------------------------------*/

178:     PetscOptionsBoolGroupBegin("-pep_largest_magnitude","compute largest eigenvalues in magnitude","PEPSetWhichEigenpairs",&flg);
179:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_LARGEST_MAGNITUDE); }
180:     PetscOptionsBoolGroup("-pep_smallest_magnitude","compute smallest eigenvalues in magnitude","PEPSetWhichEigenpairs",&flg);
181:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_SMALLEST_MAGNITUDE); }
182:     PetscOptionsBoolGroup("-pep_largest_real","compute largest real parts","PEPSetWhichEigenpairs",&flg);
183:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_LARGEST_REAL); }
184:     PetscOptionsBoolGroup("-pep_smallest_real","compute smallest real parts","PEPSetWhichEigenpairs",&flg);
185:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_SMALLEST_REAL); }
186:     PetscOptionsBoolGroup("-pep_largest_imaginary","compute largest imaginary parts","PEPSetWhichEigenpairs",&flg);
187:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_LARGEST_IMAGINARY); }
188:     PetscOptionsBoolGroup("-pep_smallest_imaginary","compute smallest imaginary parts","PEPSetWhichEigenpairs",&flg);
189:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_SMALLEST_IMAGINARY); }
190:     PetscOptionsBoolGroup("-pep_target_magnitude","compute nearest eigenvalues to target","PEPSetWhichEigenpairs",&flg);
191:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE); }
192:     PetscOptionsBoolGroup("-pep_target_real","compute eigenvalues with real parts close to target","PEPSetWhichEigenpairs",&flg);
193:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_TARGET_REAL); }
194:     PetscOptionsBoolGroupEnd("-pep_target_imaginary","compute eigenvalues with imaginary parts close to target","PEPSetWhichEigenpairs",&flg);
195:     if (flg) { PEPSetWhichEigenpairs(pep,PEP_TARGET_IMAGINARY); }

197:     PetscOptionsName("-pep_view","Print detailed information on solver used","PEPView",0);
198:     PetscOptionsName("-pep_plot_eigs","Make a plot of the computed eigenvalues","PEPSolve",0);

200:     if (pep->ops->setfromoptions) {
201:       (*pep->ops->setfromoptions)(pep);
202:     }
203:     PetscObjectProcessOptionsHandlers((PetscObject)pep);
204:   PetscOptionsEnd();

206:   if (!pep->V) { PEPGetBV(pep,&pep->V); }
207:   BVSetFromOptions(pep->V);
208:   if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
209:   RGSetFromOptions(pep->rg);
210:   if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
211:   DSSetFromOptions(pep->ds);
212:   if (!pep->st) { PEPGetST(pep,&pep->st); }
213:   STSetFromOptions(pep->st);
214:   PetscRandomSetFromOptions(pep->rand);
215:   return(0);
216: }

220: /*@
221:    PEPGetTolerances - Gets the tolerance and maximum iteration count used
222:    by the PEP convergence tests.

224:    Not Collective

226:    Input Parameter:
227: .  pep - the polynomial eigensolver context

229:    Output Parameters:
230: +  tol - the convergence tolerance
231: -  maxits - maximum number of iterations

233:    Notes:
234:    The user can specify NULL for any parameter that is not needed.

236:    Level: intermediate

238: .seealso: PEPSetTolerances()
239: @*/
240: PetscErrorCode PEPGetTolerances(PEP pep,PetscReal *tol,PetscInt *maxits)
241: {
244:   if (tol)    *tol    = pep->tol;
245:   if (maxits) *maxits = pep->max_it;
246:   return(0);
247: }

251: /*@
252:    PEPSetTolerances - Sets the tolerance and maximum iteration count used
253:    by the PEP convergence tests.

255:    Logically Collective on PEP

257:    Input Parameters:
258: +  pep - the polynomial eigensolver context
259: .  tol - the convergence tolerance
260: -  maxits - maximum number of iterations to use

262:    Options Database Keys:
263: +  -pep_tol <tol> - Sets the convergence tolerance
264: -  -pep_max_it <maxits> - Sets the maximum number of iterations allowed

266:    Notes:
267:    Use PETSC_DEFAULT for either argument to assign a reasonably good value.

269:    Level: intermediate

271: .seealso: PEPGetTolerances()
272: @*/
273: PetscErrorCode PEPSetTolerances(PEP pep,PetscReal tol,PetscInt maxits)
274: {
279:   if (tol == PETSC_DEFAULT) {
280:     pep->tol   = PETSC_DEFAULT;
281:     pep->state = PEP_STATE_INITIAL;
282:   } else {
283:     if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
284:     pep->tol = tol;
285:   }
286:   if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
287:     pep->max_it = 0;
288:     pep->state  = PEP_STATE_INITIAL;
289:   } else {
290:     if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
291:     pep->max_it = maxits;
292:   }
293:   return(0);
294: }

298: /*@
299:    PEPGetDimensions - Gets the number of eigenvalues to compute
300:    and the dimension of the subspace.

302:    Not Collective

304:    Input Parameter:
305: .  pep - the polynomial eigensolver context

307:    Output Parameters:
308: +  nev - number of eigenvalues to compute
309: .  ncv - the maximum dimension of the subspace to be used by the solver
310: -  mpd - the maximum dimension allowed for the projected problem

312:    Notes:
313:    The user can specify NULL for any parameter that is not needed.

315:    Level: intermediate

317: .seealso: PEPSetDimensions()
318: @*/
319: PetscErrorCode PEPGetDimensions(PEP pep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
320: {
323:   if (nev) *nev = pep->nev;
324:   if (ncv) *ncv = pep->ncv;
325:   if (mpd) *mpd = pep->mpd;
326:   return(0);
327: }

331: /*@
332:    PEPSetDimensions - Sets the number of eigenvalues to compute
333:    and the dimension of the subspace.

335:    Logically Collective on PEP

337:    Input Parameters:
338: +  pep - the polynomial eigensolver context
339: .  nev - number of eigenvalues to compute
340: .  ncv - the maximum dimension of the subspace to be used by the solver
341: -  mpd - the maximum dimension allowed for the projected problem

343:    Options Database Keys:
344: +  -pep_nev <nev> - Sets the number of eigenvalues
345: .  -pep_ncv <ncv> - Sets the dimension of the subspace
346: -  -pep_mpd <mpd> - Sets the maximum projected dimension

348:    Notes:
349:    Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
350:    dependent on the solution method.

352:    The parameters ncv and mpd are intimately related, so that the user is advised
353:    to set one of them at most. Normal usage is that
354:    (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
355:    (b) in cases where nev is large, the user sets mpd.

357:    The value of ncv should always be between nev and (nev+mpd), typically
358:    ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
359:    a smaller value should be used.

361:    Level: intermediate

363: .seealso: PEPGetDimensions()
364: @*/
365: PetscErrorCode PEPSetDimensions(PEP pep,PetscInt nev,PetscInt ncv,PetscInt mpd)
366: {
372:   if (nev<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
373:   pep->nev = nev;
374:   if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
375:     pep->ncv = 0;
376:   } else {
377:     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
378:     pep->ncv = ncv;
379:   }
380:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
381:     pep->mpd = 0;
382:   } else {
383:     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
384:     pep->mpd = mpd;
385:   }
386:   pep->state = PEP_STATE_INITIAL;
387:   return(0);
388: }

392: /*@
393:     PEPSetWhichEigenpairs - Specifies which portion of the spectrum is
394:     to be sought.

396:     Logically Collective on PEP

398:     Input Parameters:
399: +   pep   - eigensolver context obtained from PEPCreate()
400: -   which - the portion of the spectrum to be sought

402:     Possible values:
403:     The parameter 'which' can have one of these values

405: +     PEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
406: .     PEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
407: .     PEP_LARGEST_REAL - largest real parts
408: .     PEP_SMALLEST_REAL - smallest real parts
409: .     PEP_LARGEST_IMAGINARY - largest imaginary parts
410: .     PEP_SMALLEST_IMAGINARY - smallest imaginary parts
411: .     PEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
412: .     PEP_TARGET_REAL - eigenvalues with real part closest to target
413: -     PEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target

415:     Options Database Keys:
416: +   -pep_largest_magnitude - Sets largest eigenvalues in magnitude
417: .   -pep_smallest_magnitude - Sets smallest eigenvalues in magnitude
418: .   -pep_largest_real - Sets largest real parts
419: .   -pep_smallest_real - Sets smallest real parts
420: .   -pep_largest_imaginary - Sets largest imaginary parts
421: .   -pep_smallest_imaginary - Sets smallest imaginary parts
422: .   -pep_target_magnitude - Sets eigenvalues closest to target
423: .   -pep_target_real - Sets real parts closest to target
424: -   -pep_target_imaginary - Sets imaginary parts closest to target

426:     Notes:
427:     Not all eigensolvers implemented in PEP account for all the possible values
428:     stated above. If SLEPc is compiled for real numbers PEP_LARGEST_IMAGINARY
429:     and PEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
430:     for eigenvalue selection.

432:     Level: intermediate

434: .seealso: PEPGetWhichEigenpairs(), PEPWhich
435: @*/
436: PetscErrorCode PEPSetWhichEigenpairs(PEP pep,PEPWhich which)
437: {
441:   switch (which) {
442:     case PEP_LARGEST_MAGNITUDE:
443:     case PEP_SMALLEST_MAGNITUDE:
444:     case PEP_LARGEST_REAL:
445:     case PEP_SMALLEST_REAL:
446:     case PEP_LARGEST_IMAGINARY:
447:     case PEP_SMALLEST_IMAGINARY:
448:     case PEP_TARGET_MAGNITUDE:
449:     case PEP_TARGET_REAL:
450: #if defined(PETSC_USE_COMPLEX)
451:     case PEP_TARGET_IMAGINARY:
452: #endif
453:       if (pep->which != which) {
454:         pep->state = PEP_STATE_INITIAL;
455:         pep->which = which;
456:       }
457:       break;
458:     default:
459:       SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
460:   }
461:   return(0);
462: }

466: /*@C
467:     PEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
468:     sought.

470:     Not Collective

472:     Input Parameter:
473: .   pep - eigensolver context obtained from PEPCreate()

475:     Output Parameter:
476: .   which - the portion of the spectrum to be sought

478:     Notes:
479:     See PEPSetWhichEigenpairs() for possible values of 'which'.

481:     Level: intermediate

483: .seealso: PEPSetWhichEigenpairs(), PEPWhich
484: @*/
485: PetscErrorCode PEPGetWhichEigenpairs(PEP pep,PEPWhich *which)
486: {
490:   *which = pep->which;
491:   return(0);
492: }

496: /*@
497:    PEPSetProblemType - Specifies the type of the polynomial eigenvalue problem.

499:    Logically Collective on PEP

501:    Input Parameters:
502: +  pep  - the polynomial eigensolver context
503: -  type - a known type of polynomial eigenvalue problem

505:    Options Database Keys:
506: +  -pep_general - general problem with no particular structure
507: .  -pep_hermitian - problem whose coefficient matrices are Hermitian
508: -  -pep_gyroscopic - problem with Hamiltonian structure

510:    Notes:
511:    Allowed values for the problem type are: general (PEP_GENERAL), Hermitian
512:    (PEP_HERMITIAN), and gyroscopic (PEP_GYROSCOPIC).

514:    This function is used to instruct SLEPc to exploit certain structure in
515:    the polynomial eigenproblem. By default, no particular structure is assumed.

517:    If the problem matrices are Hermitian (symmetric in the real case) or
518:    Hermitian/skew-Hermitian then the solver can exploit this fact to perform
519:    less operations or provide better stability.

521:    Level: intermediate

523: .seealso: PEPSetOperators(), PEPSetType(), PEPGetProblemType(), PEPProblemType
524: @*/
525: PetscErrorCode PEPSetProblemType(PEP pep,PEPProblemType type)
526: {
530:   if (type!=PEP_GENERAL && type!=PEP_HERMITIAN && type!=PEP_GYROSCOPIC) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
531:   pep->problem_type = type;
532:   return(0);
533: }

537: /*@C
538:    PEPGetProblemType - Gets the problem type from the PEP object.

540:    Not Collective

542:    Input Parameter:
543: .  pep - the polynomial eigensolver context

545:    Output Parameter:
546: .  type - name of PEP problem type

548:    Level: intermediate

550: .seealso: PEPSetProblemType(), PEPProblemType
551: @*/
552: PetscErrorCode PEPGetProblemType(PEP pep,PEPProblemType *type)
553: {
557:   *type = pep->problem_type;
558:   return(0);
559: }

563: /*@
564:    PEPSetBasis - Specifies the type of polynomial basis used to describe the
565:    polynomial eigenvalue problem.

567:    Logically Collective on PEP

569:    Input Parameters:
570: +  pep   - the polynomial eigensolver context
571: -  basis - the type of polynomial basis

573:    Options Database Key:
574: .  -pep_basis <basis> - Select the basis type

576:    Notes:
577:    By default, the coefficient matrices passed via PEPSetOperators() are
578:    expressed in the monomial basis, i.e. 
579:    P(lambda) = A_0 + lambda*A_1 + lambda^2*A_2 + ... + lambda^d*A_d.
580:    Other polynomial bases may have better numerical behaviour, but the user
581:    must then pass the coefficient matrices accordingly.

583:    Level: intermediate

585: .seealso: PEPSetOperators(), PEPGetBasis(), PEPBasis
586: @*/
587: PetscErrorCode PEPSetBasis(PEP pep,PEPBasis basis)
588: {
592:   pep->basis = basis;
593:   return(0);
594: }

598: /*@C
599:    PEPGetBasis - Gets the type of polynomial basis from the PEP object.

601:    Not Collective

603:    Input Parameter:
604: .  pep - the polynomial eigensolver context

606:    Output Parameter:
607: .  basis - the polynomial basis

609:    Level: intermediate

611: .seealso: PEPSetBasis(), PEPBasis
612: @*/
613: PetscErrorCode PEPGetBasis(PEP pep,PEPBasis *basis)
614: {
618:   *basis = pep->basis;
619:   return(0);
620: }

624: /*@
625:    PEPSetTrackAll - Specifies if the solver must compute the residual of all
626:    approximate eigenpairs or not.

628:    Logically Collective on PEP

630:    Input Parameters:
631: +  pep      - the eigensolver context
632: -  trackall - whether compute all residuals or not

634:    Notes:
635:    If the user sets trackall=PETSC_TRUE then the solver explicitly computes
636:    the residual for each eigenpair approximation. Computing the residual is
637:    usually an expensive operation and solvers commonly compute the associated
638:    residual to the first unconverged eigenpair.
639:    The options '-pep_monitor_all' and '-pep_monitor_lg_all' automatically
640:    activate this option.

642:    Level: intermediate

644: .seealso: PEPGetTrackAll()
645: @*/
646: PetscErrorCode PEPSetTrackAll(PEP pep,PetscBool trackall)
647: {
651:   pep->trackall = trackall;
652:   return(0);
653: }

657: /*@
658:    PEPGetTrackAll - Returns the flag indicating whether all residual norms must
659:    be computed or not.

661:    Not Collective

663:    Input Parameter:
664: .  pep - the eigensolver context

666:    Output Parameter:
667: .  trackall - the returned flag

669:    Level: intermediate

671: .seealso: PEPSetTrackAll()
672: @*/
673: PetscErrorCode PEPGetTrackAll(PEP pep,PetscBool *trackall)
674: {
678:   *trackall = pep->trackall;
679:   return(0);
680: }

684: /*@C
685:    PEPSetConvergenceTestFunction - Sets a function to compute the error estimate
686:    used in the convergence test.

688:    Logically Collective on PEP

690:    Input Parameters:
691: +  pep     - eigensolver context obtained from PEPCreate()
692: .  func    - a pointer to the convergence test function
693: .  ctx     - [optional] context for private data for the convergence routine
694: -  destroy - [optional] destructor for the context (may be NULL;
695:              PETSC_NULL_FUNCTION in Fortran)

697:    Calling Sequence of func:
698: $   func(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)

700: +   pep    - eigensolver context obtained from PEPCreate()
701: .   eigr   - real part of the eigenvalue
702: .   eigi   - imaginary part of the eigenvalue
703: .   res    - residual norm associated to the eigenpair
704: .   errest - (output) computed error estimate
705: -   ctx    - optional context, as set by PEPSetConvergenceTest()

707:    Note:
708:    If the error estimate returned by the convergence test function is less than
709:    the tolerance, then the eigenvalue is accepted as converged.

711:    Level: advanced

713: .seealso: PEPSetConvergenceTest(), PEPSetTolerances()
714: @*/
715: PetscErrorCode PEPSetConvergenceTestFunction(PEP pep,PetscErrorCode (*func)(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
716: {

721:   if (pep->convergeddestroy) {
722:     (*pep->convergeddestroy)(pep->convergedctx);
723:   }
724:   pep->converged        = func;
725:   pep->convergeddestroy = destroy;
726:   pep->convergedctx     = ctx;
727:   if (func == PEPConvergedEigRelative) pep->conv = PEP_CONV_EIG;
728:   else if (func == PEPConvergedNormRelative) pep->conv = PEP_CONV_NORM;
729:   else if (func == PEPConvergedAbsolute) pep->conv = PEP_CONV_ABS;
730:   else pep->conv = PEP_CONV_USER;
731:   return(0);
732: }

736: /*@
737:    PEPSetConvergenceTest - Specifies how to compute the error estimate
738:    used in the convergence test.

740:    Logically Collective on PEP

742:    Input Parameters:
743: +  pep  - eigensolver context obtained from PEPCreate()
744: -  conv - the type of convergence test

746:    Options Database Keys:
747: +  -pep_conv_abs  - Sets the absolute convergence test
748: .  -pep_conv_eig  - Sets the convergence test relative to the eigenvalue
749: .  -pep_conv_norm - Sets the convergence test relative to the matrix norms
750: -  -pep_conv_user - Selects the user-defined convergence test

752:    Note:
753:    The parameter 'conv' can have one of these values
754: +     PEP_CONV_ABS  - absolute error ||r||
755: .     PEP_CONV_EIG  - error relative to the eigenvalue l, ||r||/|l|
756: .     PEP_CONV_NORM - error relative to the matrix norms
757: -     PEP_CONV_USER - function set by PEPSetConvergenceTestFunction()

759:    Level: intermediate

761: .seealso: PEPGetConvergenceTest(), PEPSetConvergenceTestFunction(), PEPConv
762: @*/
763: PetscErrorCode PEPSetConvergenceTest(PEP pep,PEPConv conv)
764: {
768:   switch (conv) {
769:     case PEP_CONV_ABS:  pep->converged = PEPConvergedAbsolute; break;
770:     case PEP_CONV_EIG:  pep->converged = PEPConvergedEigRelative; break;
771:     case PEP_CONV_NORM: pep->converged = PEPConvergedNormRelative; break;
772:     case PEP_CONV_USER: break;
773:     default:
774:       SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
775:   }
776:   pep->conv = conv;
777:   return(0);
778: }

782: /*@
783:    PEPGetConvergenceTest - Gets the method used to compute the error estimate
784:    used in the convergence test.

786:    Not Collective

788:    Input Parameters:
789: .  pep   - eigensolver context obtained from PEPCreate()

791:    Output Parameters:
792: .  conv  - the type of convergence test

794:    Level: intermediate

796: .seealso: PEPSetConvergenceTest(), PEPConv
797: @*/
798: PetscErrorCode PEPGetConvergenceTest(PEP pep,PEPConv *conv)
799: {
803:   *conv = pep->conv;
804:   return(0);
805: }

809: /*@
810:    PEPSetScale - Specifies the scaling strategy to be used.

812:    Logically Collective on PEP

814:    Input Parameters:
815: +  pep    - the eigensolver context
816: .  scale  - scaling strategy
817: .  alpha  - the scaling factor used in the scalar strategy
818: .  its    - number of iterations of the diagonal scaling algorithm
819: -  lambda - approximation to wanted eigenvalues (modulus)

821:    Options Database Keys:
822: +  -pep_scale <type> - scaling type, one of <none,scalar,diagonal,both>
823: .  -pep_scale_factor <alpha> - the scaling factor
824: .  -pep_scale_its <its> - number of iterations
825: -  -pep_scale_lambda <lambda> - approximation to eigenvalues

827:    Notes:
828:    There are two non-exclusive scaling strategies: scalar and diagonal.

830:    In the scalar strategy, scaling is applied to the eigenvalue, that is,
831:    mu = lambda/alpha is the new eigenvalue and all matrices are scaled
832:    accordingly. After solving the scaled problem, the original lambda is
833:    recovered. Parameter 'alpha' must be positive. Use PETSC_DECIDE to let
834:    the solver compute a reasonable scaling factor.

836:    In the diagonal strategy, the solver works implicitly with matrix Dr*A*Dl,
837:    where Dr and Dl are appropriate diagonal matrices. This improves the accuracy
838:    of the computed results in some cases. This option requires MATAIJ matrices.
839:    The parameter 'its' is the number of iterations performed by the method.
840:    Parameter 'lambda' must be positive. Use PETSC_DECIDE or set lambda = 1.0 if no
841:    information about eigenvalues is available.

843:    Level: intermediate

845: .seealso: PEPGetScale()
846: @*/
847: PetscErrorCode PEPSetScale(PEP pep,PEPScale scale,PetscReal alpha,PetscInt its,PetscReal lambda)
848: {
852:   pep->scale = scale;
853:   if (scale==PEP_SCALE_SCALAR || scale==PEP_SCALE_BOTH) {
855:     if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
856:       pep->sfactor = 0.0;
857:       pep->sfactor_set = PETSC_FALSE;
858:     } else {
859:       if (alpha<=0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
860:       pep->sfactor = alpha;
861:       pep->sfactor_set = PETSC_TRUE;
862:     }
863:   }
864:   if (scale==PEP_SCALE_DIAGONAL || scale==PEP_SCALE_BOTH) {
867:     if (its==PETSC_DECIDE || its==PETSC_DEFAULT) pep->sits = 5;
868:     else pep->sits = its;
869:     if (lambda==PETSC_DECIDE || lambda==PETSC_DEFAULT) pep->slambda = 1.0;
870:     else if (lambda<=0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of lambda. Must be > 0");
871:     else pep->slambda = lambda;
872:   }
873:   pep->state = PEP_STATE_INITIAL;
874:   return(0);
875: }

879: /*@
880:    PEPGetScale - Gets the scaling strategy used by the PEP object, and the
881:    associated parameters.

883:    Not Collective

885:    Input Parameter:
886: .  pep - the eigensolver context

888:    Output Parameters:
889: +  scale  - scaling strategy
890: .  alpha  - the scaling factor used in the scalar strategy
891: .  its    - number of iterations of the diagonal scaling algorithm
892: -  lambda - approximation to wanted eigenvalues (modulus)

894:    Level: intermediate

896:    Note:
897:    The user can specify NULL for any parameter that is not needed.

899: .seealso: PEPSetScale()
900: @*/
901: PetscErrorCode PEPGetScale(PEP pep,PEPScale *scale,PetscReal *alpha,PetscInt *its,PetscReal *lambda)
902: {
905:   if (scale)  *scale  = pep->scale;
906:   if (alpha)  *alpha  = pep->sfactor;
907:   if (its)    *its    = pep->sits;
908:   if (lambda) *lambda = pep->slambda;
909:   return(0);
910: }

914: /*@
915:    PEPSetExtract - Specifies the extraction strategy to be used.

917:    Logically Collective on PEP

919:    Input Parameters:
920: +  pep     - the eigensolver context
921: -  extract - extraction strategy

923:    Options Database Keys:
924: .  -pep_extract <type> - extraction type, one of <norm,residual,structured>

926:    Level: intermediate

928: .seealso: PEPGetExtract()
929: @*/
930: PetscErrorCode PEPSetExtract(PEP pep,PEPExtract extract)
931: {
935:   pep->extract = extract;
936:   return(0);
937: }

941: /*@
942:    PEPGetExtract - Gets the extraction strategy used by the PEP object.

944:    Not Collective

946:    Input Parameter:
947: .  pep - the eigensolver context

949:    Output Parameter:
950: .  extract - extraction strategy

952:    Level: intermediate

954: .seealso: PEPSetExtract()
955: @*/
956: PetscErrorCode PEPGetExtract(PEP pep,PEPExtract *extract)
957: {
960:   if (extract) *extract = pep->extract;
961:   return(0);
962: }

966: /*@
967:    PEPSetRefine - Specifies the refinement type (and options) to be used
968:    after the solve.

970:    Logically Collective on PEP

972:    Input Parameters:
973: +  pep    - the polynomial eigensolver context
974: .  refine - refinement type
975: .  npart  - number of partitions of the communicator
976: .  tol    - the convergence tolerance
977: .  its    - maximum number of refinement iterations
978: -  schur  - boolean flag to activate the Schur complement approach

980:    Options Database Keys:
981: +  -pep_refine <type> - refinement type, one of <none,simple,multiple>
982: .  -pep_refine_partitions <n> - the number of partitions
983: .  -pep_refine_tol <tol> - the tolerance
984: .  -pep_refine_its <its> - number of iterations
985: -  -pep_refine_schur - to set the Schur complement approach

987:    Notes:
988:    By default, iterative refinement is disabled, since it may be very
989:    costly. There are two possible refinement strategies: simple and multiple.
990:    The simple approach performs iterative refinement on each of the
991:    converged eigenpairs individually, whereas the multiple strategy works
992:    with the invariant pair as a whole, refining all eigenpairs simultaneously.
993:    The latter may be required for the case of multiple eigenvalues.

995:    In some cases, especially when using direct solvers within the
996:    iterative refinement method, it may be helpful for improved scalability
997:    to split the communicator in several partitions. The npart parameter
998:    indicates how many partitions to use (defaults to 1).

1000:    The tol and its parameters specify the stopping criterion. In the simple
1001:    method, refinement continues until the residual of each eigenpair is
1002:    below the tolerance (tol defaults to the PEP tol, but may be set to a
1003:    different value). In contrast, the multiple method simply performs its
1004:    refinement iterations (just one by default).

1006:    The schur flag is used to change the way in which linear systems are
1007:    solved, so that a Schur complement approach is used instead of explicitly
1008:    building the coefficient matrix.

1010:    Level: intermediate

1012: .seealso: PEPGetRefine()
1013: @*/
1014: PetscErrorCode PEPSetRefine(PEP pep,PEPRefine refine,PetscInt npart,PetscReal tol,PetscInt its,PetscBool schur)
1015: {
1017:   PetscMPIInt    size;

1026:   pep->refine = refine;
1027:   if (refine) {  /* process parameters only if not REFINE_NONE */
1028:     if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
1029:       pep->npart = 1;
1030:     } else {
1031:       MPI_Comm_size(PetscObjectComm((PetscObject)pep),&size);
1032:       if (npart<1 || npart>size) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
1033:       pep->npart = npart;
1034:     }
1035:     if (tol == PETSC_DEFAULT || tol == PETSC_DECIDE) {
1036:       pep->rtol = pep->tol;
1037:     } else {
1038:       if (tol<=0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
1039:       pep->rtol = tol;
1040:     }
1041:     if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
1042:       pep->rits = PETSC_DEFAULT;
1043:     } else {
1044:       if (its<=0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
1045:       pep->rits = its;
1046:     }
1047:     pep->schur = schur;
1048:   }
1049:   pep->state = PEP_STATE_INITIAL;
1050:   return(0);
1051: }

1055: /*@
1056:    PEPGetRefine - Gets the refinement strategy used by the PEP object, and the
1057:    associated parameters.

1059:    Not Collective

1061:    Input Parameter:
1062: .  pep - the polynomial eigensolver context

1064:    Output Parameters:
1065: +  refine - refinement type
1066: .  npart  - number of partitions of the communicator
1067: .  tol    - the convergence tolerance
1068: .  its    - maximum number of refinement iterations
1069: -  schur  - whether the Schur complement approach is being used

1071:    Level: intermediate

1073:    Note:
1074:    The user can specify NULL for any parameter that is not needed.

1076: .seealso: PEPSetRefine()
1077: @*/
1078: PetscErrorCode PEPGetRefine(PEP pep,PEPRefine *refine,PetscInt *npart,PetscReal *tol,PetscInt *its,PetscBool *schur)
1079: {
1082:   if (refine) *refine = pep->refine;
1083:   if (npart)  *npart  = pep->npart;
1084:   if (tol)    *tol    = pep->rtol;
1085:   if (its)    *its    = pep->rits;
1086:   if (schur)  *schur  = pep->schur;
1087:   return(0);
1088: }

1092: /*@C
1093:    PEPSetOptionsPrefix - Sets the prefix used for searching for all
1094:    PEP options in the database.

1096:    Logically Collective on PEP

1098:    Input Parameters:
1099: +  pep - the polynomial eigensolver context
1100: -  prefix - the prefix string to prepend to all PEP option requests

1102:    Notes:
1103:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1104:    The first character of all runtime options is AUTOMATICALLY the
1105:    hyphen.

1107:    For example, to distinguish between the runtime options for two
1108:    different PEP contexts, one could call
1109: .vb
1110:       PEPSetOptionsPrefix(pep1,"qeig1_")
1111:       PEPSetOptionsPrefix(pep2,"qeig2_")
1112: .ve

1114:    Level: advanced

1116: .seealso: PEPAppendOptionsPrefix(), PEPGetOptionsPrefix()
1117: @*/
1118: PetscErrorCode PEPSetOptionsPrefix(PEP pep,const char *prefix)
1119: {

1124:   if (!pep->st) { PEPGetST(pep,&pep->st); }
1125:   STSetOptionsPrefix(pep->st,prefix);
1126:   if (!pep->V) { PEPGetBV(pep,&pep->V); }
1127:   BVSetOptionsPrefix(pep->V,prefix);
1128:   if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
1129:   DSSetOptionsPrefix(pep->ds,prefix);
1130:   if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
1131:   RGSetOptionsPrefix(pep->rg,prefix);
1132:   PetscObjectSetOptionsPrefix((PetscObject)pep,prefix);
1133:   return(0);
1134: }

1138: /*@C
1139:    PEPAppendOptionsPrefix - Appends to the prefix used for searching for all
1140:    PEP options in the database.

1142:    Logically Collective on PEP

1144:    Input Parameters:
1145: +  pep - the polynomial eigensolver context
1146: -  prefix - the prefix string to prepend to all PEP option requests

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

1152:    Level: advanced

1154: .seealso: PEPSetOptionsPrefix(), PEPGetOptionsPrefix()
1155: @*/
1156: PetscErrorCode PEPAppendOptionsPrefix(PEP pep,const char *prefix)
1157: {
1159:   PetscBool      flg;
1160:   EPS            eps;

1164:   if (!pep->st) { PEPGetST(pep,&pep->st); }
1165:   STAppendOptionsPrefix(pep->st,prefix);
1166:   if (!pep->V) { PEPGetBV(pep,&pep->V); }
1167:   BVSetOptionsPrefix(pep->V,prefix);
1168:   if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
1169:   DSSetOptionsPrefix(pep->ds,prefix);
1170:   if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
1171:   RGSetOptionsPrefix(pep->rg,prefix);
1172:   PetscObjectAppendOptionsPrefix((PetscObject)pep,prefix);
1173:   PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&flg);
1174:   if (flg) {
1175:     PEPLinearGetEPS(pep,&eps);
1176:     EPSSetOptionsPrefix(eps,((PetscObject)pep)->prefix);
1177:     EPSAppendOptionsPrefix(eps,"pep_");
1178:   }
1179:   return(0);
1180: }

1184: /*@C
1185:    PEPGetOptionsPrefix - Gets the prefix used for searching for all
1186:    PEP options in the database.

1188:    Not Collective

1190:    Input Parameters:
1191: .  pep - the polynomial eigensolver context

1193:    Output Parameters:
1194: .  prefix - pointer to the prefix string used is returned

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

1199:    Level: advanced

1201: .seealso: PEPSetOptionsPrefix(), PEPAppendOptionsPrefix()
1202: @*/
1203: PetscErrorCode PEPGetOptionsPrefix(PEP pep,const char *prefix[])
1204: {

1210:   PetscObjectGetOptionsPrefix((PetscObject)pep,prefix);
1211:   return(0);
1212: }