1 #include <petscconvest.h> /*I "petscconvest.h" I*/
2 #include <petscdmplex.h>
3 #include <petscds.h>
4
5 #include <petsc/private/petscconvestimpl.h>
6
zero_private(PetscInt dim,PetscReal time,const PetscReal x[],PetscInt Nc,PetscScalar * u,void * ctx)7 static PetscErrorCode zero_private(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
8 {
9 PetscInt c;
10 for (c = 0; c < Nc; ++c) u[c] = 0.0;
11 return 0;
12 }
13
14 /*@
15 PetscConvEstDestroy - Destroys a PetscConvEst object
16
17 Collective on PetscConvEst
18
19 Input Parameter:
20 . ce - The PetscConvEst object
21
22 Level: beginner
23
24 .seealso: PetscConvEstCreate(), PetscConvEstGetConvRate()
25 @*/
PetscConvEstDestroy(PetscConvEst * ce)26 PetscErrorCode PetscConvEstDestroy(PetscConvEst *ce)
27 {
28 PetscErrorCode ierr;
29
30 PetscFunctionBegin;
31 if (!*ce) PetscFunctionReturn(0);
32 PetscValidHeaderSpecific((*ce),PETSC_OBJECT_CLASSID,1);
33 if (--((PetscObject)(*ce))->refct > 0) {
34 *ce = NULL;
35 PetscFunctionReturn(0);
36 }
37 ierr = PetscFree3((*ce)->initGuess, (*ce)->exactSol, (*ce)->ctxs);CHKERRQ(ierr);
38 ierr = PetscFree((*ce)->errors);CHKERRQ(ierr);
39 ierr = PetscHeaderDestroy(ce);CHKERRQ(ierr);
40 PetscFunctionReturn(0);
41 }
42
43 /*@
44 PetscConvEstSetFromOptions - Sets a PetscConvEst object from options
45
46 Collective on PetscConvEst
47
48 Input Parameters:
49 . ce - The PetscConvEst object
50
51 Level: beginner
52
53 .seealso: PetscConvEstCreate(), PetscConvEstGetConvRate()
54 @*/
PetscConvEstSetFromOptions(PetscConvEst ce)55 PetscErrorCode PetscConvEstSetFromOptions(PetscConvEst ce)
56 {
57 PetscErrorCode ierr;
58
59 PetscFunctionBegin;
60 ierr = PetscOptionsBegin(PetscObjectComm((PetscObject) ce), "", "Convergence Estimator Options", "PetscConvEst");CHKERRQ(ierr);
61 ierr = PetscOptionsInt("-convest_num_refine", "The number of refinements for the convergence check", "PetscConvEst", ce->Nr, &ce->Nr, NULL);CHKERRQ(ierr);
62 ierr = PetscOptionsReal("-convest_refine_factor", "The increase in resolution in each dimension", "PetscConvEst", ce->r, &ce->r, NULL);CHKERRQ(ierr);
63 ierr = PetscOptionsBool("-convest_monitor", "Monitor the error for each convergence check", "PetscConvEst", ce->monitor, &ce->monitor, NULL);CHKERRQ(ierr);
64 ierr = PetscOptionsEnd();
65 PetscFunctionReturn(0);
66 }
67
68 /*@
69 PetscConvEstView - Views a PetscConvEst object
70
71 Collective on PetscConvEst
72
73 Input Parameters:
74 + ce - The PetscConvEst object
75 - viewer - The PetscViewer object
76
77 Level: beginner
78
79 .seealso: PetscConvEstCreate(), PetscConvEstGetConvRate()
80 @*/
PetscConvEstView(PetscConvEst ce,PetscViewer viewer)81 PetscErrorCode PetscConvEstView(PetscConvEst ce, PetscViewer viewer)
82 {
83 PetscErrorCode ierr;
84
85 PetscFunctionBegin;
86 ierr = PetscObjectPrintClassNamePrefixType((PetscObject) ce, viewer);CHKERRQ(ierr);
87 ierr = PetscViewerASCIIPrintf(viewer, "ConvEst with %D levels\n", ce->Nr+1);CHKERRQ(ierr);
88 PetscFunctionReturn(0);
89 }
90
91 /*@
92 PetscConvEstGetSolver - Gets the solver used to produce discrete solutions
93
94 Not collective
95
96 Input Parameter:
97 . ce - The PetscConvEst object
98
99 Output Parameter:
100 . solver - The solver
101
102 Level: intermediate
103
104 .seealso: PetscConvEstSetSolver(), PetscConvEstCreate(), PetscConvEstGetConvRate()
105 @*/
PetscConvEstGetSolver(PetscConvEst ce,PetscObject * solver)106 PetscErrorCode PetscConvEstGetSolver(PetscConvEst ce, PetscObject *solver)
107 {
108 PetscFunctionBegin;
109 PetscValidHeaderSpecific(ce, PETSC_OBJECT_CLASSID, 1);
110 PetscValidPointer(solver, 2);
111 *solver = ce->solver;
112 PetscFunctionReturn(0);
113 }
114
115 /*@
116 PetscConvEstSetSolver - Sets the solver used to produce discrete solutions
117
118 Not collective
119
120 Input Parameters:
121 + ce - The PetscConvEst object
122 - solver - The solver
123
124 Level: intermediate
125
126 Note: The solver MUST have an attached DM/DS, so that we know the exact solution
127
128 .seealso: PetscConvEstGetSNES(), PetscConvEstCreate(), PetscConvEstGetConvRate()
129 @*/
PetscConvEstSetSolver(PetscConvEst ce,PetscObject solver)130 PetscErrorCode PetscConvEstSetSolver(PetscConvEst ce, PetscObject solver)
131 {
132 PetscErrorCode ierr;
133
134 PetscFunctionBegin;
135 PetscValidHeaderSpecific(ce, PETSC_OBJECT_CLASSID, 1);
136 PetscValidHeader(solver, 2);
137 ce->solver = solver;
138 ierr = (*ce->ops->setsolver)(ce, solver);CHKERRQ(ierr);
139 PetscFunctionReturn(0);
140 }
141
142 /*@
143 PetscConvEstSetUp - After the solver is specified, we create structures for estimating convergence
144
145 Collective on PetscConvEst
146
147 Input Parameters:
148 . ce - The PetscConvEst object
149
150 Level: beginner
151
152 .seealso: PetscConvEstCreate(), PetscConvEstGetConvRate()
153 @*/
PetscConvEstSetUp(PetscConvEst ce)154 PetscErrorCode PetscConvEstSetUp(PetscConvEst ce)
155 {
156 PetscInt Nf, f, Nds, s;
157 PetscErrorCode ierr;
158
159 PetscFunctionBegin;
160 ierr = DMGetNumFields(ce->idm, &Nf);CHKERRQ(ierr);
161 ce->Nf = PetscMax(Nf, 1);
162 ierr = PetscMalloc1((ce->Nr+1)*ce->Nf, &ce->errors);CHKERRQ(ierr);
163 ierr = PetscCalloc3(ce->Nf, &ce->initGuess, ce->Nf, &ce->exactSol, ce->Nf, &ce->ctxs);CHKERRQ(ierr);
164 for (f = 0; f < Nf; ++f) ce->initGuess[f] = zero_private;
165 ierr = DMGetNumDS(ce->idm, &Nds);CHKERRQ(ierr);
166 for (s = 0; s < Nds; ++s) {
167 PetscDS ds;
168 DMLabel label;
169 IS fieldIS;
170 const PetscInt *fields;
171 PetscInt dsNf;
172
173 ierr = DMGetRegionNumDS(ce->idm, s, &label, &fieldIS, &ds);CHKERRQ(ierr);
174 ierr = PetscDSGetNumFields(ds, &dsNf);CHKERRQ(ierr);
175 if (fieldIS) {ierr = ISGetIndices(fieldIS, &fields);CHKERRQ(ierr);}
176 for (f = 0; f < dsNf; ++f) {
177 const PetscInt field = fields[f];
178 ierr = PetscDSGetExactSolution(ds, field, &ce->exactSol[field], &ce->ctxs[field]);CHKERRQ(ierr);
179 }
180 if (fieldIS) {ierr = ISRestoreIndices(fieldIS, &fields);CHKERRQ(ierr);}
181 }
182 for (f = 0; f < Nf; ++f) {
183 if (!ce->exactSol[f]) SETERRQ1(PetscObjectComm((PetscObject) ce), PETSC_ERR_ARG_WRONG, "DS must contain exact solution functions in order to estimate convergence, missing for field %D", f);
184 }
185 PetscFunctionReturn(0);
186 }
187
PetscConvEstComputeInitialGuess(PetscConvEst ce,PetscInt r,DM dm,Vec u)188 PetscErrorCode PetscConvEstComputeInitialGuess(PetscConvEst ce, PetscInt r, DM dm, Vec u)
189 {
190 PetscErrorCode ierr;
191
192 PetscFunctionBegin;
193 PetscValidHeaderSpecific(ce, PETSC_OBJECT_CLASSID, 1);
194 if (dm) PetscValidHeaderSpecific(dm, DM_CLASSID, 3);
195 PetscValidHeaderSpecific(u, VEC_CLASSID, 4);
196 ierr = (*ce->ops->initguess)(ce, r, dm, u);CHKERRQ(ierr);
197 PetscFunctionReturn(0);
198 }
199
PetscConvEstComputeError(PetscConvEst ce,PetscInt r,DM dm,Vec u,PetscReal errors[])200 PetscErrorCode PetscConvEstComputeError(PetscConvEst ce, PetscInt r, DM dm, Vec u, PetscReal errors[])
201 {
202 PetscErrorCode ierr;
203
204 PetscFunctionBegin;
205 PetscValidHeaderSpecific(ce, PETSC_OBJECT_CLASSID, 1);
206 if (dm) PetscValidHeaderSpecific(dm, DM_CLASSID, 3);
207 PetscValidHeaderSpecific(u, VEC_CLASSID, 4);
208 PetscValidRealPointer(errors, 5);
209 ierr = (*ce->ops->computeerror)(ce, r, dm, u, errors);CHKERRQ(ierr);
210 PetscFunctionReturn(0);
211 }
212
213 /*@
214 PetscConvEstMonitorDefault - Monitors the convergence estimation loop
215
216 Collective on PetscConvEst
217
218 Input Parameter:
219 + ce - The PetscConvEst object
220 - r - The refinement level
221
222 Options database keys:
223 . -convest_monitor - Activate the monitor
224
225 Level: intermediate
226
227 .seealso: PetscConvEstCreate(), PetscConvEstGetConvRate(), SNESSolve(), TSSolve()
228 @*/
PetscConvEstMonitorDefault(PetscConvEst ce,PetscInt r)229 PetscErrorCode PetscConvEstMonitorDefault(PetscConvEst ce, PetscInt r)
230 {
231 MPI_Comm comm;
232 PetscInt f;
233 PetscErrorCode ierr;
234
235 PetscFunctionBegin;
236 if (ce->monitor) {
237 PetscReal *errors = &ce->errors[r*ce->Nf];
238
239 ierr = PetscObjectGetComm((PetscObject) ce, &comm);CHKERRQ(ierr);
240 ierr = PetscPrintf(comm, "L_2 Error: ");CHKERRQ(ierr);
241 if (ce->Nf > 1) {ierr = PetscPrintf(comm, "[");CHKERRQ(ierr);}
242 for (f = 0; f < ce->Nf; ++f) {
243 if (f > 0) {ierr = PetscPrintf(comm, ", ");CHKERRQ(ierr);}
244 if (errors[f] < 1.0e-11) {ierr = PetscPrintf(comm, "< 1e-11");CHKERRQ(ierr);}
245 else {ierr = PetscPrintf(comm, "%g", (double) errors[f]);CHKERRQ(ierr);}
246 }
247 if (ce->Nf > 1) {ierr = PetscPrintf(comm, "]");CHKERRQ(ierr);}
248 ierr = PetscPrintf(comm, "\n");CHKERRQ(ierr);
249 }
250 PetscFunctionReturn(0);
251 }
252
PetscConvEstSetSNES_Private(PetscConvEst ce,PetscObject solver)253 static PetscErrorCode PetscConvEstSetSNES_Private(PetscConvEst ce, PetscObject solver)
254 {
255 PetscClassId id;
256 PetscErrorCode ierr;
257
258 PetscFunctionBegin;
259 ierr = PetscObjectGetClassId(ce->solver, &id);CHKERRQ(ierr);
260 if (id != SNES_CLASSID) SETERRQ(PetscObjectComm((PetscObject) ce), PETSC_ERR_ARG_WRONG, "Solver was not a SNES");
261 ierr = SNESGetDM((SNES) ce->solver, &ce->idm);CHKERRQ(ierr);
262 PetscFunctionReturn(0);
263 }
264
PetscConvEstInitGuessSNES_Private(PetscConvEst ce,PetscInt r,DM dm,Vec u)265 static PetscErrorCode PetscConvEstInitGuessSNES_Private(PetscConvEst ce, PetscInt r, DM dm, Vec u)
266 {
267 PetscErrorCode ierr;
268
269 PetscFunctionBegin;
270 ierr = DMProjectFunction(dm, 0.0, ce->initGuess, ce->ctxs, INSERT_VALUES, u);CHKERRQ(ierr);
271 PetscFunctionReturn(0);
272 }
273
PetscConvEstComputeErrorSNES_Private(PetscConvEst ce,PetscInt r,DM dm,Vec u,PetscReal errors[])274 static PetscErrorCode PetscConvEstComputeErrorSNES_Private(PetscConvEst ce, PetscInt r, DM dm, Vec u, PetscReal errors[])
275 {
276 PetscErrorCode ierr;
277
278 PetscFunctionBegin;
279 ierr = DMComputeL2FieldDiff(dm, 0.0, ce->exactSol, ce->ctxs, u, errors);CHKERRQ(ierr);
280 PetscFunctionReturn(0);
281 }
282
PetscConvEstGetConvRateSNES_Private(PetscConvEst ce,PetscReal alpha[])283 static PetscErrorCode PetscConvEstGetConvRateSNES_Private(PetscConvEst ce, PetscReal alpha[])
284 {
285 SNES snes = (SNES) ce->solver;
286 DM *dm;
287 PetscObject disc;
288 PetscReal *x, *y, slope, intercept;
289 PetscInt *dof, Nr = ce->Nr, r, f, dim, oldlevel, oldnlev;
290 void *ctx;
291 PetscErrorCode ierr;
292
293 PetscFunctionBegin;
294 if (ce->r != 2.0) SETERRQ1(PetscObjectComm((PetscObject) ce), PETSC_ERR_SUP, "Only refinement factor 2 is currently supported (not %g)", (double) ce->r);
295 ierr = DMGetDimension(ce->idm, &dim);CHKERRQ(ierr);
296 ierr = DMGetApplicationContext(ce->idm, &ctx);CHKERRQ(ierr);
297 ierr = DMPlexSetRefinementUniform(ce->idm, PETSC_TRUE);CHKERRQ(ierr);
298 ierr = DMGetRefineLevel(ce->idm, &oldlevel);CHKERRQ(ierr);
299 ierr = PetscMalloc2((Nr+1), &dm, (Nr+1)*ce->Nf, &dof);CHKERRQ(ierr);
300 /* Loop over meshes */
301 dm[0] = ce->idm;
302 for (r = 0; r <= Nr; ++r) {
303 Vec u;
304 #if defined(PETSC_USE_LOG)
305 PetscLogStage stage;
306 #endif
307 char stageName[PETSC_MAX_PATH_LEN];
308 const char *dmname, *uname;
309
310 ierr = PetscSNPrintf(stageName, PETSC_MAX_PATH_LEN-1, "ConvEst Refinement Level %D", r);CHKERRQ(ierr);
311 ierr = PetscLogStageRegister(stageName, &stage);CHKERRQ(ierr);
312 ierr = PetscLogStagePush(stage);CHKERRQ(ierr);
313 if (r > 0) {
314 ierr = DMRefine(dm[r-1], MPI_COMM_NULL, &dm[r]);CHKERRQ(ierr);
315 ierr = DMSetCoarseDM(dm[r], dm[r-1]);CHKERRQ(ierr);
316 ierr = DMCopyTransform(ce->idm, dm[r]);CHKERRQ(ierr);
317 ierr = PetscObjectGetName((PetscObject) dm[r-1], &dmname);CHKERRQ(ierr);
318 ierr = PetscObjectSetName((PetscObject) dm[r], dmname);CHKERRQ(ierr);
319 for (f = 0; f <= ce->Nf; ++f) {
320 PetscErrorCode (*nspconstr)(DM, PetscInt, PetscInt, MatNullSpace *);
321
322 ierr = DMGetNullSpaceConstructor(dm[r-1], f, &nspconstr);CHKERRQ(ierr);
323 ierr = DMSetNullSpaceConstructor(dm[r], f, nspconstr);CHKERRQ(ierr);
324 }
325 }
326 ierr = DMViewFromOptions(dm[r], NULL, "-conv_dm_view");CHKERRQ(ierr);
327 /* Create solution */
328 ierr = DMCreateGlobalVector(dm[r], &u);CHKERRQ(ierr);
329 ierr = DMGetField(dm[r], 0, NULL, &disc);CHKERRQ(ierr);
330 ierr = PetscObjectGetName(disc, &uname);CHKERRQ(ierr);
331 ierr = PetscObjectSetName((PetscObject) u, uname);CHKERRQ(ierr);
332 /* Setup solver */
333 ierr = SNESReset(snes);CHKERRQ(ierr);
334 ierr = SNESSetDM(snes, dm[r]);CHKERRQ(ierr);
335 ierr = DMPlexSetSNESLocalFEM(dm[r], ctx, ctx, ctx);CHKERRQ(ierr);
336 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr);
337 /* Create initial guess */
338 ierr = PetscConvEstComputeInitialGuess(ce, r, dm[r], u);CHKERRQ(ierr);
339 ierr = SNESSolve(snes, NULL, u);CHKERRQ(ierr);
340 ierr = PetscLogEventBegin(ce->event, ce, 0, 0, 0);CHKERRQ(ierr);
341 ierr = PetscConvEstComputeError(ce, r, dm[r], u, &ce->errors[r*ce->Nf]);CHKERRQ(ierr);
342 ierr = PetscLogEventEnd(ce->event, ce, 0, 0, 0);CHKERRQ(ierr);
343 for (f = 0; f < ce->Nf; ++f) {
344 PetscSection s, fs;
345 PetscInt lsize;
346
347 /* Could use DMGetOutputDM() to add in Dirichlet dofs */
348 ierr = DMGetLocalSection(dm[r], &s);CHKERRQ(ierr);
349 ierr = PetscSectionGetField(s, f, &fs);CHKERRQ(ierr);
350 ierr = PetscSectionGetConstrainedStorageSize(fs, &lsize);CHKERRQ(ierr);
351 ierr = MPI_Allreduce(&lsize, &dof[r*ce->Nf+f], 1, MPIU_INT, MPI_SUM, PetscObjectComm((PetscObject) snes));CHKERRQ(ierr);
352 ierr = PetscLogEventSetDof(ce->event, f, dof[r*ce->Nf+f]);CHKERRQ(ierr);
353 ierr = PetscLogEventSetError(ce->event, f, ce->errors[r*ce->Nf+f]);CHKERRQ(ierr);
354 }
355 /* Monitor */
356 ierr = PetscConvEstMonitorDefault(ce, r);CHKERRQ(ierr);
357 if (!r) {
358 /* PCReset() does not wipe out the level structure */
359 KSP ksp;
360 PC pc;
361
362 ierr = SNESGetKSP(snes, &ksp);CHKERRQ(ierr);
363 ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);
364 ierr = PCMGGetLevels(pc, &oldnlev);CHKERRQ(ierr);
365 }
366 /* Cleanup */
367 ierr = VecDestroy(&u);CHKERRQ(ierr);
368 ierr = PetscLogStagePop();CHKERRQ(ierr);
369 }
370 for (r = 1; r <= Nr; ++r) {
371 ierr = DMDestroy(&dm[r]);CHKERRQ(ierr);
372 }
373 /* Fit convergence rate */
374 ierr = PetscMalloc2(Nr+1, &x, Nr+1, &y);CHKERRQ(ierr);
375 for (f = 0; f < ce->Nf; ++f) {
376 for (r = 0; r <= Nr; ++r) {
377 x[r] = PetscLog10Real(dof[r*ce->Nf+f]);
378 y[r] = PetscLog10Real(ce->errors[r*ce->Nf+f]);
379 }
380 ierr = PetscLinearRegression(Nr+1, x, y, &slope, &intercept);CHKERRQ(ierr);
381 /* Since h^{-dim} = N, lg err = s lg N + b = -s dim lg h + b */
382 alpha[f] = -slope * dim;
383 }
384 ierr = PetscFree2(x, y);CHKERRQ(ierr);
385 ierr = PetscFree2(dm, dof);CHKERRQ(ierr);
386 /* Restore solver */
387 ierr = SNESReset(snes);CHKERRQ(ierr);
388 {
389 /* PCReset() does not wipe out the level structure */
390 KSP ksp;
391 PC pc;
392
393 ierr = SNESGetKSP(snes, &ksp);CHKERRQ(ierr);
394 ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);
395 ierr = PCMGSetLevels(pc, oldnlev, NULL);CHKERRQ(ierr);
396 ierr = DMSetRefineLevel(ce->idm, oldlevel);CHKERRQ(ierr); /* The damn DMCoarsen() calls in PCMG can reset this */
397 }
398 ierr = SNESSetDM(snes, ce->idm);CHKERRQ(ierr);
399 ierr = DMPlexSetSNESLocalFEM(ce->idm, ctx, ctx, ctx);CHKERRQ(ierr);
400 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr);
401 PetscFunctionReturn(0);
402 }
403
404 /*@
405 PetscConvEstGetConvRate - Returns an estimate of the convergence rate for the discretization
406
407 Not collective
408
409 Input Parameter:
410 . ce - The PetscConvEst object
411
412 Output Parameter:
413 . alpha - The convergence rate for each field
414
415 Note: The convergence rate alpha is defined by
416 $ || u_\Delta - u_exact || < C \Delta^alpha
417 where u_\Delta is the discrete solution, and Delta is a measure of the discretization size. We usually use h for the
418 spatial resolution and \Delta t for the temporal resolution.
419
420 We solve a series of problems using increasing resolution (refined meshes or decreased timesteps), calculate an error
421 based upon the exact solution in the DS, and then fit the result to our model above using linear regression.
422
423 Options database keys:
424 + -snes_convergence_estimate : Execute convergence estimation inside SNESSolve() and print out the rate
425 - -ts_convergence_estimate : Execute convergence estimation inside TSSolve() and print out the rate
426
427 Level: intermediate
428
429 .seealso: PetscConvEstSetSolver(), PetscConvEstCreate(), PetscConvEstGetConvRate(), SNESSolve(), TSSolve()
430 @*/
PetscConvEstGetConvRate(PetscConvEst ce,PetscReal alpha[])431 PetscErrorCode PetscConvEstGetConvRate(PetscConvEst ce, PetscReal alpha[])
432 {
433 PetscInt f;
434 PetscErrorCode ierr;
435
436 PetscFunctionBegin;
437 if (ce->event < 0) {ierr = PetscLogEventRegister("ConvEst Error", PETSC_OBJECT_CLASSID, &ce->event);CHKERRQ(ierr);}
438 for (f = 0; f < ce->Nf; ++f) alpha[f] = 0.0;
439 ierr = (*ce->ops->getconvrate)(ce, alpha);CHKERRQ(ierr);
440 PetscFunctionReturn(0);
441 }
442
443 /*@
444 PetscConvEstRateView - Displays the convergence rate to a viewer
445
446 Collective on SNES
447
448 Parameter:
449 + snes - iterative context obtained from SNESCreate()
450 . alpha - the convergence rate for each field
451 - viewer - the viewer to display the reason
452
453 Options Database Keys:
454 . -snes_convergence_estimate - print the convergence rate
455
456 Level: developer
457
458 .seealso: PetscConvEstGetRate()
459 @*/
PetscConvEstRateView(PetscConvEst ce,const PetscReal alpha[],PetscViewer viewer)460 PetscErrorCode PetscConvEstRateView(PetscConvEst ce, const PetscReal alpha[], PetscViewer viewer)
461 {
462 PetscBool isAscii;
463 PetscErrorCode ierr;
464
465 PetscFunctionBegin;
466 ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &isAscii);CHKERRQ(ierr);
467 if (isAscii) {
468 PetscInt Nf = ce->Nf, f;
469
470 ierr = PetscViewerASCIIAddTab(viewer, ((PetscObject) ce)->tablevel);CHKERRQ(ierr);
471 ierr = PetscViewerASCIIPrintf(viewer, "L_2 convergence rate: ");CHKERRQ(ierr);
472 if (Nf > 1) {ierr = PetscViewerASCIIPrintf(viewer, "[");CHKERRQ(ierr);}
473 for (f = 0; f < Nf; ++f) {
474 if (f > 0) {ierr = PetscViewerASCIIPrintf(viewer, ", ");CHKERRQ(ierr);}
475 ierr = PetscViewerASCIIPrintf(viewer, "%#.2g", (double) alpha[f]);CHKERRQ(ierr);
476 }
477 if (Nf > 1) {ierr = PetscViewerASCIIPrintf(viewer, "]");CHKERRQ(ierr);}
478 ierr = PetscViewerASCIIPrintf(viewer, "\n");CHKERRQ(ierr);
479 ierr = PetscViewerASCIISubtractTab(viewer, ((PetscObject) ce)->tablevel);CHKERRQ(ierr);
480 }
481 PetscFunctionReturn(0);
482 }
483
484 /*@
485 PetscConvEstCreate - Create a PetscConvEst object
486
487 Collective
488
489 Input Parameter:
490 . comm - The communicator for the PetscConvEst object
491
492 Output Parameter:
493 . ce - The PetscConvEst object
494
495 Level: beginner
496
497 .seealso: PetscConvEstDestroy(), PetscConvEstGetConvRate()
498 @*/
PetscConvEstCreate(MPI_Comm comm,PetscConvEst * ce)499 PetscErrorCode PetscConvEstCreate(MPI_Comm comm, PetscConvEst *ce)
500 {
501 PetscErrorCode ierr;
502
503 PetscFunctionBegin;
504 PetscValidPointer(ce, 2);
505 ierr = PetscSysInitializePackage();CHKERRQ(ierr);
506 ierr = PetscHeaderCreate(*ce, PETSC_OBJECT_CLASSID, "PetscConvEst", "ConvergenceEstimator", "SNES", comm, PetscConvEstDestroy, PetscConvEstView);CHKERRQ(ierr);
507 (*ce)->monitor = PETSC_FALSE;
508 (*ce)->r = 2.0;
509 (*ce)->Nr = 4;
510 (*ce)->event = -1;
511 (*ce)->ops->setsolver = PetscConvEstSetSNES_Private;
512 (*ce)->ops->initguess = PetscConvEstInitGuessSNES_Private;
513 (*ce)->ops->computeerror = PetscConvEstComputeErrorSNES_Private;
514 (*ce)->ops->getconvrate = PetscConvEstGetConvRateSNES_Private;
515 PetscFunctionReturn(0);
516 }
517