Actual source code: pepbasic.c

slepc-3.5.2 2014-10-10
Report Typos and Errors
  1: /*
  2:      The basic PEP 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/pepimpl.h>      /*I "slepcpep.h" I*/

 26: PetscFunctionList PEPList = 0;
 27: PetscBool         PEPRegisterAllCalled = PETSC_FALSE;
 28: PetscClassId      PEP_CLASSID = 0;
 29: PetscLogEvent     PEP_SetUp = 0,PEP_Solve = 0,PEP_Refine = 0;

 33: /*@C
 34:    PEPView - Prints the PEP data structure.

 36:    Collective on PEP

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

 42:    Options Database Key:
 43: .  -pep_view -  Calls PEPView() at end of PEPSolve()

 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: PetscViewerASCIIOpen()
 59: @*/
 60: PetscErrorCode PEPView(PEP pep,PetscViewer viewer)
 61: {
 63:   const char     *type;
 64:   char           str[50];
 65:   PetscBool      isascii,islinear,istrivial;
 66:   PetscInt       i;

 70:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));

 74: #if defined(PETSC_USE_COMPLEX)
 75: #define HERM "hermitian"
 76: #else
 77: #define HERM "symmetric"
 78: #endif
 79:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 80:   if (isascii) {
 81:     PetscObjectPrintClassNamePrefixType((PetscObject)pep,viewer);
 82:     if (pep->ops->view) {
 83:       PetscViewerASCIIPushTab(viewer);
 84:       (*pep->ops->view)(pep,viewer);
 85:       PetscViewerASCIIPopTab(viewer);
 86:     }
 87:     if (pep->problem_type) {
 88:       switch (pep->problem_type) {
 89:         case PEP_GENERAL:    type = "general polynomial eigenvalue problem"; break;
 90:         case PEP_HERMITIAN:  type = HERM " polynomial eigenvalue problem"; break;
 91:         case PEP_GYROSCOPIC: type = "gyroscopic polynomial eigenvalue problem"; break;
 92:         default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->problem_type");
 93:       }
 94:     } else type = "not yet set";
 95:     PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type);
 96:     PetscViewerASCIIPrintf(viewer,"  polynomial represented in %s basis\n",PEPBasisTypes[pep->basis]);
 97:     switch (pep->scale) {
 98:       case PEP_SCALE_NONE:
 99:         break;
100:       case PEP_SCALE_SCALAR:
101:         PetscViewerASCIIPrintf(viewer,"  scalar balancing enabled, with scaling factor=%g\n",(double)pep->sfactor);
102:         break;
103:       case PEP_SCALE_DIAGONAL:
104:         PetscViewerASCIIPrintf(viewer,"  diagonal balancing enabled, with its=%D and lambda=%g\n",pep->sits,(double)pep->slambda);
105:         break;
106:       case PEP_SCALE_BOTH:
107:         PetscViewerASCIIPrintf(viewer,"  scalar & diagonal balancing enabled, with scaling factor=%g, its=%D and lambda=%g\n",(double)pep->sfactor,pep->sits,(double)pep->slambda);
108:         break;
109:     }
110:     PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",PEPExtractTypes[pep->extract]);
111:     PetscViewerASCIIPrintf(viewer,"  iterative refinement: %s%s\n",PEPRefineTypes[pep->refine],pep->schur?", with a Schur complement approach":"");
112:     if (pep->refine) {
113:       PetscViewerASCIIPrintf(viewer,"  refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
114:       if (pep->npart>1) {
115:         PetscViewerASCIIPrintf(viewer,"  splitting communicator in %D partitions for refinement\n",pep->npart);
116:       }
117:     }
118:     PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: ");
119:     SlepcSNPrintfScalar(str,50,pep->target,PETSC_FALSE);
120:     if (!pep->which) {
121:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
122:     } else switch (pep->which) {
123:       case PEP_TARGET_MAGNITUDE:
124:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
125:         break;
126:       case PEP_TARGET_REAL:
127:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
128:         break;
129:       case PEP_TARGET_IMAGINARY:
130:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
131:         break;
132:       case PEP_LARGEST_MAGNITUDE:
133:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
134:         break;
135:       case PEP_SMALLEST_MAGNITUDE:
136:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
137:         break;
138:       case PEP_LARGEST_REAL:
139:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
140:         break;
141:       case PEP_SMALLEST_REAL:
142:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
143:         break;
144:       case PEP_LARGEST_IMAGINARY:
145:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
146:         break;
147:       case PEP_SMALLEST_IMAGINARY:
148:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
149:         break;
150:       default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->which");
151:     }
152:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %D\n",pep->nev);
153:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",pep->ncv);
154:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %D\n",pep->mpd);
155:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",pep->max_it);
156:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)pep->tol);
157:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
158:     switch (pep->conv) {
159:     case PEP_CONV_ABS:
160:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
161:     case PEP_CONV_EIG:
162:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
163:     case PEP_CONV_NORM:
164:       PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
165:       if (pep->nrma) {
166:         PetscViewerASCIIPrintf(viewer,"  computed matrix norms: %g",(double)pep->nrma[0]);
167:         for (i=1;i<pep->nmat;i++) {
168:           PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
169:         }
170:         PetscViewerASCIIPrintf(viewer,"\n");
171:       }
172:       break;
173:     case PEP_CONV_USER:
174:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
175:     }
176:     if (pep->nini) {
177:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
178:     }
179:   } else {
180:     if (pep->ops->view) {
181:       (*pep->ops->view)(pep,viewer);
182:     }
183:   }
184:   PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
185:   if (!islinear) {
186:     PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
187:     if (!pep->V) { PEPGetBV(pep,&pep->V); }
188:     BVView(pep->V,viewer);
189:     if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
190:     RGIsTrivial(pep->rg,&istrivial);
191:     if (!istrivial) { RGView(pep->rg,viewer); }
192:     if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
193:     DSView(pep->ds,viewer);
194:     PetscViewerPopFormat(viewer);
195:     if (!pep->st) { PEPGetST(pep,&pep->st); }
196:     STView(pep->st,viewer);
197:   }
198:   return(0);
199: }

