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 : /*@
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 : /*@
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 206 : PetscErrorCode MFNViewFromOptions(MFN mfn,PetscObject obj,const char name[])
101 : {
102 206 : PetscFunctionBegin;
103 206 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
104 206 : PetscCall(PetscObjectViewFromOptions((PetscObject)mfn,obj,name));
105 206 : PetscFunctionReturn(PETSC_SUCCESS);
106 : }
107 : /*@
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 103 : PetscErrorCode MFNConvergedReasonViewFromOptions(MFN mfn)
161 : {
162 103 : PetscViewer viewer;
163 103 : PetscBool flg;
164 103 : static PetscBool incall = PETSC_FALSE;
165 103 : PetscViewerFormat format;
166 :
167 103 : PetscFunctionBegin;
168 103 : if (incall) PetscFunctionReturn(PETSC_SUCCESS);
169 103 : incall = PETSC_TRUE;
170 103 : PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)mfn),((PetscObject)mfn)->options,((PetscObject)mfn)->prefix,"-mfn_converged_reason",&viewer,&format,&flg));
171 103 : if (flg) {
172 1 : PetscCall(PetscViewerPushFormat(viewer,format));
173 1 : PetscCall(MFNConvergedReasonView(mfn,viewer));
174 1 : PetscCall(PetscViewerPopFormat(viewer));
175 1 : PetscCall(PetscViewerDestroy(&viewer));
176 : }
177 103 : incall = PETSC_FALSE;
178 103 : 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 : PetscCall(MFNInitializePackage());
206 17 : PetscCall(SlepcHeaderCreate(mfn,MFN_CLASSID,"MFN","Matrix Function","MFN",comm,MFNDestroy,MFNView));
207 :
208 17 : mfn->A = NULL;
209 17 : mfn->fn = NULL;
210 17 : mfn->max_it = PETSC_DETERMINE;
211 17 : mfn->ncv = PETSC_DETERMINE;
212 17 : mfn->tol = PETSC_DETERMINE;
213 17 : mfn->errorifnotconverged = PETSC_FALSE;
214 :
215 17 : mfn->numbermonitors = 0;
216 :
217 17 : mfn->V = NULL;
218 17 : mfn->nwork = 0;
219 17 : mfn->work = NULL;
220 17 : mfn->data = NULL;
221 :
222 17 : mfn->its = 0;
223 17 : mfn->nv = 0;
224 17 : mfn->errest = 0;
225 17 : mfn->setupcalled = 0;
226 17 : mfn->reason = MFN_CONVERGED_ITERATING;
227 :
228 17 : *outmfn = mfn;
229 17 : PetscFunctionReturn(PETSC_SUCCESS);
230 : }
231 :
232 : /*@
233 : MFNSetType - Selects the particular solver to be used in the MFN object.
234 :
235 : Logically Collective
236 :
237 : Input Parameters:
238 : + mfn - the matrix function context
239 : - type - a known method
240 :
241 : Options Database Key:
242 : . -mfn_type <method> - Sets the method; use -help for a list
243 : of available methods
244 :
245 : Notes:
246 : See "slepc/include/slepcmfn.h" for available methods. The default
247 : is MFNKRYLOV
248 :
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.
256 :
257 : Level: intermediate
258 :
259 : .seealso: MFNType
260 : @*/
261 18 : PetscErrorCode MFNSetType(MFN mfn,MFNType type)
262 : {
263 18 : PetscErrorCode (*r)(MFN);
264 18 : PetscBool match;
265 :
266 18 : PetscFunctionBegin;
267 18 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
268 18 : PetscAssertPointer(type,2);
269 :
270 18 : PetscCall(PetscObjectTypeCompare((PetscObject)mfn,type,&match));
271 18 : if (match) PetscFunctionReturn(PETSC_SUCCESS);
272 :
273 18 : PetscCall(PetscFunctionListFind(MFNList,type,&r));
274 18 : PetscCheck(r,PetscObjectComm((PetscObject)mfn),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown MFN type given: %s",type);
275 :
276 18 : PetscTryTypeMethod(mfn,destroy);
277 18 : PetscCall(PetscMemzero(mfn->ops,sizeof(struct _MFNOps)));
278 :
279 18 : mfn->setupcalled = 0;
280 18 : PetscCall(PetscObjectChangeTypeName((PetscObject)mfn,type));
281 18 : PetscCall((*r)(mfn));
282 18 : PetscFunctionReturn(PETSC_SUCCESS);
283 : }
284 :
285 : /*@
286 : MFNGetType - Gets the MFN type as a string from the MFN object.
287 :
288 : Not Collective
289 :
290 : Input Parameter:
291 : . mfn - the matrix function context
292 :
293 : Output Parameter:
294 : . type - name of MFN method
295 :
296 : Level: intermediate
297 :
298 : .seealso: MFNSetType()
299 : @*/
300 2 : PetscErrorCode MFNGetType(MFN mfn,MFNType *type)
301 : {
302 2 : PetscFunctionBegin;
303 2 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
304 2 : PetscAssertPointer(type,2);
305 2 : *type = ((PetscObject)mfn)->type_name;
306 2 : PetscFunctionReturn(PETSC_SUCCESS);
307 : }
308 :
309 : /*@C
310 : MFNRegister - Adds a method to the matrix function solver package.
311 :
312 : Not Collective
313 :
314 : Input Parameters:
315 : + name - name of a new user-defined solver
316 : - function - routine to create the solver context
317 :
318 : Notes:
319 : MFNRegister() may be called multiple times to add several user-defined solvers.
320 :
321 : Example Usage:
322 : .vb
323 : MFNRegister("my_solver",MySolverCreate);
324 : .ve
325 :
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
330 :
331 : Level: advanced
332 :
333 : .seealso: MFNRegisterAll()
334 : @*/
335 32 : PetscErrorCode MFNRegister(const char *name,PetscErrorCode (*function)(MFN))
336 : {
337 32 : PetscFunctionBegin;
338 32 : PetscCall(MFNInitializePackage());
339 32 : PetscCall(PetscFunctionListAdd(&MFNList,name,function));
340 32 : PetscFunctionReturn(PETSC_SUCCESS);
341 : }
342 :
343 : /*@C
344 : MFNMonitorRegister - Adds MFN monitor routine.
345 :
346 : Not Collective
347 :
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
353 : . create - creation routine, or NULL
354 : - destroy - destruction routine, or NULL
355 :
356 : Notes:
357 : MFNMonitorRegister() may be called multiple times to add several user-defined monitors.
358 :
359 : Example Usage:
360 : .vb
361 : MFNMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
362 : .ve
363 :
364 : Then, your monitor can be chosen with the procedural interface via
365 : $ MFNMonitorSetFromOptions(mfn,"-mfn_monitor_my_monitor","my_monitor",NULL)
366 : or at runtime via the option
367 : $ -mfn_monitor_my_monitor
368 :
369 : Level: advanced
370 :
371 : .seealso: MFNMonitorRegisterAll()
372 : @*/
373 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**))
374 : {
375 32 : char key[PETSC_MAX_PATH_LEN];
376 :
377 32 : PetscFunctionBegin;
378 32 : PetscCall(MFNInitializePackage());
379 32 : PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
380 32 : PetscCall(PetscFunctionListAdd(&MFNMonitorList,key,monitor));
381 32 : if (create) PetscCall(PetscFunctionListAdd(&MFNMonitorCreateList,key,create));
382 32 : if (destroy) PetscCall(PetscFunctionListAdd(&MFNMonitorDestroyList,key,destroy));
383 32 : PetscFunctionReturn(PETSC_SUCCESS);
384 : }
385 :
386 : /*@
387 : MFNReset - Resets the MFN context to the initial state (prior to setup)
388 : and destroys any allocated Vecs and Mats.
389 :
390 : Collective
391 :
392 : Input Parameter:
393 : . mfn - matrix function context obtained from MFNCreate()
394 :
395 : Level: advanced
396 :
397 : .seealso: MFNDestroy()
398 : @*/
399 17 : PetscErrorCode MFNReset(MFN mfn)
400 : {
401 17 : PetscFunctionBegin;
402 17 : if (mfn) PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
403 17 : if (!mfn) PetscFunctionReturn(PETSC_SUCCESS);
404 17 : PetscTryTypeMethod(mfn,reset);
405 17 : PetscCall(MatDestroy(&mfn->A));
406 17 : PetscCall(BVDestroy(&mfn->V));
407 17 : PetscCall(VecDestroyVecs(mfn->nwork,&mfn->work));
408 17 : mfn->nwork = 0;
409 17 : mfn->setupcalled = 0;
410 17 : PetscFunctionReturn(PETSC_SUCCESS);
411 : }
412 :
413 : /*@
414 : MFNDestroy - Destroys the MFN context.
415 :
416 : Collective
417 :
418 : Input Parameter:
419 : . mfn - matrix function context obtained from MFNCreate()
420 :
421 : Level: beginner
422 :
423 : .seealso: MFNCreate(), MFNSetUp(), MFNSolve()
424 : @*/
425 17 : PetscErrorCode MFNDestroy(MFN *mfn)
426 : {
427 17 : PetscFunctionBegin;
428 17 : if (!*mfn) PetscFunctionReturn(PETSC_SUCCESS);
429 17 : PetscValidHeaderSpecific(*mfn,MFN_CLASSID,1);
430 17 : if (--((PetscObject)*mfn)->refct > 0) { *mfn = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
431 17 : PetscCall(MFNReset(*mfn));
432 17 : PetscTryTypeMethod(*mfn,destroy);
433 17 : PetscCall(FNDestroy(&(*mfn)->fn));
434 17 : PetscCall(MatDestroy(&(*mfn)->AT));
435 17 : PetscCall(MFNMonitorCancel(*mfn));
436 17 : PetscCall(PetscHeaderDestroy(mfn));
437 17 : PetscFunctionReturn(PETSC_SUCCESS);
438 : }
439 :
440 : /*@
441 : MFNSetBV - Associates a basis vectors object to the matrix function solver.
442 :
443 : Collective
444 :
445 : Input Parameters:
446 : + mfn - matrix function context obtained from MFNCreate()
447 : - bv - the basis vectors object
448 :
449 : Note:
450 : Use MFNGetBV() to retrieve the basis vectors context (for example,
451 : to free it at the end of the computations).
452 :
453 : Level: advanced
454 :
455 : .seealso: MFNGetBV()
456 : @*/
457 1 : PetscErrorCode MFNSetBV(MFN mfn,BV bv)
458 : {
459 1 : PetscFunctionBegin;
460 1 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
461 1 : PetscValidHeaderSpecific(bv,BV_CLASSID,2);
462 1 : PetscCheckSameComm(mfn,1,bv,2);
463 1 : PetscCall(PetscObjectReference((PetscObject)bv));
464 1 : PetscCall(BVDestroy(&mfn->V));
465 1 : mfn->V = bv;
466 1 : PetscFunctionReturn(PETSC_SUCCESS);
467 : }
468 :
469 : /*@
470 : MFNGetBV - Obtain the basis vectors object associated to the matrix
471 : function solver.
472 :
473 : Not Collective
474 :
475 : Input Parameters:
476 : . mfn - matrix function context obtained from MFNCreate()
477 :
478 : Output Parameter:
479 : . bv - basis vectors context
480 :
481 : Level: advanced
482 :
483 : .seealso: MFNSetBV()
484 : @*/
485 16 : PetscErrorCode MFNGetBV(MFN mfn,BV *bv)
486 : {
487 16 : PetscFunctionBegin;
488 16 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
489 16 : PetscAssertPointer(bv,2);
490 16 : if (!mfn->V) {
491 16 : PetscCall(BVCreate(PetscObjectComm((PetscObject)mfn),&mfn->V));
492 16 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->V,(PetscObject)mfn,0));
493 16 : PetscCall(PetscObjectSetOptions((PetscObject)mfn->V,((PetscObject)mfn)->options));
494 : }
495 16 : *bv = mfn->V;
496 16 : PetscFunctionReturn(PETSC_SUCCESS);
497 : }
498 :
499 : /*@
500 : MFNSetFN - Specifies the function to be computed.
501 :
502 : Collective
503 :
504 : Input Parameters:
505 : + mfn - matrix function context obtained from MFNCreate()
506 : - fn - the math function object
507 :
508 : Note:
509 : Use MFNGetFN() to retrieve the math function context (for example,
510 : to free it at the end of the computations).
511 :
512 : Level: beginner
513 :
514 : .seealso: MFNGetFN()
515 : @*/
516 4 : PetscErrorCode MFNSetFN(MFN mfn,FN fn)
517 : {
518 4 : PetscFunctionBegin;
519 4 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
520 4 : PetscValidHeaderSpecific(fn,FN_CLASSID,2);
521 4 : PetscCheckSameComm(mfn,1,fn,2);
522 4 : PetscCall(PetscObjectReference((PetscObject)fn));
523 4 : PetscCall(FNDestroy(&mfn->fn));
524 4 : mfn->fn = fn;
525 4 : PetscFunctionReturn(PETSC_SUCCESS);
526 : }
527 :
528 : /*@
529 : MFNGetFN - Obtain the math function object associated to the MFN object.
530 :
531 : Not Collective
532 :
533 : Input Parameters:
534 : . mfn - matrix function context obtained from MFNCreate()
535 :
536 : Output Parameter:
537 : . fn - math function context
538 :
539 : Level: beginner
540 :
541 : .seealso: MFNSetFN()
542 : @*/
543 834 : PetscErrorCode MFNGetFN(MFN mfn,FN *fn)
544 : {
545 834 : PetscFunctionBegin;
546 834 : PetscValidHeaderSpecific(mfn,MFN_CLASSID,1);
547 834 : PetscAssertPointer(fn,2);
548 834 : if (!mfn->fn) {
549 13 : PetscCall(FNCreate(PetscObjectComm((PetscObject)mfn),&mfn->fn));
550 13 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)mfn->fn,(PetscObject)mfn,0));
551 13 : PetscCall(PetscObjectSetOptions((PetscObject)mfn->fn,((PetscObject)mfn)->options));
552 : }
553 834 : *fn = mfn->fn;
554 834 : PetscFunctionReturn(PETSC_SUCCESS);
555 : }
|
