Actual source code: mfnbasic.c
slepc-3.14.0 2020-09-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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: 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
45: @*/
46: PetscErrorCode MFNView(MFN mfn,PetscViewer viewer)
47: {
49: PetscBool isascii;
53: if (!viewer) {
54: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mfn),&viewer);
55: }
59: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
60: if (isascii) {
61: PetscObjectPrintClassNamePrefixType((PetscObject)mfn,viewer);
62: if (mfn->ops->view) {
63: PetscViewerASCIIPushTab(viewer);
64: (*mfn->ops->view)(mfn,viewer);
65: PetscViewerASCIIPopTab(viewer);
66: }
67: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",mfn->ncv);
68: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",mfn->max_it);
69: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)mfn->tol);
70: } else {
71: if (mfn->ops->view) {
72: (*mfn->ops->view)(mfn,viewer);
73: }
74: }
75: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
76: if (!mfn->V) { MFNGetFN(mfn,&mfn->fn); }
77: FNView(mfn->fn,viewer);
78: if (!mfn->V) { MFNGetBV(mfn,&mfn->V); }
79: BVView(mfn->V,viewer);
80: PetscViewerPopFormat(viewer);
81: return(0);
82: }
84: /*@C
85: MFNViewFromOptions - View from options
87: Collective on MFN
89: Input Parameters:
90: + mfn - the matrix function context
91: . obj - optional object
92: - name - command line option
94: Level: intermediate
96: .seealso: MFNView(), MFNCreate()
97: @*/
98: PetscErrorCode MFNViewFromOptions(MFN mfn,PetscObject obj,const char name[])
99: {
104: PetscObjectViewFromOptions((PetscObject)mfn,obj,name);
105: return(0);
106: }
107: /*@C
108: MFNConvergedReasonView - Displays the reason an MFN solve converged or diverged.
110: Collective on mfn
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: PetscErrorCode ierr;
132: PetscBool isAscii;
133: PetscViewerFormat format;
136: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)mfn));
137: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
138: if (isAscii) {
139: PetscViewerGetFormat(viewer,&format);
140: PetscViewerASCIIAddTab(viewer,((PetscObject)mfn)->tablevel);
141: if (mfn->reason > 0 && format != PETSC_VIEWER_FAILED) {
142: PetscViewerASCIIPrintf(viewer,"%s Matrix function solve converged due to %s; iterations %D\n",((PetscObject)mfn)->prefix?((PetscObject)mfn)->prefix:"",MFNConvergedReasons[mfn->reason],mfn->its);
143: } else if (mfn->reason <= 0) {
144: 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);
145: }
146: PetscViewerASCIISubtractTab(viewer,((PetscObject)mfn)->tablevel);
147: }
148: return(0);
149: }
151: /*@
152: MFNConvergedReasonViewFromOptions - Processes command line options to determine if/how
153: the MFN converged reason is to be viewed.
155: Collective on mfn
157: Input Parameter:
158: . mfn - the matrix function context
160: Level: developer
162: .seealso: MFNConvergedReasonView()
163: @*/
164: PetscErrorCode MFNConvergedReasonViewFromOptions(MFN mfn)
165: {
166: PetscErrorCode ierr;
167: PetscViewer viewer;
168: PetscBool flg;
169: static PetscBool incall = PETSC_FALSE;
170: PetscViewerFormat format;
173: if (incall) return(0);
174: incall = PETSC_TRUE;
175: PetscOptionsGetViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->options,((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg);
176: if (flg) {
177: PetscViewerPushFormat(viewer,format);
178: MFNConvergedReasonView(mfn,viewer);
179: PetscViewerPopFormat(viewer);
180: PetscViewerDestroy(&viewer);
181: }
182: incall = PETSC_FALSE;
183: return(0);
184: }
186: /*@
187: MFNCreate - Creates the default MFN context.
189: Collective
191: Input Parameter:
192: . comm - MPI communicator
194: Output Parameter:
195: . mfn - location to put the MFN context
197: Note:
198: The default MFN type is MFNKRYLOV
200: Level: beginner
202: .seealso: MFNSetUp(), MFNSolve(), MFNDestroy(), MFN
203: @*/
204: PetscErrorCode MFNCreate(MPI_Comm comm,MFN *outmfn)
205: {
207: MFN mfn;
211: *outmfn = 0;
212: MFNInitializePackage();
213: SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView);
215: mfn->A = NULL;
216: mfn->fn = NULL;
217: mfn->max_it = PETSC_DEFAULT;
218: mfn->ncv = PETSC_DEFAULT;
219: mfn->tol = PETSC_DEFAULT;
220: mfn->errorifnotconverged = PETSC_FALSE;
222: mfn->numbermonitors = 0;
224: mfn->V = NULL;
225: mfn->nwork = 0;
226: mfn->work = NULL;
227: mfn->data = NULL;
229: mfn->its = 0;
230: mfn->nv = 0;
231: mfn->errest = 0;
232: mfn->setupcalled = 0;
233: mfn->reason = MFN_CONVERGED_ITERATING;
235: *outmfn = mfn;
236: return(0);
237: }
239: /*@C
240: MFNSetType - Selects the particular solver to be used in the MFN object.
242: Logically Collective on mfn
244: Input Parameters:
245: + mfn - the matrix function context
246: - type - a known method
248: Options Database Key:
249: . -mfn_type <method> - Sets the method; use -help for a list
250: of available methods
252: Notes:
253: See "slepc/include/slepcmfn.h" for available methods. The default
254: is MFNKRYLOV
256: Normally, it is best to use the MFNSetFromOptions() command and
257: then set the MFN type from the options database rather than by using
258: this routine. Using the options database provides the user with
259: maximum flexibility in evaluating the different available methods.
260: The MFNSetType() routine is provided for those situations where it
261: is necessary to set the iterative solver independently of the command
262: line or options database.
264: Level: intermediate
266: .seealso: MFNType
267: @*/
268: PetscErrorCode MFNSetType(MFN mfn,MFNType type)
269: {
270: PetscErrorCode ierr,(*r)(MFN);
271: PetscBool match;
277: PetscObjectTypeCompare((PetscObject)mfn,type,&match);
278: if (match) return(0);
280: PetscFunctionListFind(MFNList,type,&r);
281: if (!r) SETERRQ1(PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown MFN type given: %s",type);
283: if (mfn->ops->destroy) { (*mfn->ops->destroy)(mfn); }
284: PetscMemzero(mfn->ops,sizeof(struct _MFNOps));
286: mfn->setupcalled = 0;
287: PetscObjectChangeTypeName((PetscObject)mfn,type);
288: (*r)(mfn);
289: return(0);
290: }
292: /*@C
293: MFNGetType - Gets the MFN type as a string from the MFN object.
295: Not Collective
297: Input Parameter:
298: . mfn - the matrix function context
300: Output Parameter:
301: . name - name of MFN method
303: Level: intermediate
305: .seealso: MFNSetType()
306: @*/
307: PetscErrorCode MFNGetType(MFN mfn,MFNType *type)
308: {
312: *type = ((PetscObject)mfn)->type_name;
313: return(0);
314: }
316: /*@C
317: MFNRegister - Adds a method to the matrix function solver package.
319: Not Collective
321: Input Parameters:
322: + name - name of a new user-defined solver
323: - function - routine to create the solver context
325: Notes:
326: MFNRegister() may be called multiple times to add several user-defined solvers.
328: Sample usage:
329: .vb
330: MFNRegister("my_solver",MySolverCreate);
331: .ve
333: Then, your solver can be chosen with the procedural interface via
334: $ MFNSetType(mfn,"my_solver")
335: or at runtime via the option
336: $ -mfn_type my_solver
338: Level: advanced
340: .seealso: MFNRegisterAll()
341: @*/
342: PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
343: {
347: MFNInitializePackage();
348: PetscFunctionListAdd(&MFNList,name,function);
349: return(0);
350: }
352: /*@
353: MFNReset - Resets the MFN context to the initial state (prior to setup)
354: and destroys any allocated Vecs and Mats.
356: Collective on mfn
358: Input Parameter:
359: . mfn - matrix function context obtained from MFNCreate()
361: Level: advanced
363: .seealso: MFNDestroy()
364: @*/
365: PetscErrorCode MFNReset(MFN mfn)
366: {
371: if (!mfn) return(0);
372: if (mfn->ops->reset) { (mfn->ops->reset)(mfn); }
373: MatDestroy(&mfn->A);
374: BVDestroy(&mfn->V);
375: VecDestroyVecs(mfn->nwork,&mfn->work);
376: mfn->nwork = 0;
377: mfn->setupcalled = 0;
378: return(0);
379: }
381: /*@
382: MFNDestroy - Destroys the MFN context.
384: Collective on mfn
386: Input Parameter:
387: . mfn - matrix function context obtained from MFNCreate()
389: Level: beginner
391: .seealso: MFNCreate(), MFNSetUp(), MFNSolve()
392: @*/
393: PetscErrorCode MFNDestroy(MFN *mfn)
394: {
398: if (!*mfn) return(0);
400: if (--((PetscObject)(*mfn))->refct > 0) { *mfn = 0; return(0); }
401: MFNReset(*mfn);
402: if ((*mfn)->ops->destroy) { (*(*mfn)->ops->destroy)(*mfn); }
403: FNDestroy(&(*mfn)->fn);
404: MatDestroy(&(*mfn)->AT);
405: MFNMonitorCancel(*mfn);
406: PetscHeaderDestroy(mfn);
407: return(0);
408: }
410: /*@
411: MFNSetBV - Associates a basis vectors object to the matrix function solver.
413: Collective on mfn
415: Input Parameters:
416: + mfn - matrix function context obtained from MFNCreate()
417: - bv - the basis vectors object
419: Note:
420: Use MFNGetBV() to retrieve the basis vectors context (for example,
421: to free it at the end of the computations).
423: Level: advanced
425: .seealso: MFNGetBV()
426: @*/
427: PetscErrorCode MFNSetBV(MFN mfn,BV bv)
428: {
435: PetscObjectReference((PetscObject)bv);
436: BVDestroy(&mfn->V);
437: mfn->V = bv;
438: PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->V);
439: return(0);
440: }
442: /*@
443: MFNGetBV - Obtain the basis vectors object associated to the matrix
444: function solver.
446: Not Collective
448: Input Parameters:
449: . mfn - matrix function context obtained from MFNCreate()
451: Output Parameter:
452: . bv - basis vectors context
454: Level: advanced
456: .seealso: MFNSetBV()
457: @*/
458: PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
459: {
465: if (!mfn->V) {
466: BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V);
467: PetscObjectIncrementTabLevel((PetscObject)mfn->V,(PetscObject)mfn,0);
468: PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->V);
469: PetscObjectSetOptions((PetscObject)mfn->V,((PetscObject)mfn)->options);
470: }
471: *bv = mfn->V;
472: return(0);
473: }
475: /*@
476: MFNSetFN - Specifies the function to be computed.
478: Collective on mfn
480: Input Parameters:
481: + mfn - matrix function context obtained from MFNCreate()
482: - fn - the math function object
484: Note:
485: Use MFNGetFN() to retrieve the math function context (for example,
486: to free it at the end of the computations).
488: Level: beginner
490: .seealso: MFNGetFN()
491: @*/
492: PetscErrorCode MFNSetFN(MFN mfn,FN fn)
493: {
500: PetscObjectReference((PetscObject)fn);
501: FNDestroy(&mfn->fn);
502: mfn->fn = fn;
503: PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->fn);
504: return(0);
505: }
507: /*@
508: MFNGetFN - Obtain the math function object associated to the MFN object.
510: Not Collective
512: Input Parameters:
513: . mfn - matrix function context obtained from MFNCreate()
515: Output Parameter:
516: . fn - math function context
518: Level: beginner
520: .seealso: MFNSetFN()
521: @*/
522: PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
523: {
529: if (!mfn->fn) {
530: FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn);
531: PetscObjectIncrementTabLevel((PetscObject)mfn->fn,(PetscObject)mfn,0);
532: PetscLogObjectParent((PetscObject)mfn,(PetscObject)mfn->fn);
533: PetscObjectSetOptions((PetscObject)mfn->fn,((PetscObject)mfn)->options);
534: }
535: *fn = mfn->fn;
536: return(0);
537: }