203: /*@
204:    PEPPrintSolution - Prints the computed eigenvalues.

206:    Collective on PEP

208:    Input Parameters:
209: +  pep - the eigensolver context
210: -  viewer - optional visualization context

212:    Options Database Key:
213: .  -pep_terse - print only minimal information

215:    Note:
216:    By default, this function prints a table with eigenvalues and associated
217:    relative errors. With -pep_terse only the eigenvalues are printed.

219:    Level: intermediate

221: .seealso: PetscViewerASCIIOpen()
222: @*/
223: PetscErrorCode PEPPrintSolution(PEP pep,PetscViewer viewer)
224: {
225:   PetscBool      terse,errok,isascii;
226:   PetscReal      error,re,im;
227:   PetscScalar    kr,ki;
228:   PetscInt       i,j;

233:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
236:   PEPCheckSolved(pep,1);
237:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
238:   if (!isascii) return(0);

240:   PetscOptionsHasName(NULL,"-pep_terse",&terse);
241:   if (terse) {
242:     if (pep->nconv<pep->nev) {
243:       PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
244:     } else {
245:       errok = PETSC_TRUE;
246:       for (i=0;i<pep->nev;i++) {
247:         PEPComputeRelativeError(pep,i,&error);
248:         errok = (errok && error<5.0*pep->tol)? PETSC_TRUE: PETSC_FALSE;
249:       }
250:       if (errok) {
251:         PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
252:         for (i=0;i<=(pep->nev-1)/8;i++) {
253:           PetscViewerASCIIPrintf(viewer,"\n     ");
254:           for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
255:             PEPGetEigenpair(pep,8*i+j,&kr,&ki,NULL,NULL);
256: #if defined(PETSC_USE_COMPLEX)
257:             re = PetscRealPart(kr);
258:             im = PetscImaginaryPart(kr);
259: #else
260:             re = kr;
261:             im = ki;
262: #endif
263:             if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
264:             if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
265:             if (im!=0.0) {
266:               PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
267:             } else {
268:               PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
269:             }
270:             if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
271:           }
272:         }
273:         PetscViewerASCIIPrintf(viewer,"\n\n");
274:       } else {
275:         PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
276:       }
277:     }
278:   } else {
279:     PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",pep->nconv);
280:     if (pep->nconv>0) {
281:       PetscViewerASCIIPrintf(viewer,
282:            "           k             ||P(k)x||/||kx||\n"
283:            "   ----------------- -------------------------\n");
284:       for (i=0;i<pep->nconv;i++) {
285:         PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
286:         PEPComputeRelativeError(pep,i,&error);
287: #if defined(PETSC_USE_COMPLEX)
288:         re = PetscRealPart(kr);
289:         im = PetscImaginaryPart(kr);
290: #else
291:         re = kr;
292:         im = ki;
293: #endif
294:         if (im!=0.0) {
295:           PetscViewerASCIIPrintf(viewer," % 9f%+9f i     %12g\n",(double)re,(double)im,(double)error);
296:         } else {
297:           PetscViewerASCIIPrintf(viewer,"   % 12f           %12g\n",(double)re,(double)error);
298:         }
299:       }
300:       PetscViewerASCIIPrintf(viewer,"\n");
301:     }
302:   }
303:   return(0);
304: }

