Actual source code: epsbasic.c

slepc-3.5.3 2014-12-19
Report Typos and Errors
  1: /*
  2:      The basic EPS routines, Create, View, etc. are here.

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

  8:    This file is part of SLEPc.

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

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

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

 24: #include <slepc-private/epsimpl.h>      /*I "slepceps.h" I*/

 26: PetscFunctionList EPSList = 0;
 27: PetscBool         EPSRegisterAllCalled = PETSC_FALSE;
 28: PetscClassId      EPS_CLASSID = 0;
 29: PetscLogEvent     EPS_SetUp = 0,EPS_Solve = 0;

 33: /*@C
 34:    EPSView - Prints the EPS data structure.

 36:    Collective on EPS

 38:    Input Parameters:
 39: +  eps - the eigenproblem solver context
 40: -  viewer - optional visualization context

 42:    Options Database Key:
 43: .  -eps_view -  Calls EPSView() at end of EPSSolve()

 45:    Note:
 46:    The available visualization contexts include
 47: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 48: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 49:          output where only the first processor opens
 50:          the file.  All other processors send their
 51:          data to the first processor to print.

 53:    The user can open an alternative visualization context with
 54:    PetscViewerASCIIOpen() - output to a specified file.

 56:    Level: beginner

 58: .seealso: STView(), PetscViewerASCIIOpen()
 59: @*/
 60: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
 61: {
 63:   const char     *type,*extr,*bal;
 64:   char           str[50];
 65:   PetscBool      isascii,ispower,isexternal,istrivial;

 69:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));

 73: #if defined(PETSC_USE_COMPLEX)
 74: #define HERM "hermitian"
 75: #else
 76: #define HERM "symmetric"
 77: #endif
 78:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 79:   if (isascii) {
 80:     PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
 81:     if (eps->ops->view) {
 82:       PetscViewerASCIIPushTab(viewer);
 83:       (*eps->ops->view)(eps,viewer);
 84:       PetscViewerASCIIPopTab(viewer);
 85:     }
 86:     if (eps->problem_type) {
 87:       switch (eps->problem_type) {
 88:         case EPS_HEP:   type = HERM " eigenvalue problem"; break;
 89:         case EPS_GHEP:  type = "generalized " HERM " eigenvalue problem"; break;
 90:         case EPS_NHEP:  type = "non-" HERM " eigenvalue problem"; break;
 91:         case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
 92:         case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
 93:         case EPS_GHIEP: type = "generalized " HERM "-indefinite eigenvalue problem"; break;
 94:         default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->problem_type");
 95:       }
 96:     } else type = "not yet set";
 97:     PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type);
 98:     if (eps->extraction) {
 99:       switch (eps->extraction) {
100:         case EPS_RITZ:             extr = "Rayleigh-Ritz"; break;
101:         case EPS_HARMONIC:         extr = "harmonic Ritz"; break;
102:         case EPS_HARMONIC_RELATIVE:extr = "relative harmonic Ritz"; break;
103:         case EPS_HARMONIC_RIGHT:   extr = "right harmonic Ritz"; break;
104:         case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
105:         case EPS_REFINED:          extr = "refined Ritz"; break;
106:         case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
107:         default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->extraction");
108:       }
109:       PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",extr);
110:     }
111:     if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
112:       switch (eps->balance) {
113:         case EPS_BALANCE_ONESIDE:   bal = "one-sided Krylov"; break;
114:         case EPS_BALANCE_TWOSIDE:   bal = "two-sided Krylov"; break;
115:         case EPS_BALANCE_USER:      bal = "user-defined matrix"; break;
116:         default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->balance");
117:       }
118:       PetscViewerASCIIPrintf(viewer,"  balancing enabled: %s",bal);
119:       if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
120:         PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
121:       }
122:       if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
123:         PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
124:       }
125:       PetscViewerASCIIPrintf(viewer,"\n");
126:     }
127:     PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: ");
128:     SlepcSNPrintfScalar(str,50,eps->target,PETSC_FALSE);
129:     if (!eps->which) {
130:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
131:     } else switch (eps->which) {
132:       case EPS_WHICH_USER:
133:         PetscViewerASCIIPrintf(viewer,"user defined\n");
134:         break;
135:       case EPS_TARGET_MAGNITUDE:
136:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
137:         break;
138:       case EPS_TARGET_REAL:
139:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
140:         break;
141:       case EPS_TARGET_IMAGINARY:
142:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
143:         break;
144:       case EPS_LARGEST_MAGNITUDE:
145:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
146:         break;
147:       case EPS_SMALLEST_MAGNITUDE:
148:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
149:         break;
150:       case EPS_LARGEST_REAL:
151:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
152:         break;
153:       case EPS_SMALLEST_REAL:
154:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
155:         break;
156:       case EPS_LARGEST_IMAGINARY:
157:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
158:         break;
159:       case EPS_SMALLEST_IMAGINARY:
160:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
161:         break;
162:       case EPS_ALL:
163:         PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
164:         break;
165:       default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->which");
166:     }
167:     if (eps->trueres) {
168:       PetscViewerASCIIPrintf(viewer,"  computing true residuals explicitly\n");
169:     }
170:     if (eps->trackall) {
171:       PetscViewerASCIIPrintf(viewer,"  computing all residuals (for tracking convergence)\n");
172:     }
173:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %D\n",eps->nev);
174:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",eps->ncv);
175:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %D\n",eps->mpd);
176:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",eps->max_it);
177:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)eps->tol);
178:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
179:     switch (eps->conv) {
180:     case EPS_CONV_ABS:
181:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
182:     case EPS_CONV_EIG:
183:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
184:     case EPS_CONV_NORM:
185:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");
186:       PetscViewerASCIIPrintf(viewer,"  computed matrix norms: norm(A)=%g",(double)eps->nrma);
187:       if (eps->isgeneralized) {
188:         PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb);
189:       }
190:       PetscViewerASCIIPrintf(viewer,"\n");
191:       break;
192:     case EPS_CONV_USER:
193:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
194:     }
195:     if (eps->nini) {
196:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
197:     }
198:     if (eps->nds) {
199:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
200:     }
201:   } else {
202:     if (eps->ops->view) {
203:       (*eps->ops->view)(eps,viewer);
204:     }
205:   }
206:   PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLZPACK,EPSTRLAN,EPSBLOPEX,EPSPRIMME,"");
207:   if (!isexternal) {
208:     PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
209:     if (!eps->V) { EPSGetBV(eps,&eps->V); }
210:     BVView(eps->V,viewer);
211:     if (eps->rg) {
212:       RGIsTrivial(eps->rg,&istrivial);
213:       if (!istrivial) { RGView(eps->rg,viewer); }
214:     }
215:     PetscObjectTypeCompare((PetscObject)eps,EPSPOWER,&ispower);
216:     if (!ispower) {
217:       if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
218:       DSView(eps->ds,viewer);
219:     }
220:     PetscViewerPopFormat(viewer);
221:   }
222:   if (!eps->st) { EPSGetST(eps,&eps->st); }
223:   STView(eps->st,viewer);
224:   return(0);
225: }

