Actual source code: mfnbasic.c

slepc-3.9.2 2018-07-02
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2018, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    Basic MFN routines
 12: */

 14: #include <slepc/private/mfnimpl.h>      /*I "slepcmfn.h" I*/

 16: PetscFunctionList MFNList = 0;
 17: PetscBool         MFNRegisterAllCalled = PETSC_FALSE;
 18: PetscClassId      MFN_CLASSID = 0;
 19: PetscLogEvent     MFN_SetUp = 0,MFN_Solve = 0;

 21: /*@C
 22:    MFNView - Prints the MFN data structure.

 24:    Collective on MFN

 26:    Input Parameters:
 27: +  mfn - the matrix function solver context
 28: -  viewer - optional visualization context

 30:    Options Database Key:
 31: .  -mfn_view -  Calls MFNView() at end of MFNSolve()

 33:    Note:
 34:    The available visualization contexts include
 35: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 36: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 37:          output where only the first processor opens
 38:          the file.  All other processors send their
 39:          data to the first processor to print.

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

 44:    Level: beginner

 46: .seealso: PetscViewerASCIIOpen()
 47: @*/
 48: PetscErrorCode MFNView(MFN mfn,PetscViewer viewer)
 49: {
 51:   PetscBool      isascii;
 52:   PetscInt       tabs;

 56:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)mfn));

 60:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 61:   if (isascii) {
 62:     PetscViewerASCIIGetTab(viewer,&tabs);
 63:     PetscViewerASCIISetTab(viewer,((PetscObject)mfn)->tablevel);
 64:     PetscObjectPrintClassNamePrefixType((PetscObject)mfn,viewer);
 65:     if (mfn->ops->view) {
 66:       PetscViewerASCIIPushTab(viewer);
 67:       (*mfn->ops->view)(mfn,viewer);
 68:       PetscViewerASCIIPopTab(viewer);
 69:     }
 70:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",mfn->ncv);
 71:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",mfn->max_it);
 72:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)mfn->tol);
 73:     PetscViewerASCIISetTab(viewer,tabs);
 74:   } else {
 75:     if (mfn->ops->view) {
 76:       (*mfn->ops->view)(mfn,viewer);
 77:     }
 78:   }
 79:   PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
 80:   if (!mfn->V) { MFNGetFN(mfn,&mfn->fn); }
 81:   FNView(mfn->fn,viewer);
 82:   if (!mfn->V) { MFNGetBV(mfn,&mfn->V); }
 83:   BVView(mfn->V,viewer);
 84:   PetscViewerPopFormat(viewer);
 85:   return(0);
 86: }

 88: /*@C
 89:    MFNReasonView - Displays the reason an MFN solve converged or diverged.

 91:    Collective on MFN

 93:    Parameter:
 94: +  mfn - the matrix function context
 95: -  viewer - the viewer to display the reason

 97:    Options Database Keys:
 98: .  -mfn_converged_reason - print reason for convergence, and number of iterations

100:    Level: intermediate

102: .seealso: MFNSetTolerances(), MFNGetIterationNumber()
103: @*/
104: PetscErrorCode MFNReasonView(MFN mfn,PetscViewer viewer)
105: {
107:   PetscBool      isAscii;

110:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
111:   if (isAscii) {
112:     PetscViewerASCIIAddTab(viewer,((PetscObject)mfn)->tablevel);
113:     if (mfn->reason > 0) {
114:       PetscViewerASCIIPrintf(viewer,"%s Matrix function solve converged due to %s; iterations %D\n",((PetscObject)mfn)->prefix?((PetscObject)mfn)->prefix:"",MFNConvergedReasons[mfn->reason],mfn->its);
115:     } else {
116:       PetscViewerASCIIPrintf(viewer,"%s Matrix function solve did not converge due to %s; iterations %D\n",((PetscObject)mfn)->prefix?((PetscObject)mfn)->prefix:"",MFNConvergedReasons[mfn->reason],mfn->its);
117:     }
118:     PetscViewerASCIISubtractTab(viewer,((PetscObject)mfn)->tablevel);
119:   }
120:   return(0);
121: }