308: /*@C
309:    PEPCreate - Creates the default PEP context.

311:    Collective on MPI_Comm

313:    Input Parameter:
314: .  comm - MPI communicator

316:    Output Parameter:
317: .  pep - location to put the PEP context

319:    Note:
320:    The default PEP type is PEPTOAR

322:    Level: beginner

324: .seealso: PEPSetUp(), PEPSolve(), PEPDestroy(), PEP
325: @*/
326: PetscErrorCode PEPCreate(MPI_Comm comm,PEP *outpep)
327: {
329:   PEP            pep;

333:   *outpep = 0;
334:   PEPInitializePackage();
335:   SlepcHeaderCreate(pep,_p_PEP,struct _PEPOps,PEP_CLASSID,"PEP","Polynomial Eigenvalue Problem","PEP",comm,PEPDestroy,PEPView);

337:   pep->max_it          = 0;
338:   pep->nev             = 1;
339:   pep->ncv             = 0;
340:   pep->mpd             = 0;
341:   pep->nini            = 0;
342:   pep->target          = 0.0;
343:   pep->tol             = PETSC_DEFAULT;
344:   pep->conv            = PEP_CONV_NORM;
345:   pep->which           = (PEPWhich)0;
346:   pep->basis           = PEP_BASIS_MONOMIAL;
347:   pep->problem_type    = (PEPProblemType)0;
348:   pep->scale           = PEP_SCALE_NONE;
349:   pep->sfactor         = 1.0;
350:   pep->sits            = 5;
351:   pep->slambda         = 1.0;
352:   pep->refine          = PEP_REFINE_NONE;
353:   pep->npart           = 1;
354:   pep->rtol            = PETSC_DEFAULT;
355:   pep->rits            = PETSC_DEFAULT;
356:   pep->schur           = PETSC_FALSE;
357:   pep->extract         = PEP_EXTRACT_NORM;
358:   pep->trackall        = PETSC_FALSE;

360:   pep->converged       = PEPConvergedNormRelative;
361:   pep->convergeddestroy= NULL;
362:   pep->convergedctx    = NULL;
363:   pep->numbermonitors  = 0;

365:   pep->st              = NULL;
366:   pep->ds              = NULL;
367:   pep->V               = NULL;
368:   pep->rg              = NULL;
369:   pep->rand            = NULL;
370:   pep->A               = NULL;
371:   pep->nmat            = 0;
372:   pep->Dl              = NULL;
373:   pep->Dr              = NULL;
374:   pep->IS              = NULL;
375:   pep->eigr            = NULL;
376:   pep->eigi            = NULL;
377:   pep->errest          = NULL;
378:   pep->perm            = NULL;
379:   pep->pbc             = NULL;
380:   pep->solvematcoeffs  = NULL;
381:   pep->nwork           = 0;
382:   pep->work            = NULL;
383:   pep->data            = NULL;

385:   pep->state           = PEP_STATE_INITIAL;
386:   pep->nconv           = 0;
387:   pep->its             = 0;
388:   pep->n               = 0;
389:   pep->nloc            = 0;
390:   pep->nrma            = NULL;
391:   pep->sfactor_set     = PETSC_FALSE;
392:   pep->reason          = PEP_CONVERGED_ITERATING;

394:   PetscNewLog(pep,&pep->sc);
395:   PetscRandomCreate(comm,&pep->rand);
396:   PetscRandomSetSeed(pep->rand,0x12345678);
397:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rand);
398:   *outpep = pep;
399:   return(0);
400: }