229: /*@
230:    EPSPrintSolution - Prints the computed eigenvalues.

232:    Collective on EPS

234:    Input Parameters:
235: +  eps - the eigensolver context
236: -  viewer - optional visualization context

238:    Options Database Key:
239: .  -eps_terse - print only minimal information

241:    Note:
242:    By default, this function prints a table with eigenvalues and associated
243:    relative errors. With -eps_terse only the eigenvalues are printed.

245:    Level: intermediate

247: .seealso: PetscViewerASCIIOpen()
248: @*/
249: PetscErrorCode EPSPrintSolution(EPS eps,PetscViewer viewer)
250: {
251:   PetscBool      terse,errok,isascii;
252:   PetscReal      error,re,im;
253:   PetscScalar    kr,ki;
254:   PetscInt       i,j;

259:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
262:   EPSCheckSolved(eps,1);
263:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
264:   if (!isascii) return(0);

266:   PetscOptionsHasName(NULL,"-eps_terse",&terse);
267:   if (terse) {
268:     if (eps->nconv<eps->nev) {
269:       PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
270:     } else {
271:       errok = PETSC_TRUE;
272:       for (i=0;i<eps->nev;i++) {
273:         EPSComputeRelativeError(eps,i,&error);
274:         errok = (errok && error<5.0*eps->tol)? PETSC_TRUE: PETSC_FALSE;
275:       }
276:       if (errok) {
277:         PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
278:         for (i=0;i<=(eps->nev-1)/8;i++) {
279:           PetscViewerASCIIPrintf(viewer,"\n     ");
280:           for (j=0;j<PetscMin(8,eps->nev-8*i);j++) {
281:             EPSGetEigenpair(eps,8*i+j,&kr,&ki,NULL,NULL);
282: #if defined(PETSC_USE_COMPLEX)
283:             re = PetscRealPart(kr);
284:             im = PetscImaginaryPart(kr);
285: #else
286:             re = kr;
287:             im = ki;
288: #endif
289:             if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
290:             if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
291:             if (im!=0.0) {
292:               PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
293:             } else {
294:               PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
295:             }
296:             if (8*i+j+1<eps->nev) { PetscViewerASCIIPrintf(viewer,", "); }
297:           }
298:         }
299:         PetscViewerASCIIPrintf(viewer,"\n\n");
300:       } else {
301:         PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",eps->nev);
302:       }
303:     }
304:   } else {
305:     PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",eps->nconv);
306:     if (eps->nconv>0) {
307:       PetscViewerASCIIPrintf(viewer,
308:            "           k          ||Ax-k%sx||/||kx||\n"
309:            "   ----------------- ------------------\n",eps->isgeneralized?"B":"");
310:       for (i=0;i<eps->nconv;i++) {
311:         EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
312:         EPSComputeRelativeError(eps,i,&error);
313: #if defined(PETSC_USE_COMPLEX)
314:         re = PetscRealPart(kr);
315:         im = PetscImaginaryPart(kr);
316: #else
317:         re = kr;
318:         im = ki;
319: #endif
320:         if (im!=0.0) {
321:           PetscViewerASCIIPrintf(viewer," % 9f%+9f i %12g\n",(double)re,(double)im,(double)error);
322:         } else {
323:           PetscViewerASCIIPrintf(viewer,"   % 12f       %12g\n",(double)re,(double)error);
324:         }
325:       }
326:       PetscViewerASCIIPrintf(viewer,"\n");
327:     }
328:   }
329:   return(0);
330: }

