Actual source code: mfnbasic.c

slepc-3.18.2 2023-01-26
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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>

 16: /* Logging support */
 17: PetscClassId      MFN_CLASSID = 0;
 18: PetscLogEvent     MFN_SetUp = 0,MFN_Solve = 0;

 20: /* List of registered MFN routines */
 21: PetscFunctionList MFNList = NULL;
 22: PetscBool         MFNRegisterAllCalled = PETSC_FALSE;

 24: /* List of registered MFN monitors */
 25: PetscFunctionList MFNMonitorList              = NULL;
 26: PetscFunctionList MFNMonitorCreateList        = NULL;
 27: PetscFunctionList MFNMonitorDestroyList       = NULL;
 28: PetscBool         MFNMonitorRegisterAllCalled = PETSC_FALSE;

 30: /*@C
 31:    MFNView - Prints the MFN data structure.

 33:    Collective on mfn

 35:    Input Parameters:
 36: +  mfn - the matrix function solver context
 37: -  viewer - optional visualization context

 39:    Options Database Key:
 40: .  -mfn_view -  Calls MFNView() at end of MFNSolve()

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

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

 53:    Level: beginner

 55: .seealso: MFNCreate()
 56: @*/
 57: PetscErrorCode MFNView(MFN mfn,PetscViewer viewer)
 58: {
 59:   PetscBool      isascii;

 62:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mfn),&viewer);

 66:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 67:   if (isascii) {
 68:     PetscObjectPrintClassNamePrefixType((PetscObject)mfn,viewer);
 69:     PetscViewerASCIIPushTab(viewer);
 70:     PetscTryTypeMethod(mfn,view,viewer);
 71:     PetscViewerASCIIPopTab(viewer);
 72:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %" PetscInt_FMT "\n",mfn->ncv);
 73:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %" PetscInt_FMT "\n",mfn->max_it);
 74:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)mfn->tol);
 75:   } else PetscTryTypeMethod(mfn,view,viewer);
 76:   PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
 77:   if (!mfn->V) MFNGetFN(mfn,&mfn->fn);
 78:   FNView(mfn->fn,viewer);
 79:   if (!mfn->V) MFNGetBV(mfn,&mfn->V);
 80:   BVView(mfn->V,viewer);
 81:   PetscViewerPopFormat(viewer);
 82:   return 0;
 83: }

 85: /*@C
 86:    MFNViewFromOptions - View from options

 88:    Collective on MFN

 90:    Input Parameters:
 91: +  mfn  - the matrix function context
 92: .  obj  - optional object
 93: -  name - command line option

 95:    Level: intermediate

 97: .seealso: MFNView(), MFNCreate()
 98: @*/
 99: PetscErrorCode MFNViewFromOptions(MFN mfn,PetscObject obj,const char name[])