404: /*@C
405:    PEPSetType - Selects the particular solver to be used in the PEP object.

407:    Logically Collective on PEP

409:    Input Parameters:
410: +  pep      - the polynomial eigensolver context
411: -  type     - a known method

413:    Options Database Key:
414: .  -pep_type <method> - Sets the method; use -help for a list
415:     of available methods

417:    Notes:
418:    See "slepc/include/slepcpep.h" for available methods. The default
419:    is PEPTOAR.

421:    Normally, it is best to use the PEPSetFromOptions() command and
422:    then set the PEP type from the options database rather than by using
423:    this routine.  Using the options database provides the user with
424:    maximum flexibility in evaluating the different available methods.
425:    The PEPSetType() routine is provided for those situations where it
426:    is necessary to set the iterative solver independently of the command
427:    line or options database.

429:    Level: intermediate

431: .seealso: PEPType
432: @*/
433: PetscErrorCode PEPSetType(PEP pep,PEPType type)
434: {
435:   PetscErrorCode ierr,(*r)(PEP);
436:   PetscBool      match;


442:   PetscObjectTypeCompare((PetscObject)pep,type,&match);
443:   if (match) return(0);

445:   PetscFunctionListFind(PEPList,type,&r);
446:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown PEP type given: %s",type);

448:   if (pep->ops->destroy) { (*pep->ops->destroy)(pep); }
449:   PetscMemzero(pep->ops,sizeof(struct _PEPOps));

451:   pep->state = PEP_STATE_INITIAL;
452:   PetscObjectChangeTypeName((PetscObject)pep,type);
453:   (*r)(pep);
454:   return(0);
455: }

459: /*@C
460:    PEPGetType - Gets the PEP type as a string from the PEP object.

462:    Not Collective

464:    Input Parameter:
465: .  pep - the eigensolver context

467:    Output Parameter:
468: .  name - name of PEP method

470:    Level: intermediate

472: .seealso: PEPSetType()
473: @*/
474: PetscErrorCode PEPGetType(PEP pep,PEPType *type)
475: {
479:   *type = ((PetscObject)pep)->type_name;
480:   return(0);
481: }

485: /*@C
486:    PEPRegister - Adds a method to the polynomial eigenproblem solver package.

488:    Not Collective

490:    Input Parameters:
491: +  name - name of a new user-defined solver
492: -  function - routine to create the solver context

494:    Notes:
495:    PEPRegister() may be called multiple times to add several user-defined solvers.

497:    Sample usage:
498: .vb
499:    PEPRegister("my_solver",MySolverCreate);
500: .ve

502:    Then, your solver can be chosen with the procedural interface via
503: $     PEPSetType(pep,"my_solver")
504:    or at runtime via the option
505: $     -pep_type my_solver

507:    Level: advanced

509: .seealso: PEPRegisterAll()
510: @*/
511: PetscErrorCode PEPRegister(const char *name,PetscErrorCode (*function)(PEP))
512: {

516:   PetscFunctionListAdd(&PEPList,name,function);
517:   return(0);
518: }

522: /*@
523:    PEPReset - Resets the PEP context to the initial state and removes any
524:    allocated objects.

526:    Collective on PEP

528:    Input Parameter:
529: .  pep - eigensolver context obtained from PEPCreate()

531:    Level: advanced

533: .seealso: PEPDestroy()
534: @*/
535: PetscErrorCode PEPReset(PEP pep)
536: {
538:   PetscInt       ncols;

542:   if (pep->ops->reset) { (pep->ops->reset)(pep); }
543:   if (pep->st) { STReset(pep->st); }
544:   if (pep->ds) { DSReset(pep->ds); }
545:   if (pep->nmat) {
546:     MatDestroyMatrices(pep->nmat,&pep->A);
547:     PetscFree3(pep->pbc,pep->solvematcoeffs,pep->nrma);
548:     pep->nmat = 0;
549:   }
550:   VecDestroy(&pep->Dl);
551:   VecDestroy(&pep->Dr);
552:   BVGetSizes(pep->V,NULL,NULL,&ncols);
553:   if (ncols) {
554:     PetscFree4(pep->eigr,pep->eigi,pep->errest,pep->perm);
555:   }
556:   BVDestroy(&pep->V);
557:   VecDestroyVecs(pep->nwork,&pep->work);
558:   pep->nwork = 0;
559:   pep->state = PEP_STATE_INITIAL;
560:   return(0);
561: }