123: /*@
124:    MFNReasonViewFromOptions - Processes command line options to determine if/how
125:    the MFN converged reason is to be viewed.

127:    Collective on MFN

129:    Input Parameters:
130: .  mfn - the matrix function context

132:    Level: developer
133: @*/
134: PetscErrorCode MFNReasonViewFromOptions(MFN mfn)
135: {
136:   PetscErrorCode    ierr;
137:   PetscViewer       viewer;
138:   PetscBool         flg;
139:   static PetscBool  incall = PETSC_FALSE;
140:   PetscViewerFormat format;

143:   if (incall) return(0);
144:   incall = PETSC_TRUE;
145:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg);
146:   if (flg) {
147:     PetscViewerPushFormat(viewer,format);
148:     MFNReasonView(mfn,viewer);
149:     PetscViewerPopFormat(viewer);
150:     PetscViewerDestroy(&viewer);
151:   }
152:   incall = PETSC_FALSE;
153:   return(0);
154: }

156: /*@
157:    MFNCreate - Creates the default MFN context.

159:    Collective on MPI_Comm

161:    Input Parameter:
162: .  comm - MPI communicator

164:    Output Parameter:
165: .  mfn - location to put the MFN context

167:    Note:
168:    The default MFN type is MFNKRYLOV

170:    Level: beginner

172: .seealso: MFNSetUp(), MFNSolve(), MFNDestroy(), MFN
173: @*/
174: PetscErrorCode MFNCreate(MPI_Comm comm,MFN *outmfn)
175: {
177:   MFN            mfn;

181:   *outmfn = 0;
182:   MFNInitializePackage();
183:   SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView);

185:   mfn->A               = NULL;
186:   mfn->fn              = NULL;
187:   mfn->max_it          = 0;
188:   mfn->ncv             = 0;
189:   mfn->tol             = PETSC_DEFAULT;
190:   mfn->errorifnotconverged = PETSC_FALSE;

192:   mfn->numbermonitors  = 0;

194:   mfn->V               = NULL;
195:   mfn->nwork           = 0;
196:   mfn->work            = NULL;
197:   mfn->data            = NULL;

199:   mfn->its             = 0;
200:   mfn->nv              = 0;
201:   mfn->errest          = 0;
202:   mfn->setupcalled     = 0;
203:   mfn->reason          = MFN_CONVERGED_ITERATING;

205:   *outmfn = mfn;
206:   return(0);
207: }

209: /*@C
210:    MFNSetType - Selects the particular solver to be used in the MFN object.

212:    Logically Collective on MFN

214:    Input Parameters:
215: +  mfn  - the matrix function context
216: -  type - a known method

218:    Options Database Key:
219: .  -mfn_type <method> - Sets the method; use -help for a list
220:     of available methods

222:    Notes:
223:    See "slepc/include/slepcmfn.h" for available methods. The default
224:    is MFNKRYLOV

226:    Normally, it is best to use the MFNSetFromOptions() command and
227:    then set the MFN type from the options database rather than by using
228:    this routine.  Using the options database provides the user with
229:    maximum flexibility in evaluating the different available methods.
230:    The MFNSetType() routine is provided for those situations where it
231:    is necessary to set the iterative solver independently of the command
232:    line or options database.

234:    Level: intermediate

236: .seealso: MFNType
237: @*/
238: PetscErrorCode MFNSetType(MFN mfn,MFNType type)
239: {
240:   PetscErrorCode ierr,(*r)(MFN);
241:   PetscBool      match;


247:   PetscObjectTypeCompare((PetscObject)mfn,type,&match);
248:   if (match) return(0);

250:   PetscFunctionListFind(MFNList,type,&r);
251:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown MFN type given: %s",type);

253:   if (mfn->ops->destroy) { (*mfn->ops->destroy)(mfn); }
254:   PetscMemzero(mfn->ops,sizeof(struct _MFNOps));

256:   mfn->setupcalled = 0;
257:   PetscObjectChangeTypeName((PetscObject)mfn,type);
258:   (*r)(mfn);
259:   return(0);
260: }

262: /*@C
263:    MFNGetType - Gets the MFN type as a string from the MFN object.

265:    Not Collective

267:    Input Parameter:
268: .  mfn - the matrix function context

270:    Output Parameter:
271: .  name - name of MFN method

273:    Level: intermediate

275: .seealso: MFNSetType()
276: @*/
277: PetscErrorCode MFNGetType(MFN mfn,MFNType *type)
278: {
282:   *type = ((PetscObject)mfn)->type_name;
283:   return(0);
284: }

286: /*@C
287:    MFNRegister - Adds a method to the matrix function solver package.

289:    Not Collective

291:    Input Parameters:
292: +  name - name of a new user-defined solver
293: -  function - routine to create the solver context

295:    Notes:
296:    MFNRegister() may be called multiple times to add several user-defined solvers.

298:    Sample usage:
299: .vb
300:     MFNRegister("my_solver",MySolverCreate);
301: .ve

303:    Then, your solver can be chosen with the procedural interface via
304: $     MFNSetType(mfn,"my_solver")
305:    or at runtime via the option
306: $     -mfn_type my_solver

308:    Level: advanced

310: .seealso: MFNRegisterAll()
311: @*/
312: PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
313: {

317:   PetscFunctionListAdd(&MFNList,name,function);
318:   return(0);
319: }

