Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
oomph::Z2ErrorEstimator Class Reference

#include <error_estimator.h>

+ Inheritance diagram for oomph::Z2ErrorEstimator:

Public Types

typedef double(* CombinedErrorEstimateFctPt )(const Vector< double > &errors)
 Function pointer to combined error estimator function. More...
 

Public Member Functions

 Z2ErrorEstimator (const unsigned &recovery_order)
 Constructor: Set order of recovery shape functions. More...
 
 Z2ErrorEstimator ()
 Constructor: Leave order of recovery shape functions open for now – they will be read out from first element of the mesh when the error estimator is applied. More...
 
 Z2ErrorEstimator (const Z2ErrorEstimator &)
 Broken copy constructor. More...
 
void operator= (const Z2ErrorEstimator &)
 Broken assignment operator. More...
 
virtual ~Z2ErrorEstimator ()
 Empty virtual destructor. More...
 
void get_element_errors (Mesh *&mesh_pt, Vector< double > &elemental_error)
 Compute the elemental error measures for a given mesh and store them in a vector. More...
 
void get_element_errors (Mesh *&mesh_pt, Vector< double > &elemental_error, DocInfo &doc_info)
 Compute the elemental error measures for a given mesh and store them in a vector. If doc_info.enable_doc(), doc FE and recovered fluxes in. More...
 
unsigned & recovery_order ()
 Access function for order of recovery polynomials. More...
 
unsigned recovery_order () const
 Access function for order of recovery polynomials (const version) More...
 
CombinedErrorEstimateFctPtcombined_error_fct_pt ()
 Access function: Pointer to combined error estimate function. More...
 
CombinedErrorEstimateFctPt combined_error_fct_pt () const
 Access function: Pointer to combined error estimate function. Const version. More...
 
void setup_patches (Mesh *&mesh_pt, std::map< Node *, Vector< ElementWithZ2ErrorEstimator * > * > &adjacent_elements_pt, Vector< Node * > &vertex_node_pt)
 Setup patches: For each vertex node pointed to by nod_pt, adjacent_elements_pt[nod_pt] contains the pointer to the vector that contains the pointers to the elements that the node is part of. More...
 
double & reference_flux_norm ()
 Access function for prescribed reference flux norm. More...
 
double reference_flux_norm () const
 Access function for prescribed reference flux norm (const. version) More...
 
double get_combined_error_estimate (const Vector< double > &compound_error)
 Return a combined error estimate from all compound errors. More...
 
- Public Member Functions inherited from oomph::ErrorEstimator
 ErrorEstimator ()
 Default empty constructor. More...
 
 ErrorEstimator (const ErrorEstimator &)
 Broken copy constructor. More...
 
void operator= (const ErrorEstimator &)
 Broken assignment operator. More...
 
virtual ~ErrorEstimator ()
 Empty virtual destructor. More...
 
void get_element_errors (Mesh *&mesh_pt, Vector< double > &elemental_error)
 Compute the elemental error-measures for a given mesh and store them in a vector. More...
 

Private Member Functions

void get_recovered_flux_in_patch (const Vector< ElementWithZ2ErrorEstimator * > &patch_el_pt, const unsigned &num_recovery_terms, const unsigned &num_flux_terms, const unsigned &dim, DenseMatrix< double > *&recovered_flux_coefficient_pt)
 Given the vector of elements that make up a patch, the number of recovery and flux terms, and the spatial dimension of the problem, compute the matrix of recovered flux coefficients and return a pointer to it. More...
 
unsigned nrecovery_terms (const unsigned &dim)
 Return number of coefficients for expansion of recovered fluxes for given spatial dimension of elements. (We use complete polynomials of the specified given order.) More...
 
void shape_rec (const Vector< double > &x, const unsigned &dim, Vector< double > &psi_r)
 Recovery shape functions as functions of the global, Eulerian coordinate x of dimension dim. The recovery shape functions are complete polynomials of the order specified by Recovery_order. More...
 
Integralintegral_rec (const unsigned &dim, const bool &is_q_mesh)
 Integation scheme associated with the recovery shape functions must be of sufficiently high order to integrate the mass matrix associated with the recovery shape functions. More...
 
void doc_flux (Mesh *mesh_pt, const unsigned &num_flux_terms, MapMatrixMixed< Node *, int, double > &rec_flux_map, const Vector< double > &elemental_error, DocInfo &doc_info)
 Doc flux and recovered flux. More...
 

Private Attributes

unsigned Recovery_order
 Order of recovery polynomials. More...
 
bool Recovery_order_from_first_element
 