334: /*@C
335:    EPSCreate - Creates the default EPS context.

337:    Collective on MPI_Comm

339:    Input Parameter:
340: .  comm - MPI communicator

342:    Output Parameter:
343: .  eps - location to put the EPS context

345:    Note:
346:    The default EPS type is EPSKRYLOVSCHUR

348:    Level: beginner

350: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
351: @*/
352: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
353: {
355:   EPS            eps;

359:   *outeps = 0;
360:   EPSInitializePackage();
361:   SlepcHeaderCreate(eps,_p_EPS,struct _EPSOps,EPS_CLASSID,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView);

363:   eps->max_it          = 0;
364:   eps->nev             = 1;
365:   eps->ncv             = 0;
366:   eps->mpd             = 0;
367:   eps->nini            = 0;
368:   eps->nds             = 0;
369:   eps->target          = 0.0;
370:   eps->tol             = PETSC_DEFAULT;
371:   eps->conv            = EPS_CONV_EIG;
372:   eps->which           = (EPSWhich)0;
373:   eps->inta            = 0.0;
374:   eps->intb            = 0.0;
375:   eps->problem_type    = (EPSProblemType)0;
376:   eps->extraction      = EPS_RITZ;
377:   eps->balance         = EPS_BALANCE_NONE;
378:   eps->balance_its     = 5;
379:   eps->balance_cutoff  = 1e-8;
380:   eps->trueres         = PETSC_FALSE;
381:   eps->trackall        = PETSC_FALSE;

383:   eps->converged       = EPSConvergedEigRelative;
384:   eps->convergeddestroy= NULL;
385:   eps->arbitrary       = NULL;
386:   eps->convergedctx    = NULL;
387:   eps->arbitraryctx    = NULL;
388:   eps->numbermonitors  = 0;

390:   eps->st              = NULL;
391:   eps->ds              = NULL;
392:   eps->V               = NULL;
393:   eps->rg              = NULL;
394:   eps->rand            = NULL;
395:   eps->D               = NULL;
396:   eps->IS              = NULL;
397:   eps->defl            = NULL;
398:   eps->eigr            = NULL;
399:   eps->eigi            = NULL;
400:   eps->errest          = NULL;
401:   eps->rr              = NULL;
402:   eps->ri              = NULL;
403:   eps->perm            = NULL;
404:   eps->nwork           = 0;
405:   eps->work            = NULL;
406:   eps->data            = NULL;

408:   eps->state           = EPS_STATE_INITIAL;
409:   eps->nconv           = 0;
410:   eps->its             = 0;
411:   eps->nloc            = 0;
412:   eps->nrma            = 0.0;
413:   eps->nrmb            = 0.0;
414:   eps->isgeneralized   = PETSC_FALSE;
415:   eps->ispositive      = PETSC_FALSE;
416:   eps->ishermitian     = PETSC_FALSE;
417:   eps->reason          = EPS_CONVERGED_ITERATING;

419:   PetscNewLog(eps,&eps->sc);
420:   PetscRandomCreate(comm,&eps->rand);
421:   PetscRandomSetSeed(eps->rand,0x12345678);
422:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rand);
423:   *outeps = eps;
424:   return(0);
425: }