565: /*@C
566:    PEPDestroy - Destroys the PEP context.

568:    Collective on PEP

570:    Input Parameter:
571: .  pep - eigensolver context obtained from PEPCreate()

573:    Level: beginner

575: .seealso: PEPCreate(), PEPSetUp(), PEPSolve()
576: @*/
577: PetscErrorCode PEPDestroy(PEP *pep)
578: {

582:   if (!*pep) return(0);
584:   if (--((PetscObject)(*pep))->refct > 0) { *pep = 0; return(0); }
585:   PEPReset(*pep);
586:   if ((*pep)->ops->destroy) { (*(*pep)->ops->destroy)(*pep); }
587:   STDestroy(&(*pep)->st);
588:   RGDestroy(&(*pep)->rg);
589:   DSDestroy(&(*pep)->ds);
590:   PetscRandomDestroy(&(*pep)->rand);
591:   PetscFree((*pep)->sc);
592:   /* just in case the initial vectors have not been used */
593:   SlepcBasisDestroy_Private(&(*pep)->nini,&(*pep)->IS);
594:   if ((*pep)->convergeddestroy) {
595:     (*(*pep)->convergeddestroy)((*pep)->convergedctx);
596:   }
597:   PEPMonitorCancel(*pep);
598:   PetscHeaderDestroy(pep);
599:   return(0);
600: }

604: /*@
605:    PEPSetBV - Associates a basis vectors object to the polynomial eigensolver.

607:    Collective on PEP

609:    Input Parameters:
610: +  pep - eigensolver context obtained from PEPCreate()
611: -  bv  - the basis vectors object

613:    Note:
614:    Use PEPGetBV() to retrieve the basis vectors context (for example,
615:    to free it at the end of the computations).

617:    Level: advanced

619: .seealso: PEPGetBV()
620: @*/
621: PetscErrorCode PEPSetBV(PEP pep,BV bv)
622: {

629:   PetscObjectReference((PetscObject)bv);
630:   BVDestroy(&pep->V);
631:   pep->V = bv;
632:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->V);
633:   return(0);
634: }

638: /*@C
639:    PEPGetBV - Obtain the basis vectors object associated to the polynomial
640:    eigensolver object.

642:    Not Collective

644:    Input Parameters:
645: .  pep - eigensolver context obtained from PEPCreate()

647:    Output Parameter:
648: .  bv - basis vectors context

650:    Level: advanced

652: .seealso: PEPSetBV()
653: @*/
654: PetscErrorCode PEPGetBV(PEP pep,BV *bv)
655: {

661:   if (!pep->V) {
662:     BVCreate(PetscObjectComm((PetscObject)pep),&pep->V);
663:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->V);
664:   }
665:   *bv = pep->V;
666:   return(0);
667: }

671: /*@
672:    PEPSetRG - Associates a region object to the polynomial eigensolver.

674:    Collective on PEP

676:    Input Parameters:
677: +  pep - eigensolver context obtained from PEPCreate()
678: -  rg  - the region object

680:    Note:
681:    Use PEPGetRG() to retrieve the region context (for example,
682:    to free it at the end of the computations).

684:    Level: advanced

686: .seealso: PEPGetRG()
687: @*/
688: PetscErrorCode PEPSetRG(PEP pep,RG rg)
689: {

696:   PetscObjectReference((PetscObject)rg);
697:   RGDestroy(&pep->rg);
698:   pep->rg = rg;
699:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rg);
700:   return(0);
701: }

705: /*@C
706:    PEPGetRG - Obtain the region object associated to the
707:    polynomial eigensolver object.

709:    Not Collective

711:    Input Parameters:
712: .  pep - eigensolver context obtained from PEPCreate()

714:    Output Parameter:
715: .  rg - region context

717:    Level: advanced

719: .seealso: PEPSetRG()
720: @*/
721: PetscErrorCode PEPGetRG(PEP pep,RG *rg)
722: {

728:   if (!pep->rg) {
729:     RGCreate(PetscObjectComm((PetscObject)pep),&pep->rg);
730:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rg);
731:   }
732:   *rg = pep->rg;
733:   return(0);
734: }