double Reference_flux_norm
 Prescribed reference flux norm. More...
 
CombinedErrorEstimateFctPt Combined_error_fct_pt
 Function pointer to combined error estimator function. More...
 

Detailed Description

Z2-error-estimator: Elements that can be used with Z2 error estimation should be derived from the base class ElementWithZ2ErrorEstimator and implement its pure virtual member functions to provide the following functionality:

As an example consider the finite element solution of the Laplace problem, $ \partial^2 u/\partial x_i^2 = 0 $. If we approximate the unknown $ u $ on a finite element mesh with $ N $ nodes as

\[ u^{[FE]}(x_k) = \sum_{j=1}^{N} U_j \ \psi_j(x_k), \]

where the $ \psi_j(x_k) $ are the (global) $ C^0 $ basis functions, the finite-element representation of the flux, $ f_i = \partial u/\partial x_i $,

\[ f_i^{[FE]} = \sum_{j=1}^{N} U_j \ \frac{\partial \psi_j}{\partial x_i} \]

is discontinuous between elements but the magnitude of the jump decreases under mesh refinement. We denote the number of flux terms by $N_{flux}$, so for a 2D (3D) Laplace problem, $N_{flux}=2 \ (3).$

The idea behind Z2 error estimation is to compute an approximation for the flux that is more accurate than the value $ f_i^{[FE]} $ obtained directly from the finite element solution. We refer to the improved approximation for the flux as the "recovered flux" and denote it by $ f_i^{[rec]} $. Since $ f_i^{[rec]} $ is more accurate than $ f_i^{[FE]} $, the difference between the two provides a measure of the error. In Z2 error estimation, the "recovered flux" is determined by projecting the discontinuous, FE-based flux $ f_i^{[FE]} $ onto a set of continuous basis functions, known as the "recovery shape functions". In principle, one could use the finite element shape functions $ \psi_j(x_k) $ as the recovery shape functions but then the determination of the recovered flux would involve the solution of a linear system that is as big as the original problem. To avoid this, the projection is performed over small patches of elements within which low-order polynomials (defined in terms of the global, Eulerian coordinates) are used as the recovery shape functions.

Z2 error estimation is known to be "optimal" for many self-adjoint problems. See, e.g., Zienkiewicz & Zhu's original papers "The superconvergent patch recovery and a posteriori error estimates." International Journal for Numerical Methods in Engineering 33 (1992) Part I: 1331-1364; Part II: 1365-1382. In non-self adjoint problems, the error estimator only analyses the "self-adjoint part" of the differential operator. However, in many cases, this still produces good error indicators since the error estimator detects under-resolved, sharp gradients in the solution.

Z2 error estimation works as follows:

  1. We combine the elements in the mesh into a number of "patches", which consist of all elements that share a common vertex node. Most elements will therefore be members of multiple patches.
  2. Within each patch $p$, we expand the "recovered flux" as

    \[ f^{[rec](p)}_i(x_k) = \sum_{j=1}^{N_{rec}} F^{(p)}_{ij} \ \psi^{[rec]}_j(x_k) \mbox{ \ \ \ for $i=1,...,N_{flux}$,} \]

    where the functions $ \psi^{[rec]}_j(x_k)$ are the recovery shape functions, which are functions of the global, Eulerian coordinates. Typically, these are chosen to be low-order polynomials. For instance, in 2D, a bilinear representation of $ f^{(p)}_i(x_0,x_1) $ involves the $N_{rec}=3$ recovery shape functions $ \psi^{[rec]}_0(x_0,x_1)=1, \ \psi^{[rec]}_1(x_0,x_1)=x_0 $ and $ \psi^{[rec]}_2(x_0,x_1)=x_1$.

    We determine the coefficients $ F^{(p)}_{ij} $ by enforcing $ f^{(p)}_i(x_k) = f^{[FE]}_i(x_k)$ in its weak form:

    \[ \int_{\mbox{Patch $p$}} \left( f^{[FE]}_i(x_k) - \sum_{j=1}^{N_{rec}} F^{(p)}_{ij} \ \psi^{[rec]}_j(x_k) \right) \psi^{[rec]}_l(x_k)\ dv = 0 \mbox{ \ \ \ \ for $l=1,...,N_{rec}$ and $i=1,...,N_{flux}$}. \]

    Once the $ F^{(p)}_{ij} $ are determined in a given patch, we determine the values of the recovered flux at all nodes that are part of the patch. We denote the value of the recovered flux at node $ k $ by $ {\cal F}^{(p)}_{ik}$.

    We repeat this procedure for every patch. For nodes that are part of multiple patches, the procedure will provide multiple, slightly different nodal values for the recovered flux. We average these values via

    \[ {\cal F}_{ij} = \frac{1}{N_p(j)} \sum_{\mbox{Node $j \in $ patch $p$}} {\cal F}^{(p)}_{ij}, \]

    where $N_p(j)$ denotes the number of patches that node $ j$ is a member of. This allows us to obtain a globally-continuous, finite-element based representation of the recovered flux as

    \[ f_i^{[rec]} = \sum_{j=1}^{N} {\cal F}_{ij}\ \psi_j, \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ (1) \]

    where the $ \psi_j $ are the (global) finite element shape functions.

  3. Since the recovered flux in (1), is based on nodal values, we can evaluate it locally within each of the $ N_e$ elements in the mesh to obtain a normalised elemental error estimate via

    \[ E_{e} = \sqrt{ \frac{ \int_{\mbox{$e$}} \left( f_i^{[rec]} - f_i^{[FE]} \right)^2 dv} {\sum_{e'=1}^{N_e} \int_{\mbox{$e'$}} \left( f_i^{[rec]} \right)^2 dv} } \mbox{\ \ \ for $e=1,...,N_e$.} \]

    In this (default) form, mesh refinement, based on this error estimate will lead to an equidistribution of the error across all elements. Usually, this is the desired behaviour. However, there are cases in which the solution evolves towards a state in which the flux tends to zero while the solution itself becomes so simple that it can be represented exactly on any finite element mesh (e.g. in spin-up problems in which the fluid motion approaches a rigid body rotation). In that case the mesh adaptation tries to equidistribute any small roundoff errors in the solution, causing strong, spatially uniform mesh refinement throughout the domain, even though the solution is already fully converged. For such cases, it is possible to replace the denominator in the above expression by a prescribed reference flux, which may be specified via the access function