429: /*@C
430:    EPSSetType - Selects the particular solver to be used in the EPS object.

432:    Logically Collective on EPS

434:    Input Parameters:
435: +  eps  - the eigensolver context
436: -  type - a known method

438:    Options Database Key:
439: .  -eps_type <method> - Sets the method; use -help for a list
440:     of available methods

442:    Notes:
443:    See "slepc/include/slepceps.h" for available methods. The default
444:    is EPSKRYLOVSCHUR.

446:    Normally, it is best to use the EPSSetFromOptions() command and
447:    then set the EPS type from the options database rather than by using
448:    this routine.  Using the options database provides the user with
449:    maximum flexibility in evaluating the different available methods.
450:    The EPSSetType() routine is provided for those situations where it
451:    is necessary to set the iterative solver independently of the command
452:    line or options database.

454:    Level: intermediate

456: .seealso: STSetType(), EPSType
457: @*/
458: PetscErrorCode EPSSetType(EPS eps,EPSType type)
459: {
460:   PetscErrorCode ierr,(*r)(EPS);
461:   PetscBool      match;


467:   PetscObjectTypeCompare((PetscObject)eps,type,&match);
468:   if (match) return(0);

470:   PetscFunctionListFind(EPSList,type,&r);
471:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);

473:   if (eps->ops->destroy) { (*eps->ops->destroy)(eps); }
474:   PetscMemzero(eps->ops,sizeof(struct _EPSOps));

476:   eps->state = EPS_STATE_INITIAL;
477:   PetscObjectChangeTypeName((PetscObject)eps,type);
478:   (*r)(eps);
479:   return(0);
480: }

484: /*@C
485:    EPSGetType - Gets the EPS type as a string from the EPS object.

487:    Not Collective

489:    Input Parameter:
490: .  eps - the eigensolver context

492:    Output Parameter:
493: .  name - name of EPS method

495:    Level: intermediate

497: .seealso: EPSSetType()
498: @*/
499: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
500: {
504:   *type = ((PetscObject)eps)->type_name;
505:   return(0);
506: }

510: /*@C
511:    EPSRegister - Adds a method to the eigenproblem solver package.

513:    Not Collective

515:    Input Parameters:
516: +  name - name of a new user-defined solver
517: -  function - routine to create the solver context

519:    Notes:
520:    EPSRegister() may be called multiple times to add several user-defined solvers.

522:    Sample usage:
523: .vb
524:    EPSRegister("my_solver",MySolverCreate);
525: .ve

527:    Then, your solver can be chosen with the procedural interface via
528: $     EPSSetType(eps,"my_solver")
529:    or at runtime via the option
530: $     -eps_type my_solver

532:    Level: advanced

534: .seealso: EPSRegisterAll()
535: @*/
536: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
537: {

541:   PetscFunctionListAdd(&EPSList,name,function);
542:   return(0);
543: }

547: /*@
548:    EPSReset - Resets the EPS context to the initial state and removes any
549:    allocated objects.

551:    Collective on EPS

553:    Input Parameter:
554: .  eps - eigensolver context obtained from EPSCreate()

556:    Level: advanced

558: .seealso: EPSDestroy()
559: @*/
560: PetscErrorCode EPSReset(EPS eps)
561: {
563:   PetscInt       ncols;

567:   if (eps->ops->reset) { (eps->ops->reset)(eps); }
568:   if (eps->st) { STReset(eps->st); }
569:   if (eps->ds) { DSReset(eps->ds); }
570:   VecDestroy(&eps->D);
571:   BVGetSizes(eps->V,NULL,NULL,&ncols);
572:   if (ncols) {
573:     PetscFree4(eps->eigr,eps->eigi,eps->errest,eps->perm);
574:     PetscFree2(eps->rr,eps->ri);
575:   }
576:   BVDestroy(&eps->V);
577:   VecDestroyVecs(eps->nwork,&eps->work);
578:   eps->nwork = 0;
579:   eps->state = EPS_STATE_INITIAL;
580:   return(0);
581: }

585: /*@C
586:    EPSDestroy - Destroys the EPS context.

588:    Collective on EPS

590:    Input Parameter:
591: .  eps - eigensolver context obtained from EPSCreate()

593:    Level: beginner

595: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
596: @*/
597: PetscErrorCode EPSDestroy(EPS *eps)
598: {

602:   if (!*eps) return(0);
604:   if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; return(0); }
605:   EPSReset(*eps);
606:   if ((*eps)->ops->destroy) { (*(*eps)->ops->destroy)(*eps); }
607:   STDestroy(&(*eps)->st);
608:   RGDestroy(&(*eps)->rg);
609:   DSDestroy(&(*eps)->ds);
610:   PetscRandomDestroy(&(*eps)->rand);
611:   PetscFree((*eps)->sc);
612:   /* just in case the initial vectors have not been used */
613:   SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl);
614:   SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS);
615:   if ((*eps)->convergeddestroy) {
616:     (*(*eps)->convergeddestroy)((*eps)->convergedctx);
617:   }
618:   EPSMonitorCancel(*eps);
619:   PetscHeaderDestroy(eps);
620:   return(0);
621: }

625: /*@
626:    EPSSetTarget - Sets the value of the target.

628:    Logically Collective on EPS

630:    Input Parameters:
631: +  eps    - eigensolver context
632: -  target - the value of the target

634:    Options Database Key:
635: .  -eps_target <scalar> - the value of the target

637:    Notes:
638:    The target is a scalar value used to determine the portion of the spectrum
639:    of interest. It is used in combination with EPSSetWhichEigenpairs().

641:    In the case of complex scalars, a complex value can be provided in the
642:    command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
643:    -eps_target 1.0+2.0i

645:    Level: beginner

647: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
648: @*/
649: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
650: {

656:   eps->target = target;
657:   if (!eps->st) { EPSGetST(eps,&eps->st); }
658:   STSetDefaultShift(eps->st,target);
659:   return(0);
660: }

664: /*@
665:    EPSGetTarget - Gets the value of the target.

667:    Not Collective

669:    Input Parameter:
670: .  eps - eigensolver context

672:    Output Parameter:
673: .  target - the value of the target

675:    Note:
676:    If the target was not set by the user, then zero is returned.

678:    Level: beginner

680: .seealso: EPSSetTarget()
681: @*/
682: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
683: {
687:   *target = eps->target;
688:   return(0);
689: }

