Actual source code: nepmon.c
slepc-3.23.0 2025-03-29
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: NEP routines related to monitors
12: */
14: #include <slepc/private/nepimpl.h>
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
22: PetscInt i,n = nep->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: NEPMonitorSet - Sets an ADDITIONAL function to be called at every
31: iteration to monitor the error estimates for each requested eigenpair.
33: Logically Collective
35: Input Parameters:
36: + nep - eigensolver context obtained from NEPCreate()
37: . monitor - pointer to function (if this is NULL, it turns off monitoring)
38: . mctx - [optional] context for private data for the
39: monitor routine (use NULL if no context is desired)
40: - monitordestroy - [optional] routine that frees monitor context (may be NULL),
41: see PetscCtxDestroyFn for the calling sequence
43: Calling sequence of monitor:
44: $ PetscErrorCode monitor(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx)
45: + nep - nonlinear eigensolver context obtained from NEPCreate()
46: . its - iteration number
47: . nconv - number of converged eigenpairs
48: . eigr - real part of the eigenvalues
49: . eigi - imaginary part of the eigenvalues
50: . errest - error estimates for each eigenpair
51: . nest - number of error estimates
52: - mctx - optional monitoring context, as set by NEPMonitorSet()
54: Options Database Keys:
55: + -nep_monitor - print only the first error estimate
56: . -nep_monitor_all - print error estimates at each iteration
57: . -nep_monitor_conv - print the eigenvalue approximations only when
58: convergence has been reached
59: . -nep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
60: approximate eigenvalue
61: . -nep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
62: approximate eigenvalues
63: . -nep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
64: - -nep_monitor_cancel - cancels all monitors that have been hardwired into
65: a code by calls to NEPMonitorSet(), but does not cancel those set via
66: the options database.
68: Notes:
69: Several different monitoring routines may be set by calling
70: NEPMonitorSet() multiple times; all will be called in the
71: order in which they were set.
73: Level: intermediate
75: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
76: @*/
77: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscCtxDestroyFn *monitordestroy)
78: {
79: PetscFunctionBegin;
81: PetscCheck(nep->numbermonitors<MAXNEPMONITORS,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
82: nep->monitor[nep->numbermonitors] = monitor;
83: nep->monitorcontext[nep->numbermonitors] = (void*)mctx;
84: nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
85: PetscFunctionReturn(PETSC_SUCCESS);
86: }
88: /*@
89: NEPMonitorCancel - Clears all monitors for a NEP object.
91: Logically Collective
93: Input Parameters:
94: . nep - eigensolver context obtained from NEPCreate()
96: Options Database Key:
97: . -nep_monitor_cancel - Cancels all monitors that have been hardwired
98: into a code by calls to NEPMonitorSet(),
99: but does not cancel those set via the options database.
101: Level: intermediate
103: .seealso: NEPMonitorSet()
104: @*/
105: PetscErrorCode NEPMonitorCancel(NEP nep)
106: {
107: PetscInt i;
109: PetscFunctionBegin;
111: for (i=0; i<nep->numbermonitors; i++) {
112: if (nep->monitordestroy[i]) PetscCall((*nep->monitordestroy[i])(&nep->monitorcontext[i]));
113: }
114: nep->numbermonitors = 0;
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: /*@C
119: NEPGetMonitorContext - Gets the monitor context, as set by
120: NEPMonitorSet() for the FIRST monitor only.
122: Not Collective
124: Input Parameter:
125: . nep - eigensolver context obtained from NEPCreate()
127: Output Parameter:
128: . ctx - monitor context
130: Level: intermediate
132: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
133: @*/
134: PetscErrorCode NEPGetMonitorContext(NEP nep,void *ctx)
135: {
136: PetscFunctionBegin;
138: *(void**)ctx = nep->monitorcontext[0];
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: /*@C
143: NEPMonitorFirst - Print the first unconverged approximate value and
144: error estimate at each iteration of the nonlinear eigensolver.
146: Collective
148: Input Parameters:
149: + nep - nonlinear eigensolver context
150: . its - iteration number
151: . nconv - number of converged eigenpairs so far
152: . eigr - real part of the eigenvalues
153: . eigi - imaginary part of the eigenvalues
154: . errest - error estimates
155: . nest - number of error estimates to display
156: - vf - viewer and format for monitoring
158: Options Database Key:
159: . -nep_monitor - activates NEPMonitorFirst()
161: Level: intermediate
163: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
164: @*/
165: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
166: {
167: PetscViewer viewer = vf->viewer;
169: PetscFunctionBegin;
172: if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix));
173: if (nconv<nest) {
174: PetscCall(PetscViewerPushFormat(viewer,vf->format));
175: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
176: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
177: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
178: #if defined(PETSC_USE_COMPLEX)
179: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv])));
180: #else
181: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]));
182: if (eigi[nconv]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]));
183: #endif
184: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
185: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
186: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
187: PetscCall(PetscViewerPopFormat(viewer));
188: }
189: PetscFunctionReturn(PETSC_SUCCESS);
190: }
192: /*@C
193: NEPMonitorAll - Print the current approximate values and
194: error estimates at each iteration of the nonlinear eigensolver.
196: Collective
198: Input Parameters:
199: + nep - nonlinear eigensolver context
200: . its - iteration number
201: . nconv - number of converged eigenpairs so far
202: . eigr - real part of the eigenvalues
203: . eigi - imaginary part of the eigenvalues
204: . errest - error estimates
205: . nest - number of error estimates to display
206: - vf - viewer and format for monitoring
208: Options Database Key:
209: . -nep_monitor_all - activates NEPMonitorAll()
211: Level: intermediate
213: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
214: @*/
215: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
216: {
217: PetscInt i;
218: PetscViewer viewer = vf->viewer;
220: PetscFunctionBegin;
223: PetscCall(PetscViewerPushFormat(viewer,vf->format));
224: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
225: if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix));
226: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
227: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
228: for (i=0;i<nest;i++) {
229: #if defined(PETSC_USE_COMPLEX)
230: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i])));
231: #else
232: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]));
233: if (eigi[i]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]));
234: #endif
235: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
236: }
237: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
238: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
239: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
240: PetscCall(PetscViewerPopFormat(viewer));
241: PetscFunctionReturn(PETSC_SUCCESS);
242: }
244: /*@C
245: NEPMonitorConverged - Print the approximate values and
246: error estimates as they converge.
248: Collective
250: Input Parameters:
251: + nep - nonlinear eigensolver context
252: . its - iteration number
253: . nconv - number of converged eigenpairs so far
254: . eigr - real part of the eigenvalues
255: . eigi - imaginary part of the eigenvalues
256: . errest - error estimates
257: . nest - number of error estimates to display
258: - vf - viewer and format for monitoring
260: Options Database Key:
261: . -nep_monitor_conv - activates NEPMonitorConverged()
263: Level: intermediate
265: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
266: @*/
267: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
268: {
269: PetscInt i;
270: PetscViewer viewer = vf->viewer;
271: SlepcConvMon ctx;
273: PetscFunctionBegin;
276: ctx = (SlepcConvMon)vf->data;
277: if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)nep)->prefix));
278: if (its==1) ctx->oldnconv = 0;
279: if (ctx->oldnconv!=nconv) {
280: PetscCall(PetscViewerPushFormat(viewer,vf->format));
281: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
282: for (i=ctx->oldnconv;i<nconv;i++) {
283: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP converged value (error) #%" PetscInt_FMT,its,i));
284: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
285: #if defined(PETSC_USE_COMPLEX)
286: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i])));
287: #else
288: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]));
289: if (eigi[i]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]));
290: #endif
291: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
292: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
293: }
294: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
295: PetscCall(PetscViewerPopFormat(viewer));
296: ctx->oldnconv = nconv;
297: }
298: PetscFunctionReturn(PETSC_SUCCESS);
299: }
301: PetscErrorCode NEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
302: {
303: SlepcConvMon mctx;
305: PetscFunctionBegin;
306: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
307: PetscCall(PetscNew(&mctx));
308: mctx->ctx = ctx;
309: (*vf)->data = (void*)mctx;
310: PetscFunctionReturn(PETSC_SUCCESS);
311: }
313: PetscErrorCode NEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
314: {
315: PetscFunctionBegin;
316: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
317: PetscCall(PetscFree((*vf)->data));
318: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
319: PetscCall(PetscFree(*vf));
320: PetscFunctionReturn(PETSC_SUCCESS);
321: }
323: /*@C
324: NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
325: approximation at each iteration of the nonlinear eigensolver.
327: Collective
329: Input Parameters:
330: + nep - nonlinear eigensolver context
331: . its - iteration number
332: . nconv - number of converged eigenpairs so far
333: . eigr - real part of the eigenvalues
334: . eigi - imaginary part of the eigenvalues
335: . errest - error estimates
336: . nest - number of error estimates to display
337: - vf - viewer and format for monitoring
339: Options Database Key:
340: . -nep_monitor draw::draw_lg - activates NEPMonitorFirstDrawLG()
342: Level: intermediate
344: .seealso: NEPMonitorSet()
345: @*/
346: PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
347: {
348: PetscViewer viewer = vf->viewer;
349: PetscDrawLG lg;
350: PetscReal x,y;
352: PetscFunctionBegin;
355: PetscCall(PetscViewerPushFormat(viewer,vf->format));
356: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
357: if (its==1) {
358: PetscCall(PetscDrawLGReset(lg));
359: PetscCall(PetscDrawLGSetDimension(lg,1));
360: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
361: }
362: if (nconv<nest) {
363: x = (PetscReal)its;
364: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
365: else y = 0.0;
366: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
367: if (its <= 20 || !(its % 5) || nep->reason) {
368: PetscCall(PetscDrawLGDraw(lg));
369: PetscCall(PetscDrawLGSave(lg));
370: }
371: }
372: PetscCall(PetscViewerPopFormat(viewer));
373: PetscFunctionReturn(PETSC_SUCCESS);
374: }
376: /*@C
377: NEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
379: Collective
381: Input Parameters:
382: + viewer - the viewer
383: . format - the viewer format
384: - ctx - an optional user context
386: Output Parameter:
387: . vf - the viewer and format context
389: Level: intermediate
391: .seealso: NEPMonitorSet()
392: @*/
393: PetscErrorCode NEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
394: {
395: PetscFunctionBegin;
396: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
397: (*vf)->data = ctx;
398: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
399: PetscFunctionReturn(PETSC_SUCCESS);
400: }
402: /*@C
403: NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
404: approximations at each iteration of the nonlinear eigensolver.
406: Collective
408: Input Parameters:
409: + nep - nonlinear eigensolver context
410: . its - iteration number
411: . nconv - number of converged eigenpairs so far
412: . eigr - real part of the eigenvalues
413: . eigi - imaginary part of the eigenvalues
414: . errest - error estimates
415: . nest - number of error estimates to display
416: - vf - viewer and format for monitoring
418: Options Database Key:
419: . -nep_monitor_all draw::draw_lg - activates NEPMonitorAllDrawLG()
421: Level: intermediate
423: .seealso: NEPMonitorSet()
424: @*/
425: PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
426: {
427: PetscViewer viewer = vf->viewer;
428: PetscDrawLG lg;
429: PetscInt i,n = PetscMin(nep->nev,255);
430: PetscReal *x,*y;
432: PetscFunctionBegin;
435: PetscCall(PetscViewerPushFormat(viewer,vf->format));
436: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
437: if (its==1) {
438: PetscCall(PetscDrawLGReset(lg));
439: PetscCall(PetscDrawLGSetDimension(lg,n));
440: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
441: }
442: PetscCall(PetscMalloc2(n,&x,n,&y));
443: for (i=0;i<n;i++) {
444: x[i] = (PetscReal)its;
445: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
446: else y[i] = 0.0;
447: }
448: PetscCall(PetscDrawLGAddPoint(lg,x,y));
449: if (its <= 20 || !(its % 5) || nep->reason) {
450: PetscCall(PetscDrawLGDraw(lg));
451: PetscCall(PetscDrawLGSave(lg));
452: }
453: PetscCall(PetscFree2(x,y));
454: PetscCall(PetscViewerPopFormat(viewer));
455: PetscFunctionReturn(PETSC_SUCCESS);
456: }
458: /*@C
459: NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
461: Collective
463: Input Parameters:
464: + viewer - the viewer
465: . format - the viewer format
466: - ctx - an optional user context
468: Output Parameter:
469: . vf - the viewer and format context
471: Level: intermediate
473: .seealso: NEPMonitorSet()
474: @*/
475: PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
476: {
477: PetscFunctionBegin;
478: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
479: (*vf)->data = ctx;
480: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
481: PetscFunctionReturn(PETSC_SUCCESS);
482: }
484: /*@C
485: NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
486: at each iteration of the nonlinear eigensolver.
488: Collective
490: Input Parameters:
491: + nep - nonlinear eigensolver context
492: . its - iteration number
493: . nconv - number of converged eigenpairs so far
494: . eigr - real part of the eigenvalues
495: . eigi - imaginary part of the eigenvalues
496: . errest - error estimates
497: . nest - number of error estimates to display
498: - vf - viewer and format for monitoring
500: Options Database Key:
501: . -nep_monitor_conv draw::draw_lg - activates NEPMonitorConvergedDrawLG()
503: Level: intermediate
505: .seealso: NEPMonitorSet()
506: @*/
507: PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
508: {
509: PetscViewer viewer = vf->viewer;
510: PetscDrawLG lg;
511: PetscReal x,y;
513: PetscFunctionBegin;
516: PetscCall(PetscViewerPushFormat(viewer,vf->format));
517: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
518: if (its==1) {
519: PetscCall(PetscDrawLGReset(lg));
520: PetscCall(PetscDrawLGSetDimension(lg,1));
521: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,nep->nev));
522: }
523: x = (PetscReal)its;
524: y = (PetscReal)nep->nconv;
525: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
526: if (its <= 20 || !(its % 5) || nep->reason) {
527: PetscCall(PetscDrawLGDraw(lg));
528: PetscCall(PetscDrawLGSave(lg));
529: }
530: PetscCall(PetscViewerPopFormat(viewer));
531: PetscFunctionReturn(PETSC_SUCCESS);
532: }
534: /*@C
535: NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
537: Collective
539: Input Parameters:
540: + viewer - the viewer
541: . format - the viewer format
542: - ctx - an optional user context
544: Output Parameter:
545: . vf - the viewer and format context
547: Level: intermediate
549: .seealso: NEPMonitorSet()
550: @*/
551: PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
552: {
553: SlepcConvMon mctx;
555: PetscFunctionBegin;
556: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
557: PetscCall(PetscNew(&mctx));
558: mctx->ctx = ctx;
559: (*vf)->data = (void*)mctx;
560: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
561: PetscFunctionReturn(PETSC_SUCCESS);
562: }