Definition at line 313 of file error_estimator.h.

Member Typedef Documentation

typedef double(* oomph::Z2ErrorEstimator::CombinedErrorEstimateFctPt)(const Vector< double > &errors)

Function pointer to combined error estimator function.

Definition at line 319 of file error_estimator.h.

Constructor & Destructor Documentation

oomph::Z2ErrorEstimator::Z2ErrorEstimator ( const unsigned &  recovery_order)
inline

Constructor: Set order of recovery shape functions.

Definition at line 322 of file error_estimator.h.

oomph::Z2ErrorEstimator::Z2ErrorEstimator ( )
inline

Constructor: Leave order of recovery shape functions open for now – they will be read out from first element of the mesh when the error estimator is applied.

Definition at line 331 of file error_estimator.h.

oomph::Z2ErrorEstimator::Z2ErrorEstimator ( const Z2ErrorEstimator )
inline

Broken copy constructor.

Definition at line 337 of file error_estimator.h.

References oomph::BrokenCopy::broken_copy().

virtual oomph::Z2ErrorEstimator::~Z2ErrorEstimator ( )
inlinevirtual

Empty virtual destructor.

Definition at line 349 of file error_estimator.h.

Member Function Documentation

CombinedErrorEstimateFctPt& oomph::Z2ErrorEstimator::combined_error_fct_pt ( )
inline

Access function: Pointer to combined error estimate function.

Definition at line 379 of file error_estimator.h.

References Combined_error_fct_pt.

CombinedErrorEstimateFctPt oomph::Z2ErrorEstimator::combined_error_fct_pt ( ) const
inline

Access function: Pointer to combined error estimate function. Const version.

Definition at line 384 of file error_estimator.h.

References Combined_error_fct_pt.

void oomph::Z2ErrorEstimator::doc_flux ( Mesh mesh_pt,
const unsigned &  num_flux_terms,
MapMatrixMixed< Node *, int, double > &  rec_flux_map,
const Vector< double > &  elemental_error,
DocInfo doc_info 
)
private
double oomph::Z2ErrorEstimator::get_combined_error_estimate ( const Vector< double > &  compound_error)

Return a combined error estimate from all compound errors.

Return a combined error estimate from all compound flux errors The default is to return the maximum of the compound flux errors which will always force refinment if any field is above the single mesh error threshold and unrefinement if both are below the lower limit. Any other fancy combinations can be selected by specifying a user-defined combined estimate by setting a function pointer.

Definition at line 422 of file error_estimator.cc.

References Combined_error_fct_pt, and i.

Referenced by get_element_errors().

void oomph::Z2ErrorEstimator::get_element_errors ( Mesh *&  mesh_pt,
Vector< double > &  elemental_error 
)
inline

Compute the elemental error measures for a given mesh and store them in a vector.

