Actual source code: qepopts.c

  1: /*
  2:       QEP 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-2012, Universitat Politecnica de Valencia, Spain

  9:    This file is part of SLEPc.
 10:       
 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/qepimpl.h>       /*I "slepcqep.h" I*/

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

 34:    Collective on QEP

 36:    Input Parameters:
 37: .  qep - the quadratic eigensolver context

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

 42:    Level: beginner
 43: @*/
 44: PetscErrorCode QEPSetFromOptions(QEP qep)
 45: {
 46:   PetscErrorCode   ierr;
 47:   char             type[256],monfilename[PETSC_MAX_PATH_LEN];
 48:   PetscBool        flg,val;
 49:   PetscReal        r;
 50:   PetscInt         i,j,k;
 51:   PetscViewer      monviewer;
 52:   SlepcConvMonitor ctx;

 56:   if (!QEPRegisterAllCalled) { QEPRegisterAll(PETSC_NULL); }
 57:   if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
 58:   PetscOptionsBegin(((PetscObject)qep)->comm,((PetscObject)qep)->prefix,"Quadratic Eigenvalue Problem (QEP) Solver Options","QEP");
 59:     PetscOptionsList("-qep_type","Quadratic Eigenvalue Problem method","QEPSetType",QEPList,(char*)(((PetscObject)qep)->type_name?((PetscObject)qep)->type_name:QEPLINEAR),type,256,&flg);
 60:     if (flg) {
 61:       QEPSetType(qep,type);
 62:     } else if (!((PetscObject)qep)->type_name) {
 63:       QEPSetType(qep,QEPLINEAR);
 64:     }

 66:     PetscOptionsBoolGroupBegin("-qep_general","general quadratic eigenvalue problem","QEPSetProblemType",&flg);
 67:     if (flg) {QEPSetProblemType(qep,QEP_GENERAL);}
 68:     PetscOptionsBoolGroup("-qep_hermitian","hermitian quadratic eigenvalue problem","QEPSetProblemType",&flg);
 69:     if (flg) {QEPSetProblemType(qep,QEP_HERMITIAN);}
 70:     PetscOptionsBoolGroupEnd("-qep_gyroscopic","gyroscopic quadratic eigenvalue problem","QEPSetProblemType",&flg);
 71:     if (flg) {QEPSetProblemType(qep,QEP_GYROSCOPIC);}

 73:     r = PETSC_IGNORE;
 74:     PetscOptionsReal("-qep_scale","Scale factor","QEPSetScaleFactor",qep->sfactor,&r,PETSC_NULL);
 75:     QEPSetScaleFactor(qep,r);

 77:     r = i = PETSC_IGNORE;
 78:     PetscOptionsInt("-qep_max_it","Maximum number of iterations","QEPSetTolerances",qep->max_it,&i,PETSC_NULL);
 79:     PetscOptionsReal("-qep_tol","Tolerance","QEPSetTolerances",qep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:qep->tol,&r,PETSC_NULL);
 80:     QEPSetTolerances(qep,r,i);
 81:     PetscOptionsBoolGroupBegin("-qep_convergence_default","Default (relative error) convergence test","QEPSetConvergenceTest",&flg);
 82:     if (flg) {QEPSetConvergenceTest(qep,QEPConvergedDefault,PETSC_NULL);}
 83:     PetscOptionsBoolGroupEnd("-qep_convergence_absolute","Absolute error convergence test","QEPSetConvergenceTest",&flg);
 84:     if (flg) {QEPSetConvergenceTest(qep,QEPConvergedAbsolute,PETSC_NULL);}

 86:     i = j = k = PETSC_IGNORE;
 87:     PetscOptionsInt("-qep_nev","Number of eigenvalues to compute","QEPSetDimensions",qep->nev,&i,PETSC_NULL);
 88:     PetscOptionsInt("-qep_ncv","Number of basis vectors","QEPSetDimensions",qep->ncv,&j,PETSC_NULL);
 89:     PetscOptionsInt("-qep_mpd","Maximum dimension of projected problem","QEPSetDimensions",qep->mpd,&k,PETSC_NULL);
 90:     QEPSetDimensions(qep,i,j,k);
 91: 
 92:     /* -----------------------------------------------------------------------*/
 93:     /*
 94:       Cancels all monitors hardwired into code before call to QEPSetFromOptions()
 95:     */
 96:     flg  = PETSC_FALSE;
 97:     PetscOptionsBool("-qep_monitor_cancel","Remove any hardwired monitor routines","QEPMonitorCancel",flg,&flg,PETSC_NULL);
 98:     if (flg) {
 99:       QEPMonitorCancel(qep);
100:     }
101:     /*
102:       Prints approximate eigenvalues and error estimates at each iteration
103:     */
104:     PetscOptionsString("-qep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
105:     if (flg) {
106:       PetscViewerASCIIOpen(((PetscObject)qep)->comm,monfilename,&monviewer);
107:       QEPMonitorSet(qep,QEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
108:     }
109:     PetscOptionsString("-qep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
110:     if (flg) {
111:       PetscNew(struct _n_SlepcConvMonitor,&ctx);
112:       PetscViewerASCIIOpen(((PetscObject)qep)->comm,monfilename,&ctx->viewer);
113:       QEPMonitorSet(qep,QEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
114:     }
115:     PetscOptionsString("-qep_monitor_all","Monitor approximate eigenvalues and error estimates","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
116:     if (flg) {
117:       PetscViewerASCIIOpen(((PetscObject)qep)->comm,monfilename,&monviewer);
118:       QEPMonitorSet(qep,QEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
119:       QEPSetTrackAll(qep,PETSC_TRUE);
120:     }
121:     flg = PETSC_FALSE;
122:     PetscOptionsBool("-qep_monitor_draw","Monitor first unconverged approximate error estimate graphically","QEPMonitorSet",flg,&flg,PETSC_NULL);
123:     if (flg) {
124:       QEPMonitorSet(qep,QEPMonitorLG,PETSC_NULL,PETSC_NULL);
125:     }
126:     flg = PETSC_FALSE;
127:     PetscOptionsBool("-qep_monitor_draw_all","Monitor error estimates graphically","QEPMonitorSet",flg,&flg,PETSC_NULL);
128:     if (flg) {
129:       QEPMonitorSet(qep,QEPMonitorLGAll,PETSC_NULL,PETSC_NULL);
130:       QEPSetTrackAll(qep,PETSC_TRUE);
131:     }
132:   /* -----------------------------------------------------------------------*/

134:     PetscOptionsBoolGroupBegin("-qep_largest_magnitude","compute largest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
135:     if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_MAGNITUDE);}
136:     PetscOptionsBoolGroup("-qep_smallest_magnitude","compute smallest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
137:     if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_MAGNITUDE);}
138:     PetscOptionsBoolGroup("-qep_largest_real","compute largest real parts","QEPSetWhichEigenpairs",&flg);
139:     if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_REAL);}
140:     PetscOptionsBoolGroup("-qep_smallest_real","compute smallest real parts","QEPSetWhichEigenpairs",&flg);
141:     if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_REAL);}
142:     PetscOptionsBoolGroup("-qep_largest_imaginary","compute largest imaginary parts","QEPSetWhichEigenpairs",&flg);
143:     if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_IMAGINARY);}
144:     PetscOptionsBoolGroupEnd("-qep_smallest_imaginary","compute smallest imaginary parts","QEPSetWhichEigenpairs",&flg);
145:     if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_IMAGINARY);}

