Actual source code: epsmon.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: EPS routines related to monitors
12: */
14: #include <slepc/private/epsimpl.h>
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode EPSMonitor(EPS eps,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
22: PetscInt i,n = eps->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*eps->monitor[i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: EPSMonitorSet - 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: + eps - eigensolver context obtained from EPSCreate()
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(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx)
45: + eps - eigensolver context obtained from EPSCreate()
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 - relative error estimates for each eigenpair
51: . nest - number of error estimates
52: - mctx - optional monitoring context, as set by EPSMonitorSet()
54: Options Database Keys:
55: + -eps_monitor - print only the first error estimate
56: . -eps_monitor_all - print error estimates at each iteration
57: . -eps_monitor_conv - print the eigenvalue approximations only when
58: convergence has been reached
59: . -eps_monitor draw::draw_lg - sets line graph monitor for the first unconverged
60: approximate eigenvalue
61: . -eps_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
62: approximate eigenvalues
63: . -eps_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
64: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
65: a code by calls to EPSMonitorSet(), but does not cancel those set via
66: the options database.
68: Notes:
69: Several different monitoring routines may be set by calling
70: EPSMonitorSet() multiple times; all will be called in the
71: order in which they were set.
73: Level: intermediate
75: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorCancel()
76: @*/
77: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscCtxDestroyFn *monitordestroy)
78: {
79: PetscFunctionBegin;
81: PetscCheck(eps->numbermonitors<MAXEPSMONITORS,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
82: eps->monitor[eps->numbermonitors] = monitor;
83: eps->monitorcontext[eps->numbermonitors] = (void*)mctx;
84: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
85: PetscFunctionReturn(PETSC_SUCCESS);
86: }
88: /*@
89: EPSMonitorCancel - Clears all monitors for an EPS object.
91: Logically Collective
93: Input Parameters:
94: . eps - eigensolver context obtained from EPSCreate()
96: Options Database Key:
97: . -eps_monitor_cancel - Cancels all monitors that have been hardwired
98: into a code by calls to EPSMonitorSet(),
99: but does not cancel those set via the options database.
101: Level: intermediate
103: .seealso: EPSMonitorSet()
104: @*/
105: PetscErrorCode EPSMonitorCancel(EPS eps)
106: {
107: PetscInt i;
109: PetscFunctionBegin;
111: for (i=0; i<eps->numbermonitors; i++) {
112: if (eps->monitordestroy[i]) PetscCall((*eps->monitordestroy[i])(&eps->monitorcontext[i]));
113: }
114: eps->numbermonitors = 0;
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: /*@C
119: EPSGetMonitorContext - Gets the monitor context, as set by
120: EPSMonitorSet() for the FIRST monitor only.
122: Not Collective
124: Input Parameter:
125: . eps - eigensolver context obtained from EPSCreate()
127: Output Parameter:
128: . ctx - monitor context
130: Level: intermediate
132: .seealso: EPSMonitorSet()
133: @*/
134: PetscErrorCode EPSGetMonitorContext(EPS eps,void *ctx)
135: {
136: PetscFunctionBegin;
138: *(void**)ctx = eps->monitorcontext[0];
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: static inline PetscErrorCode EPSMonitorPrintEval(EPS eps,PetscViewer viewer,PetscScalar er,PetscScalar ei)
143: {
144: PetscFunctionBegin;
145: if (eps->problem_type==EPS_HEP || eps->problem_type==EPS_GHEP || eps->problem_type==EPS_BSE) PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)PetscRealPart(er)));
146: else {
147: #if defined(PETSC_USE_COMPLEX)
148: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
149: #else
150: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
151: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
152: #endif
153: }
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }
157: /*@C
158: EPSMonitorFirst - Print the first unconverged approximate value and
159: error estimate at each iteration of the eigensolver.
161: Collective
163: Input Parameters:
164: + eps - eigensolver context
165: . its - iteration number
166: . nconv - number of converged eigenpairs so far
167: . eigr - real part of the eigenvalues
168: . eigi - imaginary part of the eigenvalues
169: . errest - error estimates
170: . nest - number of error estimates to display
171: - vf - viewer and format for monitoring
173: Options Database Key:
174: . -eps_monitor - activates EPSMonitorFirst()
176: Level: intermediate
178: .seealso: EPSMonitorSet(), EPSMonitorAll(), EPSMonitorConverged()
179: @*/
180: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
181: {
182: PetscScalar er,ei;
183: PetscViewer viewer = vf->viewer;
185: PetscFunctionBegin;
188: if (its==1 && ((PetscObject)eps)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix));
189: if (nconv<nest) {
190: PetscCall(PetscViewerPushFormat(viewer,vf->format));
191: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
192: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " EPS nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
193: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
194: er = eigr[nconv]; ei = eigi[nconv];
195: PetscCall(STBackTransform(eps->st,1,&er,&ei));
196: PetscCall(EPSMonitorPrintEval(eps,viewer,er,ei));
197: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
198: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
199: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
200: PetscCall(PetscViewerPopFormat(viewer));
201: }
202: PetscFunctionReturn(PETSC_SUCCESS);
203: }
205: /*@C
206: EPSMonitorAll - Print the current approximate values and
207: error estimates at each iteration of the eigensolver.
209: Collective
211: Input Parameters:
212: + eps - eigensolver context
213: . its - iteration number
214: . nconv - number of converged eigenpairs so far
215: . eigr - real part of the eigenvalues
216: . eigi - imaginary part of the eigenvalues
217: . errest - error estimates
218: . nest - number of error estimates to display
219: - vf - viewer and format for monitoring
221: Options Database Key:
222: . -eps_monitor_all - activates EPSMonitorAll()
224: Level: intermediate
226: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
227: @*/
228: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
229: {
230: PetscInt i;
231: PetscScalar er,ei;
232: PetscViewer viewer = vf->viewer;
234: PetscFunctionBegin;
237: PetscCall(PetscViewerPushFormat(viewer,vf->format));
238: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
239: if (its==1 && ((PetscObject)eps)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix));
240: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " EPS nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
241: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
242: for (i=0;i<nest;i++) {
243: er = eigr[i]; ei = eigi[i];
244: PetscCall(STBackTransform(eps->st,1,&er,&ei));
245: PetscCall(EPSMonitorPrintEval(eps,viewer,er,ei));
246: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
247: }
248: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
249: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
250: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
251: PetscCall(PetscViewerPopFormat(viewer));
252: PetscFunctionReturn(PETSC_SUCCESS);
253: }
255: /*@C
256: EPSMonitorConverged - Print the approximate values and
257: error estimates as they converge.
259: Collective
261: Input Parameters:
262: + eps - eigensolver context
263: . its - iteration number
264: . nconv - number of converged eigenpairs so far
265: . eigr - real part of the eigenvalues
266: . eigi - imaginary part of the eigenvalues
267: . errest - error estimates
268: . nest - number of error estimates to display
269: - vf - viewer and format for monitoring
271: Options Database Key:
272: . -eps_monitor_conv - activates EPSMonitorConverged()
274: Level: intermediate
276: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
277: @*/
278: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
279: {
280: PetscInt i;
281: PetscScalar er,ei;
282: PetscViewer viewer = vf->viewer;
283: SlepcConvMon ctx;
285: PetscFunctionBegin;
288: ctx = (SlepcConvMon)vf->data;
289: if (its==1 && ((PetscObject)eps)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)eps)->prefix));
290: if (its==1) ctx->oldnconv = 0;
291: if (ctx->oldnconv!=nconv) {
292: PetscCall(PetscViewerPushFormat(viewer,vf->format));
293: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
294: for (i=ctx->oldnconv;i<nconv;i++) {
295: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " EPS converged value (error) #%" PetscInt_FMT,its,i));
296: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
297: er = eigr[i]; ei = eigi[i];
298: PetscCall(STBackTransform(eps->st,1,&er,&ei));
299: PetscCall(EPSMonitorPrintEval(eps,viewer,er,ei));
300: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
301: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
302: }
303: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
304: PetscCall(PetscViewerPopFormat(viewer));
305: ctx->oldnconv = nconv;
306: }
307: PetscFunctionReturn(PETSC_SUCCESS);
308: }
310: PetscErrorCode EPSMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
311: {
312: SlepcConvMon mctx;
314: PetscFunctionBegin;
315: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
316: PetscCall(PetscNew(&mctx));
317: mctx->ctx = ctx;
318: (*vf)->data = (void*)mctx;
319: PetscFunctionReturn(PETSC_SUCCESS);
320: }
322: PetscErrorCode EPSMonitorConvergedDestroy(PetscViewerAndFormat **vf)
323: {
324: PetscFunctionBegin;
325: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
326: PetscCall(PetscFree((*vf)->data));
327: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
328: PetscCall(PetscFree(*vf));
329: PetscFunctionReturn(PETSC_SUCCESS);
330: }
332: /*@C
333: EPSMonitorFirstDrawLG - Plots the error estimate of the first unconverged
334: approximation at each iteration of the eigensolver.
336: Collective
338: Input Parameters:
339: + eps - eigensolver context
340: . its - iteration number
341: . nconv - number of converged eigenpairs so far
342: . eigr - real part of the eigenvalues
343: . eigi - imaginary part of the eigenvalues
344: . errest - error estimates
345: . nest - number of error estimates to display
346: - vf - viewer and format for monitoring
348: Options Database Key:
349: . -eps_monitor draw::draw_lg - activates EPSMonitorFirstDrawLG()
351: Level: intermediate
353: .seealso: EPSMonitorSet()
354: @*/
355: PetscErrorCode EPSMonitorFirstDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
356: {
357: PetscViewer viewer = vf->viewer;
358: PetscDrawLG lg;
359: PetscReal x,y;
361: PetscFunctionBegin;
364: PetscCall(PetscViewerPushFormat(viewer,vf->format));
365: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
366: if (its==1) {
367: PetscCall(PetscDrawLGReset(lg));
368: PetscCall(PetscDrawLGSetDimension(lg,1));
369: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0));
370: }
371: if (nconv<nest) {
372: x = (PetscReal)its;
373: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
374: else y = 0.0;
375: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
376: if (its <= 20 || !(its % 5) || eps->reason) {
377: PetscCall(PetscDrawLGDraw(lg));
378: PetscCall(PetscDrawLGSave(lg));
379: }
380: }
381: PetscCall(PetscViewerPopFormat(viewer));
382: PetscFunctionReturn(PETSC_SUCCESS);
383: }
385: /*@C
386: EPSMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
388: Collective
390: Input Parameters:
391: + viewer - the viewer
392: . format - the viewer format
393: - ctx - an optional user context
395: Output Parameter:
396: . vf - the viewer and format context
398: Level: intermediate
400: .seealso: EPSMonitorSet()
401: @*/
402: PetscErrorCode EPSMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
403: {
404: PetscFunctionBegin;
405: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
406: (*vf)->data = ctx;
407: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
408: PetscFunctionReturn(PETSC_SUCCESS);
409: }
411: /*@C
412: EPSMonitorAllDrawLG - Plots the error estimate of all unconverged
413: approximations at each iteration of the eigensolver.
415: Collective
417: Input Parameters:
418: + eps - eigensolver context
419: . its - iteration number
420: . nconv - number of converged eigenpairs so far
421: . eigr - real part of the eigenvalues
422: . eigi - imaginary part of the eigenvalues
423: . errest - error estimates
424: . nest - number of error estimates to display
425: - vf - viewer and format for monitoring
427: Options Database Key:
428: . -eps_monitor_all draw::draw_lg - activates EPSMonitorAllDrawLG()
430: Level: intermediate
432: .seealso: EPSMonitorSet()
433: @*/
434: PetscErrorCode EPSMonitorAllDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
435: {
436: PetscViewer viewer = vf->viewer;
437: PetscDrawLG lg;
438: PetscInt i,n = PetscMin(eps->nev,255);
439: PetscReal *x,*y;
441: PetscFunctionBegin;
444: PetscCall(PetscViewerPushFormat(viewer,vf->format));
445: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
446: if (its==1) {
447: PetscCall(PetscDrawLGReset(lg));
448: PetscCall(PetscDrawLGSetDimension(lg,n));
449: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0));
450: }
451: PetscCall(PetscMalloc2(n,&x,n,&y));
452: for (i=0;i<n;i++) {
453: x[i] = (PetscReal)its;
454: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
455: else y[i] = 0.0;
456: }
457: PetscCall(PetscDrawLGAddPoint(lg,x,y));
458: if (its <= 20 || !(its % 5) || eps->reason) {
459: PetscCall(PetscDrawLGDraw(lg));
460: PetscCall(PetscDrawLGSave(lg));
461: }
462: PetscCall(PetscFree2(x,y));
463: PetscCall(PetscViewerPopFormat(viewer));
464: PetscFunctionReturn(PETSC_SUCCESS);
465: }
467: /*@C
468: EPSMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
470: Collective
472: Input Parameters:
473: + viewer - the viewer
474: . format - the viewer format
475: - ctx - an optional user context
477: Output Parameter:
478: . vf - the viewer and format context
480: Level: intermediate
482: .seealso: EPSMonitorSet()
483: @*/
484: PetscErrorCode EPSMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
485: {
486: PetscFunctionBegin;
487: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
488: (*vf)->data = ctx;
489: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
490: PetscFunctionReturn(PETSC_SUCCESS);
491: }
493: /*@C
494: EPSMonitorConvergedDrawLG - Plots the number of converged eigenvalues
495: at each iteration of the eigensolver.
497: Collective
499: Input Parameters:
500: + eps - eigensolver context
501: . its - iteration number
502: . nconv - number of converged eigenpairs so far
503: . eigr - real part of the eigenvalues
504: . eigi - imaginary part of the eigenvalues
505: . errest - error estimates
506: . nest - number of error estimates to display
507: - vf - viewer and format for monitoring
509: Options Database Key:
510: . -eps_monitor_conv draw::draw_lg - activates EPSMonitorConvergedDrawLG()
512: Level: intermediate
514: .seealso: EPSMonitorSet()
515: @*/
516: PetscErrorCode EPSMonitorConvergedDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
517: {
518: PetscViewer viewer = vf->viewer;
519: PetscDrawLG lg;
520: PetscReal x,y;
522: PetscFunctionBegin;
525: PetscCall(PetscViewerPushFormat(viewer,vf->format));
526: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
527: if (its==1) {
528: PetscCall(PetscDrawLGReset(lg));
529: PetscCall(PetscDrawLGSetDimension(lg,1));
530: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,eps->nev));
531: }
532: x = (PetscReal)its;
533: y = (PetscReal)eps->nconv;
534: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
535: if (its <= 20 || !(its % 5) || eps->reason) {
536: PetscCall(PetscDrawLGDraw(lg));
537: PetscCall(PetscDrawLGSave(lg));
538: }
539: PetscCall(PetscViewerPopFormat(viewer));
540: PetscFunctionReturn(PETSC_SUCCESS);
541: }
543: /*@C
544: EPSMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
546: Collective
548: Input Parameters:
549: + viewer - the viewer
550: . format - the viewer format
551: - ctx - an optional user context
553: Output Parameter:
554: . vf - the viewer and format context
556: Level: intermediate
558: .seealso: EPSMonitorSet()
559: @*/
560: PetscErrorCode EPSMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
561: {
562: SlepcConvMon mctx;
564: PetscFunctionBegin;
565: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
566: PetscCall(PetscNew(&mctx));
567: mctx->ctx = ctx;
568: (*vf)->data = (void*)mctx;
569: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
570: PetscFunctionReturn(PETSC_SUCCESS);
571: }