Definition at line 353 of file error_estimator.h.

References oomph::DocInfo::disable_doc().

void oomph::Z2ErrorEstimator::get_element_errors ( Mesh *&  mesh_pt,
Vector< double > &  elemental_error,
DocInfo doc_info 
)
virtual

Compute the elemental error measures for a given mesh and store them in a vector. If doc_info.enable_doc(), doc FE and recovered fluxes in.

  • flux_fe*.dat
  • flux_rec*.dat

Get Vector of Z2-based error estimates for all elements in mesh. If doc_info.is_doc_enabled()=true, doc FE and recovered fluxes in

  • flux_fe*.dat
  • flux_rec*.dat

Implements oomph::ErrorEstimator.

Definition at line 901 of file error_estimator.cc.

References oomph::MPI_Helpers::communicator_pt(), oomph::Mesh::communicator_pt(), oomph::FiniteElement::dim(), doc_flux(), e, oomph::Mesh::element_pt(), oomph::ElementWithZ2ErrorEstimator::geometric_jacobian(), get_combined_error_estimate(), get_recovered_flux_in_patch(), oomph::ElementWithZ2ErrorEstimator::get_Z2_compound_flux_indices(), oomph::ElementWithZ2ErrorEstimator::get_Z2_flux(), oomph::Mesh::haloed_element_pt(), i, oomph::FiniteElement::integral_pt(), oomph::FiniteElement::interpolated_x(), oomph::DocInfo::is_doc_enabled(), oomph::GeneralisedElement::is_halo(), oomph::Mesh::is_mesh_distributed(), oomph::FiniteElement::J_eulerian(), oomph::Integral::knot(), oomph::MPI_Helpers::mpi_has_been_initialised(), oomph::ElementWithZ2ErrorEstimator::ncompound_fluxes(), oomph::Mesh::nelement(), oomph::Mesh::nnode(), oomph::FiniteElement::nnode(), oomph::Mesh::node_pt(), oomph::FiniteElement::node_pt(), oomph::ElementWithZ2ErrorEstimator::nrecovery_order(), nrecovery_terms(), oomph::ElementWithZ2ErrorEstimator::num_Z2_flux_terms(), oomph::Integral::nweight(), recovery_order(), Recovery_order, Recovery_order_from_first_element, Reference_flux_norm, s, setup_patches(), oomph::FiniteElement::shape(), shape_rec(), oomph::QuadTreeNames::W, oomph::Integral::weight(), and oomph::Node::x().

void oomph::Z2ErrorEstimator::get_recovered_flux_in_patch ( const Vector< ElementWithZ2ErrorEstimator * > &  patch_el_pt,
const unsigned &  num_recovery_terms,
const unsigned &  num_flux_terms,
const unsigned &  dim,
DenseMatrix< double > *&  recovered_flux_coefficient_pt 
)
private

Given the vector of elements that make up a patch, the number of recovery and flux terms, and the spatial dimension of the problem, compute the matrix of recovered flux coefficients and return a pointer to it.

Given the vector of elements that make up a patch, the number of recovery and flux terms, and the spatial dimension of the problem, compute the matrix of recovered flux coefficients and return a pointer to it.

Definition at line 618 of file error_estimator.cc.

References e, oomph::ElementWithZ2ErrorEstimator::geometric_jacobian(), oomph::ElementWithZ2ErrorEstimator::get_Z2_flux(), i, integral_rec(), oomph::FiniteElement::interpolated_x(), oomph::FiniteElement::J_eulerian(), oomph::Integral::knot(), oomph::DenseDoubleMatrix::lubksub(), oomph::DenseDoubleMatrix::ludecompose(), oomph::Integral::nweight(), s, shape_rec(), oomph::QuadTreeNames::W, and oomph::Integral::weight().

Referenced by get_element_errors().

Integral * oomph::Z2ErrorEstimator::integral_rec ( const unsigned &  dim,
const bool &  is_q_mesh 
)
private

Integation scheme associated with the recovery shape functions must be of sufficiently high order to integrate the mass matrix associated with the recovery shape functions.

Integation scheme associated with the recovery shape functions must be of sufficiently high order to integrate the mass matrix associated with the recovery shape functions. The argument is the dimension of the elements. The integration is performed locally over the elements, so the integration scheme does depend on the geometry of the element. The type of element is specified by the boolean which is true if elements in the patch are QElements and false if they are TElements (will need change if we ever have other element types)

Which spatial dimension are we dealing with?

Find order of recovery shape functions

Find order of recovery shape functions

Find order of recovery shape functions

Definition at line 256 of file error_estimator.cc.

References recovery_order().

