Actual source code: epsbasic.c

slepc-3.5.2 2014-10-10
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) { EPSGetRG(eps,&eps->rg); }
212:     RGIsTrivial(eps->rg,&istrivial);
213:     if (!istrivial) { RGView(eps->rg,viewer); }
214:     PetscObjectTypeCompare((PetscObject)eps,EPSPOWER,&ispower);
215:     if (!ispower) {
216:       if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
217:       DSView(eps->ds,viewer);
218:     }
219:     PetscViewerPopFormat(viewer);
220:   }
221:   if (!eps->st) { EPSGetST(eps,&eps->st); }
222:   STView(eps->st,viewer);
223:   return(0);
224: }

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

231:    Collective on EPS

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

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

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

244:    Level: intermediate

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

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

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

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

336:    Collective on MPI_Comm

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

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

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

347:    Level: beginner

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

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

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

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

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

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

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

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

431:    Logically Collective on EPS

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

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

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

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

453:    Level: intermediate

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


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

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

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

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

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

486:    Not Collective

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

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

494:    Level: intermediate

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

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

512:    Not Collective

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

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

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

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

531:    Level: advanced

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

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

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

550:    Collective on EPS

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

555:    Level: advanced

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

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

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

587:    Collective on EPS

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

592:    Level: beginner

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

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

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

627:    Logically Collective on EPS

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

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

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

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

644:    Level: beginner

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

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

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

666:    Not Collective

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

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

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

677:    Level: beginner

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

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

695:    Logically Collective on EPS

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

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

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

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

716:    Level: intermediate

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

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

738:    Not Collective

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

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

747:    Level: intermediate

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

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

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

770:    Collective on EPS

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

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

780:    Level: developer

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

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

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

805:    Not Collective

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

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

813:    Level: beginner

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

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

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

837:    Collective on EPS

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

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

847:    Level: advanced

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

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

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

871:    Not Collective

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

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

879:    Level: advanced

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

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

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

903:    Collective on EPS

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

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

913:    Level: advanced

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

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

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

937:    Not Collective

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

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

945:    Level: advanced

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

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

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

969:    Collective on EPS

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

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

979:    Level: advanced

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

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

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

1003:    Not Collective

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

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

1011:    Level: advanced

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

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

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

1036:    Not collective

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

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

1044:    Level: intermediate

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

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

1063:    Not collective

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

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

1071:    Level: intermediate

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

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

1090:    Not collective

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

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

1098:    Level: intermediate

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