147:     PetscOptionsBool("-qep_left_vectors","Compute left eigenvectors also","QEPSetLeftVectorsWanted",qep->leftvecs,&val,&flg);
148:     if (flg) {
149:       QEPSetLeftVectorsWanted(qep,val);
150:     }

152:     PetscOptionsName("-qep_view","Print detailed information on solver used","QEPView",0);
153:     PetscOptionsName("-qep_view_binary","Save the matrices associated to the eigenproblem","QEPSetFromOptions",0);
154:     PetscOptionsName("-qep_plot_eigs","Make a plot of the computed eigenvalues","QEPSolve",0);
155: 
156:     if (qep->ops->setfromoptions) {
157:       (*qep->ops->setfromoptions)(qep);
158:     }
159:     PetscObjectProcessOptionsHandlers((PetscObject)qep);
160:   PetscOptionsEnd();

162:   if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
163:   IPSetFromOptions(qep->ip);
164:   if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
165:   DSSetFromOptions(qep->ds);
166:   PetscRandomSetFromOptions(qep->rand);
167:   return(0);
168: }

172: /*@
173:    QEPGetTolerances - Gets the tolerance and maximum iteration count used
174:    by the QEP convergence tests. 

176:    Not Collective

178:    Input Parameter:
179: .  qep - the quadratic eigensolver context
180:   
181:    Output Parameters:
182: +  tol - the convergence tolerance
183: -  maxits - maximum number of iterations

185:    Notes:
186:    The user can specify PETSC_NULL for any parameter that is not needed.

188:    Level: intermediate

190: .seealso: QEPSetTolerances()
191: @*/
192: PetscErrorCode QEPGetTolerances(QEP qep,PetscReal *tol,PetscInt *maxits)
193: {
196:   if (tol)    *tol    = qep->tol;
197:   if (maxits) *maxits = qep->max_it;
198:   return(0);
199: }

203: /*@
204:    QEPSetTolerances - Sets the tolerance and maximum iteration count used
205:    by the QEP convergence tests. 

207:    Logically Collective on QEP

209:    Input Parameters:
210: +  qep - the quadratic eigensolver context
211: .  tol - the convergence tolerance
212: -  maxits - maximum number of iterations to use

214:    Options Database Keys:
215: +  -qep_tol <tol> - Sets the convergence tolerance
216: -  -qep_max_it <maxits> - Sets the maximum number of iterations allowed

218:    Notes:
219:    Use PETSC_IGNORE for an argument that need not be changed.

221:    Use PETSC_DECIDE for maxits to assign a reasonably good value, which is 
222:    dependent on the solution method.

224:    Level: intermediate

226: .seealso: QEPGetTolerances()
227: @*/
228: PetscErrorCode QEPSetTolerances(QEP qep,PetscReal tol,PetscInt maxits)
229: {
234:   if (tol != PETSC_IGNORE) {
235:     if (tol == PETSC_DEFAULT) {
236:       qep->tol = PETSC_DEFAULT;
237:     } else {
238:       if (tol < 0.0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
239:       qep->tol = tol;
240:     }
241:   }
242:   if (maxits != PETSC_IGNORE) {
243:     if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
244:       qep->max_it = 0;
245:       qep->setupcalled = 0;
246:     } else {
247:       if (maxits < 0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
248:       qep->max_it = maxits;
249:     }
250:   }
251:   return(0);
252: }

256: /*@
257:    QEPGetDimensions - Gets the number of eigenvalues to compute
258:    and the dimension of the subspace.

260:    Not Collective

262:    Input Parameter:
263: .  qep - the quadratic eigensolver context
264:   
265:    Output Parameters:
266: +  nev - number of eigenvalues to compute
267: .  ncv - the maximum dimension of the subspace to be used by the solver
268: -  mpd - the maximum dimension allowed for the projected problem

270:    Notes:
271:    The user can specify PETSC_NULL for any parameter that is not needed.

273:    Level: intermediate

275: .seealso: QEPSetDimensions()
276: @*/
277: PetscErrorCode QEPGetDimensions(QEP qep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
278: {
281:   if (nev) *nev = qep->nev;
282:   if (ncv) *ncv = qep->ncv;
283:   if (mpd) *mpd = qep->mpd;
284:   return(0);
285: }

289: /*@
290:    QEPSetDimensions - Sets the number of eigenvalues to compute
291:    and the dimension of the subspace.

293:    Logically Collective on QEP

295:    Input Parameters:
296: +  qep - the quadratic eigensolver context
297: .  nev - number of eigenvalues to compute
298: .  ncv - the maximum dimension of the subspace to be used by the solver
299: -  mpd - the maximum dimension allowed for the projected problem

301:    Options Database Keys:
302: +  -qep_nev <nev> - Sets the number of eigenvalues
303: .  -qep_ncv <ncv> - Sets the dimension of the subspace
304: -  -qep_mpd <mpd> - Sets the maximum projected dimension

306:    Notes:
307:    Use PETSC_IGNORE to retain the previous value of any parameter.

309:    Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
310:    dependent on the solution method.

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

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

321:    Level: intermediate

323: .seealso: QEPGetDimensions()
324: @*/
325: PetscErrorCode QEPSetDimensions(QEP qep,PetscInt nev,PetscInt ncv,PetscInt mpd)
326: {
332:   if (nev != PETSC_IGNORE) {
333:     if (nev<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
334:     qep->nev = nev;
335:     qep->setupcalled = 0;
336:   }
337:   if (ncv != PETSC_IGNORE) {
338:     if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
339:       qep->ncv = 0;
340:     } else {
341:       if (ncv<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
342:       qep->ncv = ncv;
343:     }
344:     qep->setupcalled = 0;
345:   }
346:   if (mpd != PETSC_IGNORE) {
347:     if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
348:       qep->mpd = 0;
349:     } else {
350:       if (mpd<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
351:       qep->mpd = mpd;
352:     }
353:   }
354:   return(0);
355: }

359: /*@
360:     QEPSetWhichEigenpairs - Specifies which portion of the spectrum is 
361:     to be sought.

363:     Logically Collective on QEP

365:     Input Parameters:
366: +   qep   - eigensolver context obtained from QEPCreate()
367: -   which - the portion of the spectrum to be sought

369:     Possible values:
370:     The parameter 'which' can have one of these values
371:     
372: +     QEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
373: .     QEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
374: .     QEP_LARGEST_REAL - largest real parts
375: .     QEP_SMALLEST_REAL - smallest real parts
376: .     QEP_LARGEST_IMAGINARY - largest imaginary parts
377: -     QEP_SMALLEST_IMAGINARY - smallest imaginary parts

379:     Options Database Keys:
380: +   -qep_largest_magnitude - Sets largest eigenvalues in magnitude
381: .   -qep_smallest_magnitude - Sets smallest eigenvalues in magnitude
382: .   -qep_largest_real - Sets largest real parts
383: .   -qep_smallest_real - Sets smallest real parts
384: .   -qep_largest_imaginary - Sets largest imaginary parts
385: -   -qep_smallest_imaginary - Sets smallest imaginary parts

387:     Notes:
388:     Not all eigensolvers implemented in QEP account for all the possible values
389:     stated above. If SLEPc is compiled for real numbers QEP_LARGEST_IMAGINARY
390:     and QEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part 
391:     for eigenvalue selection.
392:     
393:     Level: intermediate

395: .seealso: QEPGetWhichEigenpairs(), QEPWhich
396: @*/
397: PetscErrorCode QEPSetWhichEigenpairs(QEP qep,QEPWhich which)
398: {
402:   if (which!=PETSC_IGNORE) {
403:     if (which==PETSC_DECIDE || which==PETSC_DEFAULT) qep->which = (QEPWhich)0;
404:     else switch (which) {
405:       case QEP_LARGEST_MAGNITUDE:
406:       case QEP_SMALLEST_MAGNITUDE:
407:       case QEP_LARGEST_REAL:
408:       case QEP_SMALLEST_REAL:
409:       case QEP_LARGEST_IMAGINARY:
410:       case QEP_SMALLEST_IMAGINARY:
411:         if (qep->which != which) {
412:           qep->setupcalled = 0;
413:           qep->which = which;
414:         }
415:         break;
416:       default:
417:         SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
418:     }
419:   }
420:   return(0);
421: }

425: /*@C
426:     QEPGetWhichEigenpairs - Returns which portion of the spectrum is to be 
427:     sought.

429:     Not Collective

431:     Input Parameter:
432: .   qep - eigensolver context obtained from QEPCreate()

434:     Output Parameter:
435: .   which - the portion of the spectrum to be sought

437:     Notes:
438:     See QEPSetWhichEigenpairs() for possible values of 'which'.

440:     Level: intermediate

442: .seealso: QEPSetWhichEigenpairs(), QEPWhich
443: @*/
444: PetscErrorCode QEPGetWhichEigenpairs(QEP qep,QEPWhich *which)
445: {
449:   *which = qep->which;
450:   return(0);
451: }

455: /*@
456:     QEPSetLeftVectorsWanted - Specifies which eigenvectors are required.

458:     Logically Collective on QEP

460:     Input Parameters:
461: +   qep      - the quadratic eigensolver context
462: -   leftvecs - whether left eigenvectors are required or not

464:     Options Database Keys:
465: .   -qep_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'

467:     Notes:
468:     If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
469:     the algorithm that computes both right and left eigenvectors. This is
470:     usually much more costly. This option is not available in all solvers.

472:     Level: intermediate

474: .seealso: QEPGetLeftVectorsWanted()
475: @*/
476: PetscErrorCode QEPSetLeftVectorsWanted(QEP qep,PetscBool leftvecs)
477: {
481:   if (qep->leftvecs != leftvecs) {
482:     qep->leftvecs = leftvecs;
483:     qep->setupcalled = 0;
484:   }
485:   return(0);
486: }

490: /*@C
491:     QEPGetLeftVectorsWanted - Returns the flag indicating whether left 
492:     eigenvectors are required or not.

494:     Not Collective

496:     Input Parameter:
497: .   qep - the eigensolver context

499:     Output Parameter:
500: .   leftvecs - the returned flag

502:     Level: intermediate

504: .seealso: QEPSetLeftVectorsWanted()
505: @*/
506: PetscErrorCode QEPGetLeftVectorsWanted(QEP qep,PetscBool *leftvecs)
507: {
511:   *leftvecs = qep->leftvecs;
512:   return(0);
513: }

517: /*@
518:    QEPGetScaleFactor - Gets the factor used for scaling the quadratic eigenproblem.

520:    Not Collective

522:    Input Parameter:
523: .  qep - the quadratic eigensolver context
524:   
525:    Output Parameters:
526: .  alpha - the scaling factor

528:    Notes:
529:    If the user did not specify a scaling factor, then after QEPSolve() the
530:    default value is returned.

532:    Level: intermediate

534: .seealso: QEPSetScaleFactor(), QEPSolve()
535: @*/
536: PetscErrorCode QEPGetScaleFactor(QEP qep,PetscReal *alpha)
537: {
541:   *alpha = qep->sfactor;
542:   return(0);
543: }

547: /*@
548:    QEPSetScaleFactor - Sets the scaling factor to be used for scaling the
549:    quadratic problem before attempting to solve.

551:    Logically Collective on QEP

553:    Input Parameters:
554: +  qep   - the quadratic eigensolver context
555: -  alpha - the scaling factor

557:    Options Database Keys:
558: .  -qep_scale <alpha> - Sets the scaling factor

560:    Notes:
561:    For the problem (l^2*M + l*C + K)*x = 0, the effect of scaling is to work
562:    with matrices (alpha^2*M, alpha*C, K), then scale the computed eigenvalue.

564:    The default is to scale with alpha = norm(K)/norm(M).

566:    Level: intermediate

568: .seealso: QEPGetScaleFactor()
569: @*/
570: PetscErrorCode QEPSetScaleFactor(QEP qep,PetscReal alpha)
571: {
575:   if (alpha != PETSC_IGNORE) {
576:     if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
577:       qep->sfactor = 0.0;
578:     } else {
579:       if (alpha < 0.0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
580:       qep->sfactor = alpha;
581:     }
582:   }
583:   return(0);
584: }

588: /*@
589:    QEPSetProblemType - Specifies the type of the quadratic eigenvalue problem.

591:    Logically Collective on QEP

593:    Input Parameters:
594: +  qep      - the quadratic eigensolver context
595: -  type     - a known type of quadratic eigenvalue problem 

597:    Options Database Keys:
598: +  -qep_general - general problem with no particular structure
599: .  -qep_hermitian - problem whose coefficient matrices are Hermitian
600: -  -qep_gyroscopic - problem with Hamiltonian structure
601:     
602:    Notes:  
603:    Allowed values for the problem type are: general (QEP_GENERAL), Hermitian
604:    (QEP_HERMITIAN), and gyroscopic (QEP_GYROSCOPIC).

606:    This function is used to instruct SLEPc to exploit certain structure in
607:    the quadratic eigenproblem. By default, no particular structure is assumed.

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

613:    Level: intermediate

615: .seealso: QEPSetOperators(), QEPSetType(), QEPGetProblemType(), QEPProblemType
616: @*/
617: PetscErrorCode QEPSetProblemType(QEP qep,QEPProblemType type)
618: {
622:   if (type!=QEP_GENERAL && type!=QEP_HERMITIAN && type!=QEP_GYROSCOPIC)
623:     SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
624:   qep->problem_type = type;
625:   return(0);
626: }

630: /*@C
631:    QEPGetProblemType - Gets the problem type from the QEP object.

633:    Not Collective

635:    Input Parameter:
636: .  qep - the quadratic eigensolver context 

638:    Output Parameter:
639: .  type - name of QEP problem type 

641:    Level: intermediate

643: .seealso: QEPSetProblemType(), QEPProblemType
644: @*/
645: PetscErrorCode QEPGetProblemType(QEP qep,QEPProblemType *type)
646: {
650:   *type = qep->problem_type;
651:   return(0);
652: }

656: /*@C
657:     QEPSetConvergenceTest - Sets a function to compute the error estimate used in 
658:     the convergence test.

660:     Logically Collective on QEP

662:     Input Parameters:
663: +   qep  - eigensolver context obtained from QEPCreate()
664: .   func - a pointer to the convergence test function
665: -   ctx  - a context pointer (the last parameter to the convergence test function)

667:     Calling Sequence of func:
668: $   func(QEP qep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal* errest,void *ctx)

670: +   qep    - eigensolver context obtained from QEPCreate()
671: .   eigr   - real part of the eigenvalue
672: .   eigi   - imaginary part of the eigenvalue
673: .   res    - residual norm associated to the eigenpair
674: .   errest - (output) computed error estimate
675: -   ctx    - optional context, as set by QEPSetConvergenceTest()

677:     Note:
678:     If the error estimate returned by the convergence test function is less than
679:     the tolerance, then the eigenvalue is accepted as converged.

681:     Level: advanced

683: .seealso: QEPSetTolerances()
684: @*/
685: extern PetscErrorCode QEPSetConvergenceTest(QEP qep,PetscErrorCode (*func)(QEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
686: {
689:   qep->conv_func = func;
690:   qep->conv_ctx = ctx;
691:   return(0);
692: }

696: /*@
697:    QEPSetTrackAll - Specifies if the solver must compute the residual of all
698:    approximate eigenpairs or not.

700:    Logically Collective on QEP

702:    Input Parameters:
703: +  qep      - the eigensolver context
704: -  trackall - whether compute all residuals or not

706:    Notes:
707:    If the user sets trackall=PETSC_TRUE then the solver explicitly computes
708:    the residual for each eigenpair approximation. Computing the residual is
709:    usually an expensive operation and solvers commonly compute the associated
710:    residual to the first unconverged eigenpair.
711:    The options '-qep_monitor_all' and '-qep_monitor_draw_all' automatically
712:    activates this option.

714:    Level: intermediate

716: .seealso: QEPGetTrackAll()
717: @*/
718: PetscErrorCode QEPSetTrackAll(QEP qep,PetscBool trackall)
719: {
723:   qep->trackall = trackall;
724:   return(0);
725: }

729: /*@
730:    QEPGetTrackAll - Returns the flag indicating whether all residual norms must
731:    be computed or not.

733:    Not Collective

735:    Input Parameter:
736: .  qep - the eigensolver context

738:    Output Parameter:
739: .  trackall - the returned flag

741:    Level: intermediate

743: .seealso: QEPSetTrackAll()
744: @*/
745: PetscErrorCode QEPGetTrackAll(QEP qep,PetscBool *trackall)
746: {
750:   *trackall = qep->trackall;
751:   return(0);
752: }

756: /*@C
757:    QEPSetOptionsPrefix - Sets the prefix used for searching for all 
758:    QEP options in the database.

760:    Logically Collective on QEP

762:    Input Parameters:
763: +  qep - the quadratic eigensolver context
764: -  prefix - the prefix string to prepend to all QEP option requests

766:    Notes:
767:    A hyphen (-) must NOT be given at the beginning of the prefix name.
768:    The first character of all runtime options is AUTOMATICALLY the
769:    hyphen.

771:    For example, to distinguish between the runtime options for two
772:    different QEP contexts, one could call
773: .vb
774:       QEPSetOptionsPrefix(qep1,"qeig1_")
775:       QEPSetOptionsPrefix(qep2,"qeig2_")
776: .ve

778:    Level: advanced

780: .seealso: QEPAppendOptionsPrefix(), QEPGetOptionsPrefix()
781: @*/
782: PetscErrorCode QEPSetOptionsPrefix(QEP qep,const char *prefix)
783: {

788:   if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
789:   IPSetOptionsPrefix(qep->ip,prefix);
790:   if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
791:   DSSetOptionsPrefix(qep->ds,prefix);
792:   PetscObjectSetOptionsPrefix((PetscObject)qep,prefix);
793:   return(0);
794: }
795: 
798: /*@C
799:    QEPAppendOptionsPrefix - Appends to the prefix used for searching for all 
800:    QEP options in the database.

802:    Logically Collective on QEP

804:    Input Parameters:
805: +  qep - the quadratic eigensolver context
806: -  prefix - the prefix string to prepend to all QEP option requests

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

812:    Level: advanced

814: .seealso: QEPSetOptionsPrefix(), QEPGetOptionsPrefix()
815: @*/
816: PetscErrorCode QEPAppendOptionsPrefix(QEP qep,const char *prefix)
817: {
819:   PetscBool      flg;
820:   EPS            eps;

824:   if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
825:   IPSetOptionsPrefix(qep->ip,prefix);
826:   if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
827:   DSSetOptionsPrefix(qep->ds,prefix);
828:   PetscObjectAppendOptionsPrefix((PetscObject)qep,prefix);
829:   PetscObjectTypeCompare((PetscObject)qep,QEPLINEAR,&flg);
830:   if (flg) {
831:     QEPLinearGetEPS(qep,&eps);
832:     EPSSetOptionsPrefix(eps,((PetscObject)qep)->prefix);
833:     EPSAppendOptionsPrefix(eps,"qep_");
834:   }
835:   return(0);
836: }

840: /*@C
841:    QEPGetOptionsPrefix - Gets the prefix used for searching for all 
842:    QEP options in the database.

844:    Not Collective

846:    Input Parameters:
847: .  qep - the quadratic eigensolver context

849:    Output Parameters:
850: .  prefix - pointer to the prefix string used is returned

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

855:    Level: advanced

857: .seealso: QEPSetOptionsPrefix(), QEPAppendOptionsPrefix()
858: @*/
859: PetscErrorCode QEPGetOptionsPrefix(QEP qep,const char *prefix[])
860: {

866:   PetscObjectGetOptionsPrefix((PetscObject)qep,prefix);
867:   return(0);
868: }