Referenced by get_recovered_flux_in_patch().

unsigned oomph::Z2ErrorEstimator::nrecovery_terms ( const unsigned &  dim)
private

Return number of coefficients for expansion of recovered fluxes for given spatial dimension of elements. (We use complete polynomials of the specified given order.)

Number of coefficients for expansion of recovered fluxes for given spatial dimension of elements. Use complete polynomial of given order for recovery

Definition at line 762 of file error_estimator.cc.

References Recovery_order, and oomph::Global_string_for_annotation::string().

Referenced by get_element_errors().

void oomph::Z2ErrorEstimator::operator= ( const Z2ErrorEstimator )
inline

Broken assignment operator.

Definition at line 343 of file error_estimator.h.

References oomph::BrokenCopy::broken_assign().

unsigned& oomph::Z2ErrorEstimator::recovery_order ( )
inline

Access function for order of recovery polynomials.

Definition at line 373 of file error_estimator.h.

References Recovery_order.

Referenced by get_element_errors(), integral_rec(), and shape_rec().

unsigned oomph::Z2ErrorEstimator::recovery_order ( ) const
inline

Access function for order of recovery polynomials (const version)

Definition at line 376 of file error_estimator.h.

References Recovery_order.

double& oomph::Z2ErrorEstimator::reference_flux_norm ( )
inline

Access function for prescribed reference flux norm.

Definition at line 397 of file error_estimator.h.

References Reference_flux_norm.

double oomph::Z2ErrorEstimator::reference_flux_norm ( ) const
inline

Access function for prescribed reference flux norm (const. version)

Definition at line 400 of file error_estimator.h.

References Reference_flux_norm.

void oomph::Z2ErrorEstimator::setup_patches ( Mesh *&  mesh_pt,
std::map< Node *, Vector< ElementWithZ2ErrorEstimator * > * > &  adjacent_elements_pt,
Vector< Node * > &  vertex_node_pt 
)

Setup patches: For each vertex node pointed to by nod_pt, adjacent_elements_pt[nod_pt] contains the pointer to the vector that contains the pointers to the elements that the node is part of.

Setup patches: For each vertex node pointed to by nod_pt, adjacent_elements_pt[nod_pt] contains the pointer to the vector that contains the pointers to the elements that the node is part of. Also returns a Vector of vertex nodes for use in get_element_errors.

Definition at line 463 of file error_estimator.cc.

References e, oomph::Mesh::element_pt(), oomph::Mesh::nelement(), oomph::FiniteElement::nnode(), oomph::FiniteElement::node_pt(), oomph::ElementWithZ2ErrorEstimator::nrecovery_order(), oomph::ElementWithZ2ErrorEstimator::nvertex_node(), oomph::oomph_info, and oomph::ElementWithZ2ErrorEstimator::vertex_node_pt().

Referenced by get_element_errors().

void oomph::Z2ErrorEstimator::shape_rec ( const Vector< double > &  x,
const unsigned &  dim,
Vector< double > &  psi_r 
)
private

Recovery shape functions as functions of the global, Eulerian coordinate x of dimension dim. The recovery shape functions are complete polynomials of the order specified by Recovery_order.

Recovery shape functions as functions of the global, Eulerian coordinate x of dimension dim. The recovery shape functions are complete polynomials of the order specified by Recovery_order.

Which spatial dimension are we dealing with?

Find order of recovery shape functions

Find order of recovery shape functions

Find order of recovery shape functions

Definition at line 50 of file error_estimator.cc.

References recovery_order().

Referenced by get_element_errors(), and get_recovered_flux_in_patch().

Member Data Documentation

CombinedErrorEstimateFctPt oomph::Z2ErrorEstimator::Combined_error_fct_pt
private

Function pointer to combined error estimator function.

Definition at line 456 of file error_estimator.h.

Referenced by combined_error_fct_pt(), and get_combined_error_estimate().

unsigned oomph::Z2ErrorEstimator::Recovery_order
private

Order of recovery polynomials.

Definition at line 439 of file error_estimator.h.

Referenced by get_element_errors(), nrecovery_terms(), and recovery_order().

bool oomph::Z2ErrorEstimator::Recovery_order_from_first_element
private

Bool to indicate if recovery order is to be read out from first element in mesh or set globally

Definition at line 443 of file error_estimator.h.

Referenced by get_element_errors().

double oomph::Z2ErrorEstimator::Reference_flux_norm
private

Prescribed reference flux norm.

Definition at line 453 of file error_estimator.h.

Referenced by get_element_errors(), and reference_flux_norm().


The documentation for this class was generated from the following files: