Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
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 : */
13 :
14 : #include <slepc/private/mfnimpl.h> /*I "slepcmfn.h" I*/
15 :
16 : /* Logging support */
17 : PetscClassId MFN_CLASSID = 0;
18 : PetscLogEvent MFN_SetUp = 0,MFN_Solve = 0;
19 :
20 : /* List of registered MFN routines */
21 : PetscFunctionList MFNList = NULL;
22 : PetscBool MFNRegisterAllCalled = PETSC_FALSE;
23 :
24 : /* List of registered MFN monitors */
25 : PetscFunctionList MFNMonitorList = NULL;
26 : PetscFunctionList MFNMonitorCreateList = NULL;
27 : PetscFunctionList MFNMonitorDestroyList = NULL;
28 : PetscBool MFNMonitorRegisterAllCalled = PETSC_FALSE;
29 :
30 : /*@C
31 : MFNView - Prints the MFN data structure.
32 :
33 : Collective
34 :
35 : Input Parameters:
36 : + mfn - the matrix function solver context
37 : - viewer - optional visualization context
38 :
39 : Options Database Key:
40 : . -mfn_view - Calls MFNView() at end of MFNSolve()
41 :
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.
49 :
50 : The user can open an alternative visualization context with
51 : PetscViewerASCIIOpen() - output to a specified file.
52 :
53 : Level: beginner
54 :
55 : .seealso: MFNCreate()
56 : @*/
57 2 : PetscErrorCode MFNView(MFN mfn,PetscViewer viewer)
58 : {
59 2 : PetscBool isascii;
60 :
61 2 : PetscFunctionBegin;
62 2 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
63 2 : if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)mfn),&viewer));
64 2 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
65 2 : PetscCheckSameComm(mfn,1,viewer,2);
66 :
67 2 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
68 2 : if (isascii) {
69 2 : PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)mfn,viewer));
70 2 : PetscCall(PetscViewerASCIIPushTab(viewer));
71 2 : PetscTryTypeMethod(mfn,view,viewer);
72 2 : PetscCall(PetscViewerASCIIPopTab(viewer));
73 2 : PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",mfn->ncv));
74 2 : PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",mfn->max_it));
75 2 : PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)mfn->tol));
76 0 : } else PetscTryTypeMethod(mfn,view,viewer);
77 2 : PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
78 2 : if (!mfn->V) PetscCall(MFNGetFN(mfn,&mfn->fn));
79 2 : PetscCall(FNView(mfn->fn,viewer));
80 2 : if (!mfn->V) PetscCall(MFNGetBV(mfn,&mfn->V));
81 2 : PetscCall(BVView(mfn->V,viewer));
82 2 : PetscCall(PetscViewerPopFormat(viewer));
83 2 : PetscFunctionReturn(PETSC_SUCCESS);
84 : }
85 :
86 : /*@C
87 : MFNViewFromOptions - View from options
88 :
89 : Collective
90 :
91 : Input Parameters:
92 : + mfn - the matrix function context
93 : . obj - optional object
94 : - name - command line option
95 :
96 : Level: intermediate
97 :
98 : .seealso: MFNView(), MFNCreate()
99 : @*/
100 202 : PetscErrorCode MFNViewFromOptions(MFN mfn,PetscObject obj,const char name[])
101 : {
102 202 : PetscFunctionBegin;
103 202 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
104 202 : PetscCall(PetscObjectViewFromOptions((PetscObject)mfn,obj,name));
105 202 : PetscFunctionReturn(PETSC_SUCCESS);
106 : }
107 : /*@C
108 : MFNConvergedReasonView - Displays the reason an MFN solve converged or diverged.
109 :
110 : Collective
111 :
112 : Input Parameters:
113 : + mfn - the matrix function context
114 : - viewer - the viewer to display the reason
115 :
116 : Options Database Keys:
117 : . -mfn_converged_reason - print reason for convergence, and number of iterations
118 :
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
124 :
125 : Level: intermediate
126 :
127 : .seealso: MFNSetTolerances(), MFNGetIterationNumber(), MFNConvergedReasonViewFromOptions()
128 : @*/
129 1 : PetscErrorCode MFNConvergedReasonView(MFN mfn,PetscViewer viewer)
130 : {
131 1 : PetscBool isAscii;
132 1 : PetscViewerFormat format;
133 :
134 1 : PetscFunctionBegin;
135 1 : if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)mfn));
136 1 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
137 1 : if (isAscii) {
138 1 : PetscCall(PetscViewerGetFormat(viewer,&format));
139 1 : PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)mfn)->tablevel));
140 1 : 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 0 : 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 1 : PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)mfn)->tablevel));
143 : }
144 1 : PetscFunctionReturn(PETSC_SUCCESS);
145 : }
146 :
147 : /*@
148 : MFNConvergedReasonViewFromOptions - Processes command line options to determine if/how
149 : the MFN converged reason is to be viewed.
150 :
151 : Collective
152 :
153 : Input Parameter:
154 : . mfn - the matrix function context
155 :
156 : Level: developer
157 :
158 : .seealso: MFNConvergedReasonView()
159 : @*/
160 101 : PetscErrorCode MFNConvergedReasonViewFromOptions(MFN mfn)
161 : {
162 101 : PetscViewer viewer;
163 101 : PetscBool flg;
164 101 : static PetscBool incall = PETSC_FALSE;
165 101 : PetscViewerFormat format;
166 :
167 101 : PetscFunctionBegin;
168 101 : if (incall) PetscFunctionReturn(PETSC_SUCCESS);
169 101 : incall = PETSC_TRUE;
170 101 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->options,((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg));
171 101 : if (flg) {
172 1 : PetscCall(PetscViewerPushFormat(viewer,format));
173 1 : PetscCall(MFNConvergedReasonView(mfn,viewer));
174 1 : PetscCall(PetscViewerPopFormat(viewer));
175 1 : PetscCall(PetscOptionsRestoreViewer(&viewer));
176 : }
177 101 : incall = PETSC_FALSE;
178 101 : PetscFunctionReturn(PETSC_SUCCESS);
179 : }
180 :
181 : /*@
182 : MFNCreate - Creates the default MFN context.
183 :
184 : Collective
185 :
186 : Input Parameter:
187 : . comm - MPI communicator
188 :
189 : Output Parameter:
190 : . outmfn - location to put the MFN context
191 :
192 : Note:
193 : The default MFN type is MFNKRYLOV
194 :
195 : Level: beginner
196 :
197 : .seealso: MFNSetUp(), MFNSolve(), MFNDestroy(), MFN
198 : @*/
199 17 : PetscErrorCode MFNCreate(MPI_Comm comm,MFN *outmfn)
200 : {
201 17 : MFN mfn;
202 :
203 17 : PetscFunctionBegin;
204 17 : PetscAssertPointer(outmfn,2);
205 17 : *outmfn = NULL;
206 17 : PetscCall(MFNInitializePackage());
207 17 : PetscCall(SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView));
208 :
209 17 : mfn->A = NULL;
210 17 : mfn->fn = NULL;
211 17 : mfn->max_it = PETSC_DEFAULT;
212 17 : mfn->ncv = PETSC_DEFAULT;
213 17 : mfn->tol = PETSC_DEFAULT;
214 17 : mfn->errorifnotconverged = PETSC_FALSE;
215 :
216 17 : mfn->numbermonitors = 0;
217 :
218 17 : mfn->V = NULL;
219 17 : mfn->nwork = 0;
220 17 : mfn->work = NULL;
221 17 : mfn->data = NULL;
222 :
223 17 : mfn->its = 0;
224 17 : mfn->nv = 0;
225 17 : mfn->errest = 0;
226 17 : mfn->setupcalled = 0;
227 17 : mfn->reason = MFN_CONVERGED_ITERATING;
228 :
229 17 : *outmfn = mfn;
230 17 : PetscFunctionReturn(PETSC_SUCCESS);
231 : }
232 :
233 : /*@C
234 : MFNSetType - Selects the particular solver to be used in the MFN object.
235 :
236 : Logically Collective
237 :
238 : Input Parameters:
239 : + mfn - the matrix function context
240 : - type - a known method
241 :
242 : Options Database Key:
243 : . -mfn_type <method> - Sets the method; use -help for a list
244 : of available methods
245 :
246 : Notes:
247 : See "slepc/include/slepcmfn.h" for available methods. The default
248 : is MFNKRYLOV
249 :
250 : Normally, it is best to use the MFNSetFromOptions() command and
251 : then set the MFN type from the options database rather than by using
252 : this routine. Using the options database provides the user with
253 : maximum flexibility in evaluating the different available methods.
254 : The MFNSetType() routine is provided for those situations where it
255 : is necessary to set the iterative solver independently of the command
256 : line or options database.
257 :
258 : Level: intermediate
259 :
260 : .seealso: MFNType
261 : @*/
262 18 : PetscErrorCode MFNSetType(MFN mfn,MFNType type)
263 : {
264 18 : PetscErrorCode (*r)(MFN);
265 18 : PetscBool match;
266 :
267 18 : PetscFunctionBegin;
268 18 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
269 18 : PetscAssertPointer(type,2);
270 :
271 18 : PetscCall(PetscObjectTypeCompare((PetscObject)mfn,type,&match));
272 18 : if (match) PetscFunctionReturn(PETSC_SUCCESS);
273 :
274 18 : PetscCall(PetscFunctionListFind(MFNList,type,&r));
275 18 : PetscCheck(r,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown MFN type given: %s",type);
276 :
277 18 : PetscTryTypeMethod(mfn,destroy);
278 18 : PetscCall(PetscMemzero(mfn->ops,sizeof(struct _MFNOps)));
279 :
280 18 : mfn->setupcalled = 0;
281 18 : PetscCall(PetscObjectChangeTypeName((PetscObject)mfn,type));
282 18 : PetscCall((*r)(mfn));
283 18 : PetscFunctionReturn(PETSC_SUCCESS);
284 : }
285 :
286 : /*@C
287 : MFNGetType - Gets the MFN type as a string from the MFN object.
288 :
289 : Not Collective
290 :
291 : Input Parameter:
292 : . mfn - the matrix function context
293 :
294 : Output Parameter:
295 : . type - name of MFN method
296 :
297 : Level: intermediate
298 :
299 : .seealso: MFNSetType()
300 : @*/
301 2 : PetscErrorCode MFNGetType(MFN mfn,MFNType *type)
302 : {
303 2 : PetscFunctionBegin;
304 2 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
305 2 : PetscAssertPointer(type,2);
306 2 : *type = ((PetscObject)mfn)->type_name;
307 2 : PetscFunctionReturn(PETSC_SUCCESS);
308 : }
309 :
310 : /*@C
311 : MFNRegister - Adds a method to the matrix function solver package.
312 :
313 : Not Collective
314 :
315 : Input Parameters:
316 : + name - name of a new user-defined solver
317 : - function - routine to create the solver context
318 :
319 : Notes:
320 : MFNRegister() may be called multiple times to add several user-defined solvers.
321 :
322 : Example Usage:
323 : .vb
324 : MFNRegister("my_solver",MySolverCreate);
325 : .ve
326 :
327 : Then, your solver can be chosen with the procedural interface via
328 : $ MFNSetType(mfn,"my_solver")
329 : or at runtime via the option
330 : $ -mfn_type my_solver
331 :
332 : Level: advanced
333 :
334 : .seealso: MFNRegisterAll()
335 : @*/
336 32 : PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
337 : {
338 32 : PetscFunctionBegin;
339 32 : PetscCall(MFNInitializePackage());
340 32 : PetscCall(PetscFunctionListAdd(&MFNList,name,function));
341 32 : PetscFunctionReturn(PETSC_SUCCESS);
342 : }
343 :
344 : /*@C
345 : MFNMonitorRegister - Adds MFN monitor routine.
346 :
347 : Not Collective
348 :
349 : Input Parameters:
350 : + name - name of a new monitor routine
351 : . vtype - a PetscViewerType for the output
352 : . format - a PetscViewerFormat for the output
353 : . monitor - monitor routine
354 : . create - creation routine, or NULL
355 : - destroy - destruction routine, or NULL
356 :
357 : Notes:
358 : MFNMonitorRegister() may be called multiple times to add several user-defined monitors.
359 :
360 : Example Usage:
361 : .vb
362 : MFNMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
363 : .ve
364 :
365 : Then, your monitor can be chosen with the procedural interface via
366 : $ MFNMonitorSetFromOptions(mfn,"-mfn_monitor_my_monitor","my_monitor",NULL)
367 : or at runtime via the option
368 : $ -mfn_monitor_my_monitor
369 :
370 : Level: advanced
371 :
372 : .seealso: MFNMonitorRegisterAll()
373 : @*/
374 32 : PetscErrorCode MFNMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(MFN,PetscInt,PetscReal,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
375 : {
376 32 : char key[PETSC_MAX_PATH_LEN];
377 :
378 32 : PetscFunctionBegin;
379 32 : PetscCall(MFNInitializePackage());
380 32 : PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
381 32 : PetscCall(PetscFunctionListAdd(&MFNMonitorList,key,monitor));
382 32 : if (create) PetscCall(PetscFunctionListAdd(&MFNMonitorCreateList,key,create));
383 32 : if (destroy) PetscCall(PetscFunctionListAdd(&MFNMonitorDestroyList,key,destroy));
384 32 : PetscFunctionReturn(PETSC_SUCCESS);
385 : }
386 :
387 : /*@
388 : MFNReset - Resets the MFN context to the initial state (prior to setup)
389 : and destroys any allocated Vecs and Mats.
390 :
391 : Collective
392 :
393 : Input Parameter:
394 : . mfn - matrix function context obtained from MFNCreate()
395 :
396 : Level: advanced
397 :
398 : .seealso: MFNDestroy()
399 : @*/
400 17 : PetscErrorCode MFNReset(MFN mfn)
401 : {
402 17 : PetscFunctionBegin;
403 17 : if (mfn) PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
404 17 : if (!mfn) PetscFunctionReturn(PETSC_SUCCESS);
405 17 : PetscTryTypeMethod(mfn,reset);
406 17 : PetscCall(MatDestroy(&mfn->A));
407 17 : PetscCall(BVDestroy(&mfn->V));
408 17 : PetscCall(VecDestroyVecs(mfn->nwork,&mfn->work));
409 17 : mfn->nwork = 0;
410 17 : mfn->setupcalled = 0;
411 17 : PetscFunctionReturn(PETSC_SUCCESS);
412 : }
413 :
414 : /*@C
415 : MFNDestroy - Destroys the MFN context.
416 :
417 : Collective
418 :
419 : Input Parameter:
420 : . mfn - matrix function context obtained from MFNCreate()
421 :
422 : Level: beginner
423 :
424 : .seealso: MFNCreate(), MFNSetUp(), MFNSolve()
425 : @*/
426 17 : PetscErrorCode MFNDestroy(MFN *mfn)
427 : {
428 17 : PetscFunctionBegin;
429 17 : if (!*mfn) PetscFunctionReturn(PETSC_SUCCESS);
430 17 : PetscValidHeaderSpecific(*mfn,MFN_CLASSID,1);
431 17 : if (--((PetscObject)*mfn)->refct > 0) { *mfn = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
432 17 : PetscCall(MFNReset(*mfn));
433 17 : PetscTryTypeMethod(*mfn,destroy);
434 17 : PetscCall(FNDestroy(&(*mfn)->fn));
435 17 : PetscCall(MatDestroy(&(*mfn)->AT));
436 17 : PetscCall(MFNMonitorCancel(*mfn));
437 17 : PetscCall(PetscHeaderDestroy(mfn));
438 17 : PetscFunctionReturn(PETSC_SUCCESS);
439 : }
440 :
441 : /*@
442 : MFNSetBV - Associates a basis vectors object to the matrix function solver.
443 :
444 : Collective
445 :
446 : Input Parameters:
447 : + mfn - matrix function context obtained from MFNCreate()
448 : - bv - the basis vectors object
449 :
450 : Note:
451 : Use MFNGetBV() to retrieve the basis vectors context (for example,
452 : to free it at the end of the computations).
453 :
454 : Level: advanced
455 :
456 : .seealso: MFNGetBV()
457 : @*/
458 1 : PetscErrorCode MFNSetBV(MFN mfn,BV bv)
459 : {
460 1 : PetscFunctionBegin;
461 1 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
462 1 : PetscValidHeaderSpecific(bv,BV_CLASSID,2);
463 1 : PetscCheckSameComm(mfn,1,bv,2);
464 1 : PetscCall(PetscObjectReference((PetscObject)bv));
465 1 : PetscCall(BVDestroy(&mfn->V));
466 1 : mfn->V = bv;
467 1 : PetscFunctionReturn(PETSC_SUCCESS);
468 : }
469 :
470 : /*@
471 : MFNGetBV - Obtain the basis vectors object associated to the matrix
472 : function solver.
473 :
474 : Not Collective
475 :
476 : Input Parameters:
477 : . mfn - matrix function context obtained from MFNCreate()
478 :
479 : Output Parameter:
480 : . bv - basis vectors context
481 :
482 : Level: advanced
483 :
484 : .seealso: MFNSetBV()
485 : @*/
486 16 : PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
487 : {
488 16 : PetscFunctionBegin;
489 16 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
490 16 : PetscAssertPointer(bv,2);
491 16 : if (!mfn->V) {
492 16 : PetscCall(BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V));
493 16 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->V,(PetscObject)mfn,0));
494 16 : PetscCall(PetscObjectSetOptions((PetscObject)mfn->V,((PetscObject)mfn)->options));
495 : }
496 16 : *bv = mfn->V;
497 16 : PetscFunctionReturn(PETSC_SUCCESS);
498 : }
499 :
500 : /*@
501 : MFNSetFN - Specifies the function to be computed.
502 :
503 : Collective
504 :
505 : Input Parameters:
506 : + mfn - matrix function context obtained from MFNCreate()
507 : - fn - the math function object
508 :
509 : Note:
510 : Use MFNGetFN() to retrieve the math function context (for example,
511 : to free it at the end of the computations).
512 :
513 : Level: beginner
514 :
515 : .seealso: MFNGetFN()
516 : @*/
517 4 : PetscErrorCode MFNSetFN(MFN mfn,FN fn)
518 : {
519 4 : PetscFunctionBegin;
520 4 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
521 4 : PetscValidHeaderSpecific(fn,FN_CLASSID,2);
522 4 : PetscCheckSameComm(mfn,1,fn,2);
523 4 : PetscCall(PetscObjectReference((PetscObject)fn));
524 4 : PetscCall(FNDestroy(&mfn->fn));
525 4 : mfn->fn = fn;
526 4 : PetscFunctionReturn(PETSC_SUCCESS);
527 : }
528 :
529 : /*@
530 : MFNGetFN - Obtain the math function object associated to the MFN object.
531 :
532 : Not Collective
533 :
534 : Input Parameters:
535 : . mfn - matrix function context obtained from MFNCreate()
536 :
537 : Output Parameter:
538 : . fn - math function context
539 :
540 : Level: beginner
541 :
542 : .seealso: MFNSetFN()
543 : @*/
544 370 : PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
545 : {
546 370 : PetscFunctionBegin;
547 370 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
548 370 : PetscAssertPointer(fn,2);
549 370 : if (!mfn->fn) {
550 13 : PetscCall(FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn));
551 13 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->fn,(PetscObject)mfn,0));
552 13 : PetscCall(PetscObjectSetOptions((PetscObject)mfn->fn,((PetscObject)mfn)->options));
553 : }
554 370 : *fn = mfn->fn;
555 370 : PetscFunctionReturn(PETSC_SUCCESS);
556 : }
|