321: /*@
322:    MFNReset - Resets the MFN context to the initial state (prior to setup)
323:    and destroys any allocated Vecs and Mats.

325:    Collective on MFN

327:    Input Parameter:
328: .  mfn - matrix function context obtained from MFNCreate()

330:    Level: advanced

332: .seealso: MFNDestroy()
333: @*/
334: PetscErrorCode MFNReset(MFN mfn)
335: {

340:   if (!mfn) return(0);
341:   if (mfn->ops->reset) { (mfn->ops->reset)(mfn); }
342:   MatDestroy(&mfn->A);
343:   BVDestroy(&mfn->V);
344:   VecDestroyVecs(mfn->nwork,&mfn->work);
345:   mfn->nwork = 0;
346:   mfn->setupcalled = 0;
347:   return(0);
348: }

350: /*@
351:    MFNDestroy - Destroys the MFN context.

353:    Collective on MFN

355:    Input Parameter:
356: .  mfn - matrix function context obtained from MFNCreate()

358:    Level: beginner

360: .seealso: MFNCreate(), MFNSetUp(), MFNSolve()
361: @*/
362: PetscErrorCode MFNDestroy(MFN *mfn)
363: {

367:   if (!*mfn) return(0);
369:   if (--((PetscObject)(*mfn))->refct > 0) { *mfn = 0; return(0); }
370:   MFNReset(*mfn);
371:   if ((*mfn)->ops->destroy) { (*(*mfn)->ops->destroy)(*mfn); }
372:   FNDestroy(&(*mfn)->fn);
373:   MFNMonitorCancel(*mfn);
374:   PetscHeaderDestroy(mfn);
375:   return(0);
376: }

378: /*@
379:    MFNSetBV - Associates a basis vectors object to the matrix function solver.

381:    Collective on MFN

383:    Input Parameters:
384: +  mfn - matrix function context obtained from MFNCreate()
385: -  bv  - the basis vectors object

387:    Note:
388:    Use MFNGetBV() to retrieve the basis vectors context (for example,
389:    to free it at the end of the computations).

391:    Level: advanced

393: .seealso: MFNGetBV()
394: @*/
395: PetscErrorCode MFNSetBV(MFN mfn,BV bv)
396: {

403:   PetscObjectReference((PetscObject)bv);
404:   BVDestroy(&mfn->V);
405:   mfn->V = bv;
406:   PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->V);
407:   return(0);
408: }

410: /*@
411:    MFNGetBV - Obtain the basis vectors object associated to the matrix
412:    function solver.

414:    Not Collective

416:    Input Parameters:
417: .  mfn - matrix function context obtained from MFNCreate()

419:    Output Parameter:
420: .  bv - basis vectors context

422:    Level: advanced

424: .seealso: MFNSetBV()
425: @*/
426: PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
427: {

433:   if (!mfn->V) {
434:     BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V);
435:     PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->V);
436:   }
437:   *bv = mfn->V;
438:   return(0);
439: }

441: /*@
442:    MFNSetFN - Specifies the function to be computed.

444:    Collective on MFN

446:    Input Parameters:
447: +  mfn - matrix function context obtained from MFNCreate()
448: -  fn  - the math function object

450:    Note:
451:    Use MFNGetFN() to retrieve the math function context (for example,
452:    to free it at the end of the computations).

454:    Level: beginner

456: .seealso: MFNGetFN()
457: @*/
458: PetscErrorCode MFNSetFN(MFN mfn,FN fn)
459: {

466:   PetscObjectReference((PetscObject)fn);
467:   FNDestroy(&mfn->fn);
468:   mfn->fn = fn;
469:   PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->fn);
470:   return(0);
471: }

473: /*@
474:    MFNGetFN - Obtain the math function object associated to the MFN object.

476:    Not Collective

478:    Input Parameters:
479: .  mfn - matrix function context obtained from MFNCreate()

481:    Output Parameter:
482: .  fn - math function context

484:    Level: beginner

486: .seealso: MFNSetFN()
487: @*/
488: PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
489: {

495:   if (!mfn->fn) {
496:     FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn);
497:     PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->fn);
498:   }
499:   *fn = mfn->fn;
500:   return(0);
501: }