693: /*@
694:    EPSSetInterval - Defines the computational interval for spectrum slicing.

696:    Logically Collective on EPS

698:    Input Parameters:
699: +  eps  - eigensolver context
700: .  inta - left end of the interval
701: -  intb - right end of the interval

703:    Options Database Key:
704: .  -eps_interval <a,b> - set [a,b] as the interval of interest

706:    Notes:
707:    Spectrum slicing is a technique employed for computing all eigenvalues of
708:    symmetric eigenproblems in a given interval. This function provides the
709:    interval to be considered. It must be used in combination with EPS_ALL, see
710:    EPSSetWhichEigenpairs().

712:    In the command-line option, two values must be provided. For an open interval,
713:    one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
714:    An open interval in the programmatic interface can be specified with
715:    PETSC_MAX_REAL and -PETSC_MAX_REAL.

717:    Level: intermediate

719: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
720: @*/
721: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
722: {
727:   if (inta>=intb) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
728:   eps->inta = inta;
729:   eps->intb = intb;
730:   eps->state = EPS_STATE_INITIAL;
731:   return(0);
732: }

736: /*@
737:    EPSGetInterval - Gets the computational interval for spectrum slicing.

739:    Not Collective

741:    Input Parameter:
742: .  eps - eigensolver context

744:    Output Parameters:
745: +  inta - left end of the interval
746: -  intb - right end of the interval

748:    Level: intermediate

750:    Note:
751:    If the interval was not set by the user, then zeros are returned.

753: .seealso: EPSSetInterval()
754: @*/
755: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
756: {
761:   if (inta) *inta = eps->inta;
762:   if (intb) *intb = eps->intb;
763:   return(0);
764: }

768: /*@
769:    EPSSetST - Associates a spectral transformation object to the eigensolver.

771:    Collective on EPS

773:    Input Parameters:
774: +  eps - eigensolver context obtained from EPSCreate()
775: -  st   - the spectral transformation object

777:    Note:
778:    Use EPSGetST() to retrieve the spectral transformation context (for example,
779:    to free it at the end of the computations).

781:    Level: developer

783: .seealso: EPSGetST()
784: @*/
785: PetscErrorCode EPSSetST(EPS eps,ST st)
786: {

793:   PetscObjectReference((PetscObject)st);
794:   STDestroy(&eps->st);
795:   eps->st = st;
796:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
797:   return(0);
798: }

802: /*@C
803:    EPSGetST - Obtain the spectral transformation (ST) object associated
804:    to the eigensolver object.

806:    Not Collective

808:    Input Parameters:
809: .  eps - eigensolver context obtained from EPSCreate()

811:    Output Parameter:
812: .  st - spectral transformation context

814:    Level: beginner

816: .seealso: EPSSetST()
817: @*/
818: PetscErrorCode EPSGetST(EPS eps,ST *st)
819: {

825:   if (!eps->st) {
826:     STCreate(PetscObjectComm((PetscObject)eps),&eps->st);
827:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
828:   }
829:   *st = eps->st;
830:   return(0);
831: }

835: /*@
836:    EPSSetBV - Associates a basis vectors object to the eigensolver.

838:    Collective on EPS

840:    Input Parameters:
841: +  eps - eigensolver context obtained from EPSCreate()
842: -  V   - the basis vectors object

844:    Note:
845:    Use EPSGetBV() to retrieve the basis vectors context (for example,
846:    to free them at the end of the computations).

848:    Level: advanced

850: .seealso: EPSGetBV()
851: @*/
852: PetscErrorCode EPSSetBV(EPS eps,BV V)
853: {

860:   PetscObjectReference((PetscObject)V);
861:   BVDestroy(&eps->V);
862:   eps->V = V;
863:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
864:   return(0);
865: }

869: /*@C
870:    EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.

872:    Not Collective

874:    Input Parameters:
875: .  eps - eigensolver context obtained from EPSCreate()

877:    Output Parameter:
878: .  V - basis vectors context

880:    Level: advanced

882: .seealso: EPSSetBV()
883: @*/
884: PetscErrorCode EPSGetBV(EPS eps,BV *V)
885: {

891:   if (!eps->V) {
892:     BVCreate(PetscObjectComm((PetscObject)eps),&eps->V);
893:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
894:   }
895:   *V = eps->V;
896:   return(0);
897: }