738: /*@
739:    PEPSetDS - Associates a direct solver object to the polynomial eigensolver.

741:    Collective on PEP

743:    Input Parameters:
744: +  pep - eigensolver context obtained from PEPCreate()
745: -  ds  - the direct solver object

747:    Note:
748:    Use PEPGetDS() to retrieve the direct solver context (for example,
749:    to free it at the end of the computations).

751:    Level: advanced

753: .seealso: PEPGetDS()
754: @*/
755: PetscErrorCode PEPSetDS(PEP pep,DS ds)
756: {

763:   PetscObjectReference((PetscObject)ds);
764:   DSDestroy(&pep->ds);
765:   pep->ds = ds;
766:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->ds);
767:   return(0);
768: }

772: /*@C
773:    PEPGetDS - Obtain the direct solver object associated to the
774:    polynomial eigensolver object.

776:    Not Collective

778:    Input Parameters:
779: .  pep - eigensolver context obtained from PEPCreate()

781:    Output Parameter:
782: .  ds - direct solver context

784:    Level: advanced

786: .seealso: PEPSetDS()
787: @*/
788: PetscErrorCode PEPGetDS(PEP pep,DS *ds)
789: {

795:   if (!pep->ds) {
796:     DSCreate(PetscObjectComm((PetscObject)pep),&pep->ds);
797:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->ds);
798:   }
799:   *ds = pep->ds;
800:   return(0);
801: }

805: /*@
806:    PEPSetST - Associates a spectral transformation object to the eigensolver.

808:    Collective on PEP

810:    Input Parameters:
811: +  pep - eigensolver context obtained from PEPCreate()
812: -  st   - the spectral transformation object

814:    Note:
815:    Use PEPGetST() to retrieve the spectral transformation context (for example,
816:    to free it at the end of the computations).

818:    Level: developer

820: .seealso: PEPGetST()
821: @*/
822: PetscErrorCode PEPSetST(PEP pep,ST st)
823: {

830:   PetscObjectReference((PetscObject)st);
831:   STDestroy(&pep->st);
832:   pep->st = st;
833:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->st);
834:   return(0);
835: }

839: /*@C
840:    PEPGetST - Obtain the spectral transformation (ST) object associated
841:    to the eigensolver object.

843:    Not Collective

845:    Input Parameters:
846: .  pep - eigensolver context obtained from PEPCreate()

848:    Output Parameter:
849: .  st - spectral transformation context

851:    Level: beginner

853: .seealso: PEPSetST()
854: @*/
855: PetscErrorCode PEPGetST(PEP pep,ST *st)
856: {

862:   if (!pep->st) {
863:     STCreate(PetscObjectComm((PetscObject)pep),&pep->st);
864:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->st);
865:   }
866:   *st = pep->st;
867:   return(0);
868: }

872: /*@
873:    PEPSetTarget - Sets the value of the target.

875:    Logically Collective on PEP

877:    Input Parameters:
878: +  pep    - eigensolver context
879: -  target - the value of the target

881:    Options Database Key:
882: .  -pep_target <scalar> - the value of the target

884:    Notes:
885:    The target is a scalar value used to determine the portion of the spectrum
886:    of interest. It is used in combination with PEPSetWhichEigenpairs().

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

892:    Level: beginner

894: .seealso: PEPGetTarget(), PEPSetWhichEigenpairs()
895: @*/
896: PetscErrorCode PEPSetTarget(PEP pep,PetscScalar target)
897: {

903:   pep->target = target;
904:   if (!pep->st) { PEPGetST(pep,&pep->st); }
905:   STSetDefaultShift(pep->st,target);
906:   return(0);
907: }

911: /*@
912:    PEPGetTarget - Gets the value of the target.

914:    Not Collective

916:    Input Parameter:
917: .  pep - eigensolver context

919:    Output Parameter:
920: .  target - the value of the target

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

925:    Level: beginner

927: .seealso: PEPSetTarget()
928: @*/
929: PetscErrorCode PEPGetTarget(PEP pep,PetscScalar* target)
930: {
934:   *target = pep->target;
935:   return(0);
936: }