100: {
102:   PetscObjectViewFromOptions((PetscObject)mfn,obj,name);
103:   return 0;
104: }
105: /*@C
106:    MFNConvergedReasonView - Displays the reason an MFN solve converged or diverged.

108:    Collective on mfn

110:    Input Parameters:
111: +  mfn - the matrix function context
112: -  viewer - the viewer to display the reason

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

117:    Note:
118:    To change the format of the output call PetscViewerPushFormat(viewer,format) before
119:    this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
120:    display a reason if it fails. The latter can be set in the command line with
121:    -mfn_converged_reason ::failed

123:    Level: intermediate

125: .seealso: MFNSetTolerances(), MFNGetIterationNumber(), MFNConvergedReasonViewFromOptions()
126: @*/
127: PetscErrorCode MFNConvergedReasonView(MFN mfn,PetscViewer viewer)
128: {
129:   PetscBool         isAscii;
130:   PetscViewerFormat format;

132:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)mfn));
133:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
134:   if (isAscii) {
135:     PetscViewerGetFormat(viewer,&format);
136:     PetscViewerASCIIAddTab(viewer,((PetscObject)mfn)->tablevel);
137:     if (mfn->reason > 0 && format != PETSC_VIEWER_FAILED) PetscViewerASCIIPrintf(viewer,"%s Matrix function solve converged due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)mfn)->prefix?((PetscObject)mfn)->prefix:"",MFNConvergedReasons[mfn->reason],mfn->its);
138:     else if (mfn->reason <= 0) PetscViewerASCIIPrintf(viewer,"%s Matrix function solve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)mfn)->prefix?((PetscObject)mfn)->prefix:"",MFNConvergedReasons[mfn->reason],mfn->its);
139:     PetscViewerASCIISubtractTab(viewer,((PetscObject)mfn)->tablevel);
140:   }
141:   return 0;
142: }

144: /*@
145:    MFNConvergedReasonViewFromOptions - Processes command line options to determine if/how
146:    the MFN converged reason is to be viewed.

148:    Collective on mfn

150:    Input Parameter:
151: .  mfn - the matrix function context

153:    Level: developer

155: .seealso: MFNConvergedReasonView()
156: @*/
157: PetscErrorCode MFNConvergedReasonViewFromOptions(MFN mfn)
158: {
159:   PetscViewer       viewer;
160:   PetscBool         flg;
161:   static PetscBool  incall = PETSC_FALSE;
162:   PetscViewerFormat format;

164:   if (incall) return 0;
165:   incall = PETSC_TRUE;
166:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->options,((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg);
167:   if (flg) {
168:     PetscViewerPushFormat(viewer,format);
169:     MFNConvergedReasonView(mfn,viewer);
170:     PetscViewerPopFormat(viewer);
171:     PetscViewerDestroy(&viewer);
172:   }
173:   incall = PETSC_FALSE;
174:   return 0;
175: }

177: /*@
178:    MFNCreate - Creates the default MFN context.

180:    Collective

182:    Input Parameter:
183: .  comm - MPI communicator

185:    Output Parameter:
186: .  outmfn - location to put the MFN context

188:    Note:
189:    The default MFN type is MFNKRYLOV

191:    Level: beginner

193: .seealso: MFNSetUp(), MFNSolve(), MFNDestroy(), MFN
194: @*/
195: PetscErrorCode MFNCreate(MPI_Comm comm,MFN *outmfn)
196: {
197:   MFN            mfn;

200:   *outmfn = NULL;
201:   MFNInitializePackage();
202:   SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView);

204:   mfn->A               = NULL;
205:   mfn->fn              = NULL;
206:   mfn->max_it          = PETSC_DEFAULT;
207:   mfn->ncv             = PETSC_DEFAULT;
208:   mfn->tol             = PETSC_DEFAULT;
209:   mfn->errorifnotconverged = PETSC_FALSE;

211:   mfn->numbermonitors  = 0;

213:   mfn->V               = NULL;
214:   mfn->nwork           = 0;
215:   mfn->work            = NULL;
216:   mfn->data            = NULL;

218:   mfn->its             = 0;
219:   mfn->nv              = 0;
220:   mfn->errest          = 0;
221:   mfn->setupcalled     = 0;
222:   mfn->reason          = MFN_CONVERGED_ITERATING;

224:   *outmfn = mfn;
225:   return 0;
226: }

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

231:    Logically Collective on mfn

233:    Input Parameters:
234: +  mfn  - the matrix function context
235: -  type - a known method

237:    Options Database Key:
238: .  -mfn_type <method> - Sets the method; use -help for a list
239:     of available methods

241:    Notes:
242:    See "slepc/include/slepcmfn.h" for available methods. The default
243:    is MFNKRYLOV

245:    Normally, it is best to use the MFNSetFromOptions() command and
246:    then set the MFN type from the options database rather than by using
247:    this routine.  Using the options database provides the user with
248:    maximum flexibility in evaluating the different available methods.
249:    The MFNSetType() routine is provided for those situations where it
250:    is necessary to set the iterative solver independently of the command
251:    line or options database.

253:    Level: intermediate

255: .seealso: MFNType
256: @*/
257: PetscErrorCode MFNSetType(MFN mfn,MFNType type)
258: {
259:   PetscErrorCode (*r)(MFN);
260:   PetscBool      match;


265:   PetscObjectTypeCompare((PetscObject)mfn,type,&match);
266:   if (match) return 0;

268:   PetscFunctionListFind(MFNList,type,&r);

271:   PetscTryTypeMethod(mfn,destroy);
272:   PetscMemzero(mfn->ops,sizeof(struct _MFNOps));

274:   mfn->setupcalled = 0;
275:   PetscObjectChangeTypeName((PetscObject)mfn,type);
276:   (*r)(mfn);
277:   return 0;
278: }

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

283:    Not Collective

285:    Input Parameter:
286: .  mfn - the matrix function context

288:    Output Parameter:
289: .  type - name of MFN method

291:    Level: intermediate

293: .seealso: MFNSetType()
294: @*/
295: PetscErrorCode MFNGetType(MFN mfn,MFNType *type)
296: {
299:   *type = ((PetscObject)mfn)->type_name;
300:   return 0;
301: }

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

306:    Not Collective

308:    Input Parameters:
309: +  name - name of a new user-defined solver
310: -  function - routine to create the solver context

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

315:    Sample usage:
316: .vb
317:     MFNRegister("my_solver",MySolverCreate);
318: .ve

320:    Then, your solver can be chosen with the procedural interface via
321: $     MFNSetType(mfn,"my_solver")
322:    or at runtime via the option
323: $     -mfn_type my_solver

325:    Level: advanced

327: .seealso: MFNRegisterAll()
328: @*/
329: PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
330: {
331:   MFNInitializePackage();
332:   PetscFunctionListAdd(&MFNList,name,function);
333:   return 0;
334: }

336: /*@C
337:    MFNMonitorRegister - Adds MFN monitor routine.

339:    Not Collective

341:    Input Parameters:
342: +  name    - name of a new monitor routine
343: .  vtype   - a PetscViewerType for the output
344: .  format  - a PetscViewerFormat for the output
345: .  monitor - monitor routine
346: .  create  - creation routine, or NULL
347: -  destroy - destruction routine, or NULL

349:    Notes:
350:    MFNMonitorRegister() may be called multiple times to add several user-defined monitors.

352:    Sample usage:
353: .vb
354:    MFNMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
355: .ve

357:    Then, your monitor can be chosen with the procedural interface via
358: $      MFNMonitorSetFromOptions(mfn,"-mfn_monitor_my_monitor","my_monitor",NULL)
359:    or at runtime via the option
360: $      -mfn_monitor_my_monitor

362:    Level: advanced

364: .seealso: MFNMonitorRegisterAll()
365: @*/
366: PetscErrorCode MFNMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(MFN,PetscInt,PetscReal,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
367: {
368:   char           key[PETSC_MAX_PATH_LEN];

370:   MFNInitializePackage();
371:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
372:   PetscFunctionListAdd(&MFNMonitorList,key,monitor);
373:   if (create)  PetscFunctionListAdd(&MFNMonitorCreateList,key,create);
374:   if (destroy) PetscFunctionListAdd(&MFNMonitorDestroyList,key,destroy);
375:   return 0;
376: }

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

382:    Collective on mfn

384:    Input Parameter:
385: .  mfn - matrix function context obtained from MFNCreate()

387:    Level: advanced

389: .seealso: MFNDestroy()
390: @*/
391: PetscErrorCode MFNReset(MFN mfn)
392: {
394:   if (!mfn) return 0;
395:   PetscTryTypeMethod(mfn,reset);
396:   MatDestroy(&mfn->A);
397:   BVDestroy(&mfn->V);
398:   VecDestroyVecs(mfn->nwork,&mfn->work);
399:   mfn->nwork = 0;
400:   mfn->setupcalled = 0;
401:   return 0;
402: }

404: /*@C
405:    MFNDestroy - Destroys the MFN context.

407:    Collective on mfn

409:    Input Parameter:
410: .  mfn - matrix function context obtained from MFNCreate()

412:    Level: beginner

414: .seealso: MFNCreate(), MFNSetUp(), MFNSolve()
415: @*/
416: PetscErrorCode MFNDestroy(MFN *mfn)
417: {
418:   if (!*mfn) return 0;
420:   if (--((PetscObject)(*mfn))->refct > 0) { *mfn = NULL; return 0; }
421:   MFNReset(*mfn);
422:   PetscTryTypeMethod(*mfn,destroy);
423:   FNDestroy(&(*mfn)->fn);
424:   MatDestroy(&(*mfn)->AT);
425:   MFNMonitorCancel(*mfn);
426:   PetscHeaderDestroy(mfn);
427:   return 0;
428: }

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

433:    Collective on mfn

435:    Input Parameters:
436: +  mfn - matrix function context obtained from MFNCreate()
437: -  bv  - the basis vectors object

439:    Note:
440:    Use MFNGetBV() to retrieve the basis vectors context (for example,
441:    to free it at the end of the computations).

443:    Level: advanced

445: .seealso: MFNGetBV()
446: @*/
447: PetscErrorCode MFNSetBV(MFN mfn,BV bv)
448: {
452:   PetscObjectReference((PetscObject)bv);
453:   BVDestroy(&mfn->V);
454:   mfn->V = bv;
455:   return 0;
456: }

458: /*@
459:    MFNGetBV - Obtain the basis vectors object associated to the matrix
460:    function solver.

462:    Not Collective

464:    Input Parameters:
465: .  mfn - matrix function context obtained from MFNCreate()

467:    Output Parameter:
468: .  bv - basis vectors context

470:    Level: advanced

472: .seealso: MFNSetBV()
473: @*/
474: PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
475: {
478:   if (!mfn->V) {
479:     BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V);
480:     PetscObjectIncrementTabLevel((PetscObject)mfn->V,(PetscObject)mfn,0);
481:     PetscObjectSetOptions((PetscObject)mfn->V,((PetscObject)mfn)->options);
482:   }
483:   *bv = mfn->V;
484:   return 0;
485: }

487: /*@
488:    MFNSetFN - Specifies the function to be computed.

490:    Collective on mfn

492:    Input Parameters:
493: +  mfn - matrix function context obtained from MFNCreate()
494: -  fn  - the math function object

496:    Note:
497:    Use MFNGetFN() to retrieve the math function context (for example,
498:    to free it at the end of the computations).

500:    Level: beginner

502: .seealso: MFNGetFN()
503: @*/
504: PetscErrorCode MFNSetFN(MFN mfn,FN fn)
505: {
509:   PetscObjectReference((PetscObject)fn);
510:   FNDestroy(&mfn->fn);
511:   mfn->fn = fn;
512:   return 0;
513: }

515: /*@
516:    MFNGetFN - Obtain the math function object associated to the MFN object.

518:    Not Collective

520:    Input Parameters:
521: .  mfn - matrix function context obtained from MFNCreate()

523:    Output Parameter:
524: .  fn - math function context

526:    Level: beginner

528: .seealso: MFNSetFN()
529: @*/
530: PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
531: {
534:   if (!mfn->fn) {
535:     FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn);
536:     PetscObjectIncrementTabLevel((PetscObject)mfn->fn,(PetscObject)mfn,0);
537:     PetscObjectSetOptions((PetscObject)mfn->fn,((PetscObject)mfn)->options);
538:   }
539:   *fn = mfn->fn;
540:   return 0;
541: }