901: /*@
902:    EPSSetRG - Associates a region object to the eigensolver.

904:    Collective on EPS

906:    Input Parameters:
907: +  eps - eigensolver context obtained from EPSCreate()
908: -  rg  - the region object

910:    Note:
911:    Use EPSGetRG() to retrieve the region context (for example,
912:    to free it at the end of the computations).

914:    Level: advanced

916: .seealso: EPSGetRG()
917: @*/
918: PetscErrorCode EPSSetRG(EPS eps,RG rg)
919: {

926:   PetscObjectReference((PetscObject)rg);
927:   RGDestroy(&eps->rg);
928:   eps->rg = rg;
929:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
930:   return(0);
931: }

935: /*@C
936:    EPSGetRG - Obtain the region object associated to the eigensolver.

938:    Not Collective

940:    Input Parameters:
941: .  eps - eigensolver context obtained from EPSCreate()

943:    Output Parameter:
944: .  rg - region context

946:    Level: advanced

948: .seealso: EPSSetRG()
949: @*/
950: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
951: {

957:   if (!eps->rg) {
958:     RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg);
959:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
960:   }
961:   *rg = eps->rg;
962:   return(0);
963: }

967: /*@
968:    EPSSetDS - Associates a direct solver object to the eigensolver.

970:    Collective on EPS

972:    Input Parameters:
973: +  eps - eigensolver context obtained from EPSCreate()
974: -  ds  - the direct solver object

976:    Note:
977:    Use EPSGetDS() to retrieve the direct solver context (for example,
978:    to free it at the end of the computations).

980:    Level: advanced

982: .seealso: EPSGetDS()
983: @*/
984: PetscErrorCode EPSSetDS(EPS eps,DS ds)
985: {

992:   PetscObjectReference((PetscObject)ds);
993:   DSDestroy(&eps->ds);
994:   eps->ds = ds;
995:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
996:   return(0);
997: }

1001: /*@C
1002:    EPSGetDS - Obtain the direct solver object associated to the eigensolver object.

1004:    Not Collective

1006:    Input Parameters:
1007: .  eps - eigensolver context obtained from EPSCreate()

1009:    Output Parameter:
1010: .  ds - direct solver context

1012:    Level: advanced

1014: .seealso: EPSSetDS()
1015: @*/
1016: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
1017: {

1023:   if (!eps->ds) {
1024:     DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds);
1025:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
1026:   }
1027:   *ds = eps->ds;
1028:   return(0);
1029: }

1033: /*@
1034:    EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
1035:    eigenvalue problem.

1037:    Not collective

1039:    Input Parameter:
1040: .  eps - the eigenproblem solver context

1042:    Output Parameter:
1043: .  is - the answer

1045:    Level: intermediate

1047: .seealso: EPSIsHermitian(), EPSIsPositive()
1048: @*/
1049: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
1050: {
1054:   *is = eps->isgeneralized;
1055:   return(0);
1056: }

1060: /*@
1061:    EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
1062:    eigenvalue problem.

1064:    Not collective

1066:    Input Parameter:
1067: .  eps - the eigenproblem solver context

1069:    Output Parameter:
1070: .  is - the answer

1072:    Level: intermediate

1074: .seealso: EPSIsGeneralized(), EPSIsPositive()
1075: @*/
1076: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
1077: {
1081:   *is = eps->ishermitian;
1082:   return(0);
1083: }

1087: /*@
1088:    EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
1089:    problem type that requires a positive (semi-) definite matrix B.

1091:    Not collective

1093:    Input Parameter:
1094: .  eps - the eigenproblem solver context

1096:    Output Parameter:
1097: .  is - the answer

1099:    Level: intermediate

1101: .seealso: EPSIsGeneralized(), EPSIsHermitian()
1102: @*/
1103: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
1104: {
1108:   *is = eps->ispositive;
1109:   return(0);
1110: }