Actual source code: mfnbasic.c
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: /*@
31: MFNView - Prints the MFN data structure.
33: Collective
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;
61: PetscFunctionBegin;
63: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mfn),&viewer));
65: PetscCheckSameComm(mfn,1,viewer,2);
67: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
68: if (isascii) {
69: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mfn,viewer));
70: PetscCall(PetscViewerASCIIPushTab(viewer));
71: PetscTryTypeMethod(mfn,view,viewer);
72: PetscCall(PetscViewerASCIIPopTab(viewer));
73: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",mfn->ncv));
74: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",mfn->max_it));
75: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)mfn->tol));
76: } else PetscTryTypeMethod(mfn,view,viewer);
77: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
78: if (!mfn->V) PetscCall(MFNGetFN(mfn,&mfn->fn));
79: PetscCall(FNView(mfn->fn,viewer));
80: if (!mfn->V) PetscCall(MFNGetBV(mfn,&mfn->V));
81: PetscCall(BVView(mfn->V,viewer));
82: PetscCall(PetscViewerPopFormat(viewer));
83: PetscFunctionReturn(PETSC_SUCCESS);
84: }
86: /*@
87: MFNViewFromOptions - View from options
89: Collective
91: Input Parameters:
92: + mfn - the matrix function context
93: . obj - optional object
94: - name - command line option
96: Level: intermediate
98: .seealso: `MFNView()`, `MFNCreate()`
99: @*/
100: PetscErrorCode MFNViewFromOptions(MFN mfn,PetscObject obj,const char name[])
101: {
102: PetscFunctionBegin;
104: PetscCall(PetscObjectViewFromOptions((PetscObject)mfn,obj,name));
105: PetscFunctionReturn(PETSC_SUCCESS);
106: }
107: /*@
108: MFNConvergedReasonView - Displays the reason an MFN solve converged or diverged.
110: Collective
112: Input Parameters:
113: + mfn - the matrix function context
114: - viewer - the viewer to display the reason
116: Options Database Keys:
117: . -mfn_converged_reason - print reason for convergence, and number of iterations
119: Note:
120: To change the format of the output call PetscViewerPushFormat(viewer,format) before
121: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
122: display a reason if it fails. The latter can be set in the command line with
123: -mfn_converged_reason ::failed
125: Level: intermediate
127: .seealso: `MFNSetTolerances()`, `MFNGetIterationNumber()`, `MFNConvergedReasonViewFromOptions()`
128: @*/
129: PetscErrorCode MFNConvergedReasonView(MFN mfn,PetscViewer viewer)
130: {
131: PetscBool isAscii;
132: PetscViewerFormat format;
134: PetscFunctionBegin;
135: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)mfn));
136: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
137: if (isAscii) {
138: PetscCall(PetscViewerGetFormat(viewer,&format));
139: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)mfn)->tablevel));
140: if (mfn->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(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));
141: else if (mfn->reason <= 0) PetscCall(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));
142: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)mfn)->tablevel));
143: }
144: PetscFunctionReturn(PETSC_SUCCESS);
145: }
147: /*@
148: MFNConvergedReasonViewFromOptions - Processes command line options to determine if/how
149: the MFN converged reason is to be viewed.
151: Collective
153: Input Parameter:
154: . mfn - the matrix function context
156: Level: developer
158: .seealso: `MFNConvergedReasonView()`
159: @*/
160: PetscErrorCode MFNConvergedReasonViewFromOptions(MFN mfn)
161: {
162: PetscViewer viewer;
163: PetscBool flg;
164: static PetscBool incall = PETSC_FALSE;
165: PetscViewerFormat format;
167: PetscFunctionBegin;
168: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
169: incall = PETSC_TRUE;
170: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->options,((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg));
171: if (flg) {
172: PetscCall(PetscViewerPushFormat(viewer,format));
173: PetscCall(MFNConvergedReasonView(mfn,viewer));
174: PetscCall(PetscViewerPopFormat(viewer));
175: PetscCall(PetscViewerDestroy(&viewer));
176: }
177: incall = PETSC_FALSE;
178: PetscFunctionReturn(PETSC_SUCCESS);
179: }
181: /*@
182: MFNCreate - Creates the default MFN context.
184: Collective
186: Input Parameter:
187: . comm - MPI communicator
189: Output Parameter:
190: . outmfn - location to put the MFN context
192: Note:
193: The default MFN type is MFNKRYLOV
195: Level: beginner
197: .seealso: `MFNSetUp()`, `MFNSolve()`, `MFNDestroy()`, `MFN`
198: @*/
199: PetscErrorCode MFNCreate(MPI_Comm comm,MFN *outmfn)
200: {
201: MFN mfn;
203: PetscFunctionBegin;
204: PetscAssertPointer(outmfn,2);
205: PetscCall(MFNInitializePackage());
206: PetscCall(SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView));
208: mfn->A = NULL;
209: mfn->fn = NULL;
210: mfn->max_it = PETSC_DETERMINE;
211: mfn->ncv = PETSC_DETERMINE;
212: mfn->tol = PETSC_DETERMINE;
213: mfn->errorifnotconverged = PETSC_FALSE;
215: mfn->numbermonitors = 0;
217: mfn->V = NULL;
218: mfn->nwork = 0;
219: mfn->work = NULL;
220: mfn->data = NULL;
222: mfn->its = 0;
223: mfn->nv = 0;
224: mfn->errest = 0;
225: mfn->setupcalled = 0;
226: mfn->reason = MFN_CONVERGED_ITERATING;
228: *outmfn = mfn;
229: PetscFunctionReturn(PETSC_SUCCESS);
230: }
232: /*@
233: MFNSetType - Selects the particular solver to be used in the MFN object.
235: Logically Collective
237: Input Parameters:
238: + mfn - the matrix function context
239: - type - a known method
241: Options Database Key:
242: . -mfn_type <method> - Sets the method; use -help for a list
243: of available methods
245: Notes:
246: See "slepc/include/slepcmfn.h" for available methods. The default
247: is MFNKRYLOV
249: Normally, it is best to use the MFNSetFromOptions() command and
250: then set the MFN type from the options database rather than by using
251: this routine. Using the options database provides the user with
252: maximum flexibility in evaluating the different available methods.
253: The MFNSetType() routine is provided for those situations where it
254: is necessary to set the iterative solver independently of the command
255: line or options database.
257: Level: intermediate
259: .seealso: `MFNType`
260: @*/
261: PetscErrorCode MFNSetType(MFN mfn,MFNType type)
262: {
263: PetscErrorCode (*r)(MFN);
264: PetscBool match;
266: PetscFunctionBegin;
268: PetscAssertPointer(type,2);
270: PetscCall(PetscObjectTypeCompare((PetscObject)mfn,type,&match));
271: if (match) PetscFunctionReturn(PETSC_SUCCESS);
273: PetscCall(PetscFunctionListFind(MFNList,type,&r));
274: PetscCheck(r,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown MFN type given: %s",type);
276: PetscTryTypeMethod(mfn,destroy);
277: PetscCall(PetscMemzero(mfn->ops,sizeof(struct _MFNOps)));
279: mfn->setupcalled = 0;
280: PetscCall(PetscObjectChangeTypeName((PetscObject)mfn,type));
281: PetscCall((*r)(mfn));
282: PetscFunctionReturn(PETSC_SUCCESS);
283: }
285: /*@
286: MFNGetType - Gets the MFN type as a string from the MFN object.
288: Not Collective
290: Input Parameter:
291: . mfn - the matrix function context
293: Output Parameter:
294: . type - name of MFN method
296: Level: intermediate
298: .seealso: `MFNSetType()`
299: @*/
300: PetscErrorCode MFNGetType(MFN mfn,MFNType *type)
301: {
302: PetscFunctionBegin;
304: PetscAssertPointer(type,2);
305: *type = ((PetscObject)mfn)->type_name;
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
309: /*@C
310: MFNRegister - Adds a method to the matrix function solver package.
312: Not Collective
314: Input Parameters:
315: + name - name of a new user-defined solver
316: - function - routine to create the solver context
318: Notes:
319: MFNRegister() may be called multiple times to add several user-defined solvers.
321: Example Usage:
322: .vb
323: MFNRegister("my_solver",MySolverCreate);
324: .ve
326: Then, your solver can be chosen with the procedural interface via
327: $ MFNSetType(mfn,"my_solver")
328: or at runtime via the option
329: $ -mfn_type my_solver
331: Level: advanced
333: .seealso: `MFNRegisterAll()`
334: @*/
335: PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
336: {
337: PetscFunctionBegin;
338: PetscCall(MFNInitializePackage());
339: PetscCall(PetscFunctionListAdd(&MFNList,name,function));
340: PetscFunctionReturn(PETSC_SUCCESS);
341: }
343: /*@C
344: MFNMonitorRegister - Registers an MFN monitor routine that may be accessed with MFNMonitorSetFromOptions().
346: Not Collective
348: Input Parameters:
349: + name - name of a new monitor routine
350: . vtype - a PetscViewerType for the output
351: . format - a PetscViewerFormat for the output
352: . monitor - monitor routine, see MFNMonitorRegisterFn
353: . create - creation routine, or NULL
354: - destroy - destruction routine, or NULL
356: Notes:
357: MFNMonitorRegister() may be called multiple times to add several user-defined monitors.
359: The calling sequence for the given function matches the calling sequence of MFNMonitorFn
360: functions passed to MFNMonitorSet() with the additional requirement that its final argument
361: be a PetscViewerAndFormat.
363: Example Usage:
364: .vb
365: MFNMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
366: .ve
368: Then, your monitor can be chosen with the procedural interface via
369: $ MFNMonitorSetFromOptions(mfn,"-mfn_monitor_my_monitor","my_monitor",NULL)
370: or at runtime via the option
371: $ -mfn_monitor_my_monitor
373: Level: advanced
375: .seealso: `MFNMonitorSet()`, `MFNMonitorRegisterAll()`
376: @*/
377: PetscErrorCode MFNMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,MFNMonitorRegisterFn *monitor,MFNMonitorRegisterCreateFn *create,MFNMonitorRegisterDestroyFn *destroy)
378: {
379: char key[PETSC_MAX_PATH_LEN];
381: PetscFunctionBegin;
382: PetscCall(MFNInitializePackage());
383: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
384: PetscCall(PetscFunctionListAdd(&MFNMonitorList,key,monitor));
385: if (create) PetscCall(PetscFunctionListAdd(&MFNMonitorCreateList,key,create));
386: if (destroy) PetscCall(PetscFunctionListAdd(&MFNMonitorDestroyList,key,destroy));
387: PetscFunctionReturn(PETSC_SUCCESS);
388: }
390: /*@
391: MFNReset - Resets the MFN context to the initial state (prior to setup)
392: and destroys any allocated Vecs and Mats.
394: Collective
396: Input Parameter:
397: . mfn - matrix function context obtained from MFNCreate()
399: Level: advanced
401: .seealso: `MFNDestroy()`
402: @*/
403: PetscErrorCode MFNReset(MFN mfn)
404: {
405: PetscFunctionBegin;
407: if (!mfn) PetscFunctionReturn(PETSC_SUCCESS);
408: PetscTryTypeMethod(mfn,reset);
409: PetscCall(MatDestroy(&mfn->A));
410: PetscCall(BVDestroy(&mfn->V));
411: PetscCall(VecDestroyVecs(mfn->nwork,&mfn->work));
412: mfn->nwork = 0;
413: mfn->setupcalled = 0;
414: PetscFunctionReturn(PETSC_SUCCESS);
415: }
417: /*@
418: MFNDestroy - Destroys the MFN context.
420: Collective
422: Input Parameter:
423: . mfn - matrix function context obtained from MFNCreate()
425: Level: beginner
427: .seealso: `MFNCreate()`, `MFNSetUp()`, `MFNSolve()`
428: @*/
429: PetscErrorCode MFNDestroy(MFN *mfn)
430: {
431: PetscFunctionBegin;
432: if (!*mfn) PetscFunctionReturn(PETSC_SUCCESS);
434: if (--((PetscObject)*mfn)->refct > 0) { *mfn = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
435: PetscCall(MFNReset(*mfn));
436: PetscTryTypeMethod(*mfn,destroy);
437: PetscCall(FNDestroy(&(*mfn)->fn));
438: PetscCall(MatDestroy(&(*mfn)->AT));
439: PetscCall(MFNMonitorCancel(*mfn));
440: PetscCall(PetscHeaderDestroy(mfn));
441: PetscFunctionReturn(PETSC_SUCCESS);
442: }
444: /*@
445: MFNSetBV - Associates a basis vectors object to the matrix function solver.
447: Collective
449: Input Parameters:
450: + mfn - matrix function context obtained from MFNCreate()
451: - bv - the basis vectors object
453: Note:
454: Use MFNGetBV() to retrieve the basis vectors context (for example,
455: to free it at the end of the computations).
457: Level: advanced
459: .seealso: `MFNGetBV()`
460: @*/
461: PetscErrorCode MFNSetBV(MFN mfn,BV bv)
462: {
463: PetscFunctionBegin;
466: PetscCheckSameComm(mfn,1,bv,2);
467: PetscCall(PetscObjectReference((PetscObject)bv));
468: PetscCall(BVDestroy(&mfn->V));
469: mfn->V = bv;
470: PetscFunctionReturn(PETSC_SUCCESS);
471: }
473: /*@
474: MFNGetBV - Obtain the basis vectors object associated to the matrix
475: function solver.
477: Not Collective
479: Input Parameters:
480: . mfn - matrix function context obtained from MFNCreate()
482: Output Parameter:
483: . bv - basis vectors context
485: Level: advanced
487: .seealso: `MFNSetBV()`
488: @*/
489: PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
490: {
491: PetscFunctionBegin;
493: PetscAssertPointer(bv,2);
494: if (!mfn->V) {
495: PetscCall(BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V));
496: PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->V,(PetscObject)mfn,0));
497: PetscCall(PetscObjectSetOptions((PetscObject)mfn->V,((PetscObject)mfn)->options));
498: }
499: *bv = mfn->V;
500: PetscFunctionReturn(PETSC_SUCCESS);
501: }
503: /*@
504: MFNSetFN - Specifies the function to be computed.
506: Collective
508: Input Parameters:
509: + mfn - matrix function context obtained from MFNCreate()
510: - fn - the math function object
512: Note:
513: Use MFNGetFN() to retrieve the math function context (for example,
514: to free it at the end of the computations).
516: Level: beginner
518: .seealso: `MFNGetFN()`
519: @*/
520: PetscErrorCode MFNSetFN(MFN mfn,FN fn)
521: {
522: PetscFunctionBegin;
525: PetscCheckSameComm(mfn,1,fn,2);
526: PetscCall(PetscObjectReference((PetscObject)fn));
527: PetscCall(FNDestroy(&mfn->fn));
528: mfn->fn = fn;
529: PetscFunctionReturn(PETSC_SUCCESS);
530: }
532: /*@
533: MFNGetFN - Obtain the math function object associated to the MFN object.
535: Not Collective
537: Input Parameters:
538: . mfn - matrix function context obtained from MFNCreate()
540: Output Parameter:
541: . fn - math function context
543: Level: beginner
545: .seealso: `MFNSetFN()`
546: @*/
547: PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
548: {
549: PetscFunctionBegin;
551: PetscAssertPointer(fn,2);
552: if (!mfn->fn) {
553: PetscCall(FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn));
554: PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->fn,(PetscObject)mfn,0));
555: PetscCall(PetscObjectSetOptions((PetscObject)mfn->fn,((PetscObject)mfn)->options));
556: }
557: *fn = mfn->fn;
558: PetscFunctionReturn(PETSC_SUCCESS);
559: }