petsc-3.14.6 2021-03-30
Report Typos and Errors

PetscConvEstGetConvRate

Returns an estimate of the convergence rate for the discretization

Synopsis

#include "petscconvest.h" 
PetscErrorCode PetscConvEstGetConvRate(PetscConvEst ce, PetscReal alpha[])
Not collective

Input Parameter

ce - The PetscConvEst object

Output Parameter

alpha - The convergence rate for each field

Note: The convergence rate alpha is defined by

|| u_\Delta - u_exact || < C \Delta^alpha
where u_\Delta is the discrete solution, and Delta is a measure of the discretization size. We usually use h for the spatial resolution and \Delta t for the temporal resolution.

We solve a series of problems using increasing resolution (refined meshes or decreased timesteps), calculate an error based upon the exact solution in the DS, and then fit the result to our model above using linear regression.

Options database keys

-snes_convergence_estimate : Execute convergence estimation inside SNESSolve() and print out the rate
-ts_convergence_estimate : Execute convergence estimation inside TSSolve() and print out the rate

See Also

PetscConvEstSetSolver(), PetscConvEstCreate(), PetscConvEstGetConvRate(), SNESSolve(), TSSolve()

Level

intermediate

Location

src/snes/utils/convest.c
Index of all SNES routines
Table of Contents for all manual pages
Index of all manual pages