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: Notes:
43: The available visualization contexts include
44: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
45: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
46: first process opens the file; all other processes send their data to the
47: first one to print
49: The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
50: to output to a specified file.
52: Level: beginner
54: .seealso: [](ch:mfn), `MFNCreate()`, `MFNViewFromOptions()`
55: @*/
56: PetscErrorCode MFNView(MFN mfn,PetscViewer viewer)
57: {
58: PetscBool isascii;
60: PetscFunctionBegin;
62: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mfn),&viewer));
64: PetscCheckSameComm(mfn,1,viewer,2);
66: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
67: if (isascii) {
68: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mfn,viewer));
69: PetscCall(PetscViewerASCIIPushTab(viewer));
70: PetscTryTypeMethod(mfn,view,viewer);
71: PetscCall(PetscViewerASCIIPopTab(viewer));
72: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",mfn->ncv));
73: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",mfn->max_it));
74: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)mfn->tol));
75: } else PetscTryTypeMethod(mfn,view,viewer);
76: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
77: if (!mfn->V) PetscCall(MFNGetFN(mfn,&mfn->fn));
78: PetscCall(FNView(mfn->fn,viewer));
79: if (!mfn->V) PetscCall(MFNGetBV(mfn,&mfn->V));
80: PetscCall(BVView(mfn->V,viewer));
81: PetscCall(PetscViewerPopFormat(viewer));
82: PetscFunctionReturn(PETSC_SUCCESS);
83: }
85: /*@
86: MFNViewFromOptions - View (print) an `MFN` object based on values in the options database.
88: Collective
90: Input Parameters:
91: + mfn - the matrix function solver context
92: . obj - optional object that provides the options prefix used to query the options database
93: - name - command line option
95: Level: intermediate
97: .seealso: [](ch:mfn), `MFNView()`, `MFNCreate()`, `PetscObjectViewFromOptions()`
98: @*/
99: PetscErrorCode MFNViewFromOptions(MFN mfn,PetscObject obj,const char name[])
100: {
101: PetscFunctionBegin;
103: PetscCall(PetscObjectViewFromOptions((PetscObject)mfn,obj,name));
104: PetscFunctionReturn(PETSC_SUCCESS);
105: }
106: /*@
107: MFNConvergedReasonView - Displays the reason an `MFN` solve converged or diverged.
109: Collective
111: Input Parameters:
112: + mfn - the matrix function solver context
113: - viewer - the viewer to display the reason
115: Options Database Key:
116: . -mfn_converged_reason - print reason for convergence/divergence, and number of iterations
118: Notes:
119: Use `MFNConvergedReasonViewFromOptions()` to display the reason based on values
120: in the options database.
122: To change the format of the output call `PetscViewerPushFormat()` before this
123: call. Use `PETSC_VIEWER_DEFAULT` for the default, or `PETSC_VIEWER_FAILED` to only
124: display a reason if it fails. The latter can be set in the command line with
125: `-mfn_converged_reason ::failed`.
127: Level: intermediate
129: .seealso: [](ch:mfn), `MFNSetTolerances()`, `MFNGetIterationNumber()`, `MFNConvergedReasonViewFromOptions()`
130: @*/
131: PetscErrorCode MFNConvergedReasonView(MFN mfn,PetscViewer viewer)
132: {
133: PetscBool isAscii;
134: PetscViewerFormat format;
136: PetscFunctionBegin;
137: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)mfn));
138: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
139: if (isAscii) {
140: PetscCall(PetscViewerGetFormat(viewer,&format));
141: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)mfn)->tablevel));
142: 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));
143: 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));
144: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)mfn)->tablevel));
145: }
146: PetscFunctionReturn(PETSC_SUCCESS);
147: }
149: /*@
150: MFNConvergedReasonViewFromOptions - Processes command line options to determine if/how
151: the `MFN` converged reason is to be viewed.
153: Collective
155: Input Parameter:
156: . mfn - the matrix function solver context
158: Level: intermediate
160: .seealso: [](ch:mfn), `MFNConvergedReasonView()`
161: @*/
162: PetscErrorCode MFNConvergedReasonViewFromOptions(MFN mfn)
163: {
164: PetscViewer viewer;
165: PetscBool flg;
166: static PetscBool incall = PETSC_FALSE;
167: PetscViewerFormat format;
169: PetscFunctionBegin;
170: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
171: incall = PETSC_TRUE;
172: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->options,((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg));
173: if (flg) {
174: PetscCall(PetscViewerPushFormat(viewer,format));
175: PetscCall(MFNConvergedReasonView(mfn,viewer));
176: PetscCall(PetscViewerPopFormat(viewer));
177: PetscCall(PetscViewerDestroy(&viewer));
178: }
179: incall = PETSC_FALSE;
180: PetscFunctionReturn(PETSC_SUCCESS);
181: }
183: /*@
184: MFNCreate - Creates the `MFN` context.
186: Collective
188: Input Parameter:
189: . comm - MPI communicator
191: Output Parameter:
192: . outmfn - location to put the `MFN` context
194: Note:
195: The default `MFN` type is `MFNKRYLOV`.
197: Level: beginner
199: .seealso: [](ch:mfn), `MFNSetUp()`, `MFNSolve()`, `MFNDestroy()`, `MFN`
200: @*/
201: PetscErrorCode MFNCreate(MPI_Comm comm,MFN *outmfn)
202: {
203: MFN mfn;
205: PetscFunctionBegin;
206: PetscAssertPointer(outmfn,2);
207: PetscCall(MFNInitializePackage());
208: PetscCall(SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView));
210: mfn->A = NULL;
211: mfn->fn = NULL;
212: mfn->max_it = PETSC_DETERMINE;
213: mfn->ncv = PETSC_DETERMINE;
214: mfn->tol = PETSC_DETERMINE;
215: mfn->errorifnotconverged = PETSC_FALSE;
217: mfn->numbermonitors = 0;
219: mfn->V = NULL;
220: mfn->nwork = 0;
221: mfn->work = NULL;
222: mfn->data = NULL;
224: mfn->its = 0;
225: mfn->nv = 0;
226: mfn->errest = 0;
227: mfn->setupcalled = 0;
228: mfn->reason = MFN_CONVERGED_ITERATING;
230: *outmfn = mfn;
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /*@
235: MFNSetType - Selects the particular solver to be used in the `MFN` object.
237: Logically Collective
239: Input Parameters:
240: + mfn - the matrix function solver context
241: - type - a known method
243: Options Database Key:
244: . -mfn_type \<type\> - sets the method; use `-help` for a list of available methods
246: Notes:
247: See `MFNType` for available methods. The default 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: [](ch:mfn), `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 solver context
293: Output Parameter:
294: . type - name of `MFN` method
296: Level: intermediate
298: .seealso: [](ch:mfn), `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: Note:
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: .vb
328: MFNSetType(mfn,"my_solver")
329: .ve
330: or at runtime via the option `-mfn_type my_solver`.
332: Level: advanced
334: .seealso: [](ch:mfn), `MFNRegisterAll()`
335: @*/
336: PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
337: {
338: PetscFunctionBegin;
339: PetscCall(MFNInitializePackage());
340: PetscCall(PetscFunctionListAdd(&MFNList,name,function));
341: PetscFunctionReturn(PETSC_SUCCESS);
342: }
344: /*@C
345: MFNMonitorRegister - Registers an `MFN` monitor routine that may be accessed with
346: `MFNMonitorSetFromOptions()`.
348: Not Collective
350: Input Parameters:
351: + name - name of a new monitor routine
352: . vtype - a `PetscViewerType` for the output
353: . format - a `PetscViewerFormat` for the output
354: . monitor - monitor routine, see `MFNMonitorRegisterFn`
355: . create - creation routine, or `NULL`
356: - destroy - destruction routine, or `NULL`
358: Notes:
359: `MFNMonitorRegister()` may be called multiple times to add several user-defined monitors.
361: The calling sequence for the given function matches the calling sequence of `MFNMonitorFn`
362: functions passed to `MFNMonitorSet()` with the additional requirement that its final argument
363: be a `PetscViewerAndFormat`.
365: Example Usage:
366: .vb
367: MFNMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
368: .ve
370: Then, your monitor can be chosen with the procedural interface via
371: .vb
372: MFNMonitorSetFromOptions(mfn,"-mfn_monitor_my_monitor","my_monitor",NULL);
373: .ve
374: or at runtime via the option `-mfn_monitor_my_monitor`.
376: Level: advanced
378: .seealso: [](ch:mfn), `MFNMonitorSet()`, `MFNMonitorRegisterAll()`, `MFNMonitorSetFromOptions()`
379: @*/
380: PetscErrorCode MFNMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,MFNMonitorRegisterFn *monitor,MFNMonitorRegisterCreateFn *create,MFNMonitorRegisterDestroyFn *destroy)
381: {
382: char key[PETSC_MAX_PATH_LEN];
384: PetscFunctionBegin;
385: PetscCall(MFNInitializePackage());
386: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
387: PetscCall(PetscFunctionListAdd(&MFNMonitorList,key,monitor));
388: if (create) PetscCall(PetscFunctionListAdd(&MFNMonitorCreateList,key,create));
389: if (destroy) PetscCall(PetscFunctionListAdd(&MFNMonitorDestroyList,key,destroy));
390: PetscFunctionReturn(PETSC_SUCCESS);
391: }
393: /*@
394: MFNReset - Resets the `MFN` context to the initial state (prior to setup)
395: and destroys any allocated `Vec`s and `Mat`s.
397: Collective
399: Input Parameter:
400: . mfn - the matrix function solver context
402: Level: advanced
404: .seealso: [](ch:mfn), `MFNDestroy()`
405: @*/
406: PetscErrorCode MFNReset(MFN mfn)
407: {
408: PetscFunctionBegin;
410: if (!mfn) PetscFunctionReturn(PETSC_SUCCESS);
411: PetscTryTypeMethod(mfn,reset);
412: PetscCall(MatDestroy(&mfn->A));
413: PetscCall(BVDestroy(&mfn->V));
414: PetscCall(VecDestroyVecs(mfn->nwork,&mfn->work));
415: mfn->nwork = 0;
416: mfn->setupcalled = 0;
417: PetscFunctionReturn(PETSC_SUCCESS);
418: }
420: /*@
421: MFNDestroy - Destroys the `MFN` context.
423: Collective
425: Input Parameter:
426: . mfn - the matrix function solver context
428: Level: beginner
430: .seealso: [](ch:mfn), `MFNCreate()`, `MFNSetUp()`, `MFNSolve()`
431: @*/
432: PetscErrorCode MFNDestroy(MFN *mfn)
433: {
434: PetscFunctionBegin;
435: if (!*mfn) PetscFunctionReturn(PETSC_SUCCESS);
437: if (--((PetscObject)*mfn)->refct > 0) { *mfn = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
438: PetscCall(MFNReset(*mfn));
439: PetscTryTypeMethod(*mfn,destroy);
440: PetscCall(FNDestroy(&(*mfn)->fn));
441: PetscCall(MatDestroy(&(*mfn)->AT));
442: PetscCall(MFNMonitorCancel(*mfn));
443: PetscCall(PetscHeaderDestroy(mfn));
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: /*@
448: MFNSetBV - Associates a basis vectors object to the matrix function solver.
450: Collective
452: Input Parameters:
453: + mfn - the matrix function solver context
454: - bv - the basis vectors object
456: Note:
457: Use `MFNGetBV()` to retrieve the basis vectors context (for example,
458: to free it at the end of the computations).
460: Level: advanced
462: .seealso: [](ch:mfn), `MFNGetBV()`
463: @*/
464: PetscErrorCode MFNSetBV(MFN mfn,BV bv)
465: {
466: PetscFunctionBegin;
469: PetscCheckSameComm(mfn,1,bv,2);
470: PetscCall(PetscObjectReference((PetscObject)bv));
471: PetscCall(BVDestroy(&mfn->V));
472: mfn->V = bv;
473: PetscFunctionReturn(PETSC_SUCCESS);
474: }
476: /*@
477: MFNGetBV - Obtain the basis vectors object associated to the matrix
478: function solver.
480: Not Collective
482: Input Parameter:
483: . mfn - the matrix function solver context
485: Output Parameter:
486: . bv - basis vectors context
488: Level: advanced
490: .seealso: [](ch:mfn), `MFNSetBV()`
491: @*/
492: PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
493: {
494: PetscFunctionBegin;
496: PetscAssertPointer(bv,2);
497: if (!mfn->V) {
498: PetscCall(BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V));
499: PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->V,(PetscObject)mfn,0));
500: PetscCall(PetscObjectSetOptions((PetscObject)mfn->V,((PetscObject)mfn)->options));
501: }
502: *bv = mfn->V;
503: PetscFunctionReturn(PETSC_SUCCESS);
504: }
506: /*@
507: MFNSetFN - Specifies the function to be computed.
509: Collective
511: Input Parameters:
512: + mfn - the matrix function solver context
513: - fn - the math function object
515: Notes:
516: At a later time, use `MFNGetFN()` to retrieve the math function context
517: (for example, to free it at the end of the computations).
519: This function is not called in normal usage. Instead, it is easier to
520: extract the internal `FN` object with `MFNGetFN()` and modify it.
522: Level: beginner
524: .seealso: [](ch:mfn), `MFNGetFN()`
525: @*/
526: PetscErrorCode MFNSetFN(MFN mfn,FN fn)
527: {
528: PetscFunctionBegin;
531: PetscCheckSameComm(mfn,1,fn,2);
532: PetscCall(PetscObjectReference((PetscObject)fn));
533: PetscCall(FNDestroy(&mfn->fn));
534: mfn->fn = fn;
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: /*@
539: MFNGetFN - Obtain the math function object associated to the `MFN` object.
541: Not Collective
543: Input Parameter:
544: . mfn - the matrix function solver context
546: Output Parameter:
547: . fn - math function context
549: Note:
550: This is the usual way to specify the function that needs to be applied
551: to a given vector in `MFNSolve()`.
553: Level: beginner
555: .seealso: [](ch:mfn), `MFNSetFN()`. `MFNSolve()`
556: @*/
557: PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
558: {
559: PetscFunctionBegin;
561: PetscAssertPointer(fn,2);
562: if (!mfn->fn) {
563: PetscCall(FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn));
564: PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->fn,(PetscObject)mfn,0));
565: PetscCall(PetscObjectSetOptions((PetscObject)mfn->fn,((PetscObject)mfn)->options));
566: }
567: *fn = mfn->fn;
568: PetscFunctionReturn(PETSC_SUCCESS);
569: }