1 #include <petscconvest.h> /*I "petscconvest.h" I*/ 2 #include <petscdmplex.h> 3 #include <petscds.h> 4 5 #include <petsc/private/petscconvestimpl.h> 6 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 @*/ 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 @*/ 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 @*/ 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 @*/ 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 @*/ 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 @*/ 154 PetscErrorCode PetscConvEstSetUp(PetscConvEst ce) 155 { 156 PetscDS ds; 157 PetscInt Nf, f; 158 PetscErrorCode ierr; 159 160 PetscFunctionBegin; 161 ierr = DMGetDS(ce->idm, &ds);CHKERRQ(ierr); 162 ierr = PetscDSGetNumFields(ds, &Nf);CHKERRQ(ierr); 163 ce->Nf = PetscMax(Nf, 1); 164 ierr = PetscMalloc1((ce->Nr+1)*ce->Nf, &ce->errors);CHKERRQ(ierr); 165 ierr = PetscMalloc3(ce->Nf, &ce->initGuess, ce->Nf, &ce->exactSol, ce->Nf, &ce->ctxs);CHKERRQ(ierr); 166 for (f = 0; f < Nf; ++f) ce->initGuess[f] = zero_private; 167 for (f = 0; f < Nf; ++f) { 168 ierr = PetscDSGetExactSolution(ds, f, &ce->exactSol[f], &ce->ctxs[f]);CHKERRQ(ierr); 169 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); 170 } 171 PetscFunctionReturn(0); 172 } 173 174 PetscErrorCode PetscConvEstComputeInitialGuess(PetscConvEst ce, PetscInt r, DM dm, Vec u) 175 { 176 PetscErrorCode ierr; 177 178 PetscFunctionBegin; 179 PetscValidHeaderSpecific(ce, PETSC_OBJECT_CLASSID, 1); 180 if (dm) PetscValidHeaderSpecific(dm, DM_CLASSID, 3); 181 PetscValidHeaderSpecific(u, VEC_CLASSID, 4); 182 ierr = (*ce->ops->initguess)(ce, r, dm, u);CHKERRQ(ierr); 183 PetscFunctionReturn(0); 184 } 185 186 PetscErrorCode PetscConvEstComputeError(PetscConvEst ce, PetscInt r, DM dm, Vec u, PetscReal errors[]) 187 { 188 PetscErrorCode ierr; 189 190 PetscFunctionBegin; 191 PetscValidHeaderSpecific(ce, PETSC_OBJECT_CLASSID, 1); 192 if (dm) PetscValidHeaderSpecific(dm, DM_CLASSID, 3); 193 PetscValidHeaderSpecific(u, VEC_CLASSID, 4); 194 PetscValidRealPointer(errors, 5); 195 ierr = (*ce->ops->computeerror)(ce, r, dm, u, errors);CHKERRQ(ierr); 196 PetscFunctionReturn(0); 197 } 198 199 static PetscErrorCode PetscConvEstMonitor_Private(PetscConvEst ce, PetscInt r) 200 { 201 MPI_Comm comm; 202 PetscInt f; 203 PetscErrorCode ierr; 204 205 PetscFunctionBegin; 206 if (ce->monitor) { 207 PetscReal *errors = &ce->errors[r*ce->Nf]; 208 209 ierr = PetscObjectGetComm((PetscObject) ce, &comm);CHKERRQ(ierr); 210 ierr = PetscPrintf(comm, "L_2 Error: ");CHKERRQ(ierr); 211 if (ce->Nf > 1) {ierr = PetscPrintf(comm, "[");CHKERRQ(ierr);} 212 for (f = 0; f < ce->Nf; ++f) { 213 if (f > 0) {ierr = PetscPrintf(comm, ", ");CHKERRQ(ierr);} 214 if (errors[f] < 1.0e-11) {ierr = PetscPrintf(comm, "< 1e-11");CHKERRQ(ierr);} 215 else {ierr = PetscPrintf(comm, "%g", (double) errors[f]);CHKERRQ(ierr);} 216 } 217 if (ce->Nf > 1) {ierr = PetscPrintf(comm, "]");CHKERRQ(ierr);} 218 ierr = PetscPrintf(comm, "\n");CHKERRQ(ierr); 219 } 220 PetscFunctionReturn(0); 221 } 222 223 static PetscErrorCode PetscConvEstSetSNES_Private(PetscConvEst ce, PetscObject solver) 224 { 225 PetscClassId id; 226 PetscErrorCode ierr; 227 228 PetscFunctionBegin; 229 ierr = PetscObjectGetClassId(ce->solver, &id);CHKERRQ(ierr); 230 if (id != SNES_CLASSID) SETERRQ(PetscObjectComm((PetscObject) ce), PETSC_ERR_ARG_WRONG, "Solver was not a SNES"); 231 ierr = SNESGetDM((SNES) ce->solver, &ce->idm);CHKERRQ(ierr); 232 PetscFunctionReturn(0); 233 } 234 235 static PetscErrorCode PetscConvEstInitGuessSNES_Private(PetscConvEst ce, PetscInt r, DM dm, Vec u) 236 { 237 PetscErrorCode ierr; 238 239 PetscFunctionBegin; 240 ierr = DMProjectFunction(dm, 0.0, ce->initGuess, ce->ctxs, INSERT_VALUES, u);CHKERRQ(ierr); 241 PetscFunctionReturn(0); 242 } 243 244 static PetscErrorCode PetscConvEstComputeErrorSNES_Private(PetscConvEst ce, PetscInt r, DM dm, Vec u, PetscReal errors[]) 245 { 246 PetscErrorCode ierr; 247 248 PetscFunctionBegin; 249 ierr = DMComputeL2FieldDiff(dm, 0.0, ce->exactSol, ce->ctxs, u, errors);CHKERRQ(ierr); 250 PetscFunctionReturn(0); 251 } 252 253 static PetscErrorCode PetscConvEstGetConvRateSNES_Private(PetscConvEst ce, PetscReal alpha[]) 254 { 255 SNES snes = (SNES) ce->solver; 256 DM *dm; 257 PetscObject disc; 258 PetscReal *x, *y, slope, intercept; 259 PetscInt *dof, Nr = ce->Nr, r, f, dim, oldlevel, oldnlev; 260 void *ctx; 261 PetscErrorCode ierr; 262 263 PetscFunctionBegin; 264 if (ce->r != 2.0) SETERRQ1(PetscObjectComm((PetscObject) ce), PETSC_ERR_SUP, "Only refinement factor 2 is currently supported (not %g)", (double) ce->r); 265 ierr = DMGetDimension(ce->idm, &dim);CHKERRQ(ierr); 266 ierr = DMGetApplicationContext(ce->idm, &ctx);CHKERRQ(ierr); 267 ierr = DMPlexSetRefinementUniform(ce->idm, PETSC_TRUE);CHKERRQ(ierr); 268 ierr = DMGetRefineLevel(ce->idm, &oldlevel);CHKERRQ(ierr); 269 ierr = PetscMalloc2((Nr+1), &dm, (Nr+1)*ce->Nf, &dof);CHKERRQ(ierr); 270 /* Loop over meshes */ 271 dm[0] = ce->idm; 272 for (r = 0; r <= Nr; ++r) { 273 Vec u; 274 #if defined(PETSC_USE_LOG) 275 PetscLogStage stage; 276 #endif 277 char stageName[PETSC_MAX_PATH_LEN]; 278 const char *dmname, *uname; 279 280 ierr = PetscSNPrintf(stageName, PETSC_MAX_PATH_LEN-1, "ConvEst Refinement Level %D", r);CHKERRQ(ierr); 281 ierr = PetscLogStageRegister(stageName, &stage);CHKERRQ(ierr); 282 ierr = PetscLogStagePush(stage);CHKERRQ(ierr); 283 if (r > 0) { 284 ierr = DMRefine(dm[r-1], MPI_COMM_NULL, &dm[r]);CHKERRQ(ierr); 285 ierr = DMSetCoarseDM(dm[r], dm[r-1]);CHKERRQ(ierr); 286 ierr = DMCopyDisc(ce->idm, dm[r]);CHKERRQ(ierr); 287 ierr = DMCopyTransform(ce->idm, dm[r]);CHKERRQ(ierr); 288 ierr = PetscObjectGetName((PetscObject) dm[r-1], &dmname);CHKERRQ(ierr); 289 ierr = PetscObjectSetName((PetscObject) dm[r], dmname);CHKERRQ(ierr); 290 for (f = 0; f <= ce->Nf; ++f) { 291 PetscErrorCode (*nspconstr)(DM, PetscInt, MatNullSpace *); 292 293 ierr = DMGetNullSpaceConstructor(dm[r-1], f, &nspconstr);CHKERRQ(ierr); 294 ierr = DMSetNullSpaceConstructor(dm[r], f, nspconstr);CHKERRQ(ierr); 295 } 296 } 297 ierr = DMViewFromOptions(dm[r], NULL, "-conv_dm_view");CHKERRQ(ierr); 298 /* Create solution */ 299 ierr = DMCreateGlobalVector(dm[r], &u);CHKERRQ(ierr); 300 ierr = DMGetField(dm[r], 0, NULL, &disc);CHKERRQ(ierr); 301 ierr = PetscObjectGetName(disc, &uname);CHKERRQ(ierr); 302 ierr = PetscObjectSetName((PetscObject) u, uname);CHKERRQ(ierr); 303 /* Setup solver */ 304 ierr = SNESReset(snes);CHKERRQ(ierr); 305 ierr = SNESSetDM(snes, dm[r]);CHKERRQ(ierr); 306 ierr = DMPlexSetSNESLocalFEM(dm[r], ctx, ctx, ctx);CHKERRQ(ierr); 307 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr); 308 /* Create initial guess */ 309 ierr = PetscConvEstComputeInitialGuess(ce, r, dm[r], u);CHKERRQ(ierr); 310 ierr = SNESSolve(snes, NULL, u);CHKERRQ(ierr); 311 ierr = PetscLogEventBegin(ce->event, ce, 0, 0, 0);CHKERRQ(ierr); 312 ierr = PetscConvEstComputeError(ce, r, dm[r], u, &ce->errors[r*ce->Nf]);CHKERRQ(ierr); 313 ierr = PetscLogEventEnd(ce->event, ce, 0, 0, 0);CHKERRQ(ierr); 314 for (f = 0; f < ce->Nf; ++f) { 315 PetscSection s, fs; 316 PetscInt lsize; 317 318 /* Could use DMGetOutputDM() to add in Dirichlet dofs */ 319 ierr = DMGetLocalSection(dm[r], &s);CHKERRQ(ierr); 320 ierr = PetscSectionGetField(s, f, &fs);CHKERRQ(ierr); 321 ierr = PetscSectionGetConstrainedStorageSize(fs, &lsize);CHKERRQ(ierr); 322 ierr = MPI_Allreduce(&lsize, &dof[r*ce->Nf+f], 1, MPIU_INT, MPI_SUM, PetscObjectComm((PetscObject) snes));CHKERRQ(ierr); 323 ierr = PetscLogEventSetDof(ce->event, f, dof[r*ce->Nf+f]);CHKERRQ(ierr); 324 ierr = PetscLogEventSetError(ce->event, f, ce->errors[r*ce->Nf+f]);CHKERRQ(ierr); 325 } 326 /* Monitor */ 327 ierr = PetscConvEstMonitor_Private(ce, r);CHKERRQ(ierr); 328 if (!r) { 329 /* PCReset() does not wipe out the level structure */ 330 KSP ksp; 331 PC pc; 332 333 ierr = SNESGetKSP(snes, &ksp);CHKERRQ(ierr); 334 ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr); 335 ierr = PCMGGetLevels(pc, &oldnlev);CHKERRQ(ierr); 336 } 337 /* Cleanup */ 338 ierr = VecDestroy(&u);CHKERRQ(ierr); 339 ierr = PetscLogStagePop();CHKERRQ(ierr); 340 } 341 for (r = 1; r <= Nr; ++r) { 342 ierr = DMDestroy(&dm[r]);CHKERRQ(ierr); 343 } 344 /* Fit convergence rate */ 345 ierr = PetscMalloc2(Nr+1, &x, Nr+1, &y);CHKERRQ(ierr); 346 for (f = 0; f < ce->Nf; ++f) { 347 for (r = 0; r <= Nr; ++r) { 348 x[r] = PetscLog10Real(dof[r*ce->Nf+f]); 349 y[r] = PetscLog10Real(ce->errors[r*ce->Nf+f]); 350 } 351 ierr = PetscLinearRegression(Nr+1, x, y, &slope, &intercept);CHKERRQ(ierr); 352 /* Since h^{-dim} = N, lg err = s lg N + b = -s dim lg h + b */ 353 alpha[f] = -slope * dim; 354 } 355 ierr = PetscFree2(x, y);CHKERRQ(ierr); 356 ierr = PetscFree2(dm, dof);CHKERRQ(ierr); 357 /* Restore solver */ 358 ierr = SNESReset(snes);CHKERRQ(ierr); 359 { 360 /* PCReset() does not wipe out the level structure */ 361 KSP ksp; 362 PC pc; 363 364 ierr = SNESGetKSP(snes, &ksp);CHKERRQ(ierr); 365 ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr); 366 ierr = PCMGSetLevels(pc, oldnlev, NULL);CHKERRQ(ierr); 367 ierr = DMSetRefineLevel(ce->idm, oldlevel);CHKERRQ(ierr); /* The damn DMCoarsen() calls in PCMG can reset this */ 368 } 369 ierr = SNESSetDM(snes, ce->idm);CHKERRQ(ierr); 370 ierr = DMPlexSetSNESLocalFEM(ce->idm, ctx, ctx, ctx);CHKERRQ(ierr); 371 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr); 372 PetscFunctionReturn(0); 373 } 374 375 /*@ 376 PetscConvEstGetConvRate - Returns an estimate of the convergence rate for the discretization 377 378 Not collective 379 380 Input Parameter: 381 . ce - The PetscConvEst object 382 383 Output Parameter: 384 . alpha - The convergence rate for each field 385 386 Note: The convergence rate alpha is defined by 387 $ || u_\Delta - u_exact || < C \Delta^alpha 388 where u_\Delta is the discrete solution, and Delta is a measure of the discretization size. We usually use h for the 389 spatial resolution and \Delta t for the temporal resolution. 390 391 We solve a series of problems using increasing resolution (refined meshes or decreased timesteps), calculate an error 392 based upon the exact solution in the DS, and then fit the result to our model above using linear regression. 393 394 Options database keys: 395 + -snes_convergence_estimate : Execute convergence estimation inside SNESSolve() and print out the rate 396 - -ts_convergence_estimate : Execute convergence estimation inside TSSolve() and print out the rate 397 398 Level: intermediate 399 400 .seealso: PetscConvEstSetSolver(), PetscConvEstCreate(), PetscConvEstGetConvRate(), SNESSolve(), TSSolve() 401 @*/ 402 PetscErrorCode PetscConvEstGetConvRate(PetscConvEst ce, PetscReal alpha[]) 403 { 404 PetscInt f; 405 PetscErrorCode ierr; 406 407 PetscFunctionBegin; 408 if (ce->event < 0) {ierr = PetscLogEventRegister("ConvEst Error", PETSC_OBJECT_CLASSID, &ce->event);CHKERRQ(ierr);} 409 for (f = 0; f < ce->Nf; ++f) alpha[f] = 0.0; 410 ierr = (*ce->ops->getconvrate)(ce, alpha);CHKERRQ(ierr); 411 PetscFunctionReturn(0); 412 } 413 414 /*@ 415 PetscConvEstRateView - Displays the convergence rate to a viewer 416 417 Collective on SNES 418 419 Parameter: 420 + snes - iterative context obtained from SNESCreate() 421 . alpha - the convergence rate for each field 422 - viewer - the viewer to display the reason 423 424 Options Database Keys: 425 . -snes_convergence_estimate - print the convergence rate 426 427 Level: developer 428 429 .seealso: PetscConvEstGetRate() 430 @*/ 431 PetscErrorCode PetscConvEstRateView(PetscConvEst ce, const PetscReal alpha[], PetscViewer viewer) 432 { 433 PetscBool isAscii; 434 PetscErrorCode ierr; 435 436 PetscFunctionBegin; 437 ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &isAscii);CHKERRQ(ierr); 438 if (isAscii) { 439 PetscInt Nf = ce->Nf, f; 440 441 ierr = PetscViewerASCIIAddTab(viewer, ((PetscObject) ce)->tablevel);CHKERRQ(ierr); 442 ierr = PetscViewerASCIIPrintf(viewer, "L_2 convergence rate: ");CHKERRQ(ierr); 443 if (Nf > 1) {ierr = PetscViewerASCIIPrintf(viewer, "[");CHKERRQ(ierr);} 444 for (f = 0; f < Nf; ++f) { 445 if (f > 0) {ierr = PetscViewerASCIIPrintf(viewer, ", ");CHKERRQ(ierr);} 446 ierr = PetscViewerASCIIPrintf(viewer, "%#.2g", (double) alpha[f]);CHKERRQ(ierr); 447 } 448 if (Nf > 1) {ierr = PetscViewerASCIIPrintf(viewer, "]");CHKERRQ(ierr);} 449 ierr = PetscViewerASCIIPrintf(viewer, "\n");CHKERRQ(ierr); 450 ierr = PetscViewerASCIISubtractTab(viewer, ((PetscObject) ce)->tablevel);CHKERRQ(ierr); 451 } 452 PetscFunctionReturn(0); 453 } 454 455 /*@ 456 PetscConvEstCreate - Create a PetscConvEst object 457 458 Collective 459 460 Input Parameter: 461 . comm - The communicator for the PetscConvEst object 462 463 Output Parameter: 464 . ce - The PetscConvEst object 465 466 Level: beginner 467 468 .seealso: PetscConvEstDestroy(), PetscConvEstGetConvRate() 469 @*/ 470 PetscErrorCode PetscConvEstCreate(MPI_Comm comm, PetscConvEst *ce) 471 { 472 PetscErrorCode ierr; 473 474 PetscFunctionBegin; 475 PetscValidPointer(ce, 2); 476 ierr = PetscSysInitializePackage();CHKERRQ(ierr); 477 ierr = PetscHeaderCreate(*ce, PETSC_OBJECT_CLASSID, "PetscConvEst", "ConvergenceEstimator", "SNES", comm, PetscConvEstDestroy, PetscConvEstView);CHKERRQ(ierr); 478 (*ce)->monitor = PETSC_FALSE; 479 (*ce)->r = 2.0; 480 (*ce)->Nr = 4; 481 (*ce)->event = -1; 482 (*ce)->ops->setsolver = PetscConvEstSetSNES_Private; 483 (*ce)->ops->initguess = PetscConvEstInitGuessSNES_Private; 484 (*ce)->ops->computeerror = PetscConvEstComputeErrorSNES_Private; 485 (*ce)->ops->getconvrate = PetscConvEstGetConvRateSNES_Private; 486 PetscFunctionReturn(0); 487 } 488