Calculates the Jacobian of the local residual for finite volume spatial discretizations using a finite difference method.
More...
#include <fvbasefdlocallinearizer.hh>
|
| FvBaseFdLocalLinearizer () |
|
| ~FvBaseFdLocalLinearizer () |
|
void | init (Simulator &simulator) |
| Initialize the local Jacobian object. More...
|
|
void | linearize (const Element &element) |
| Compute an element's local Jacobian matrix and evaluate its residual. More...
|
|
void | linearize (ElementContext &elemCtx) |
| Compute an element's local Jacobian matrix and evaluate its residual. More...
|
|
Scalar | numericEpsilon (const ElementContext &elemCtx, int dofIdx, int pvIdx) const |
| Returns the epsilon value which is added and removed from the current solution. More...
|
|
LocalResidual & | localResidual () |
| Return reference to the local residual. More...
|
|
const LocalResidual & | localResidual () const |
| Return reference to the local residual. More...
|
|
const MatrixBlock & | jacobian (int domainScvIdx, int rangeScvIdx) const |
| Returns the local Jacobian matrix of the residual of a sub-control volume. More...
|
|
const MatrixBlock & | jacobianStorage (int dofIdx) const |
| Returns the local Jacobian matrix the storage term of a sub-control volume. More...
|
|
const VectorBlock & | residual (int dofIdx) const |
| Returns the local residual of a sub-control volume. More...
|
|
const VectorBlock & | residualStorage (int dofIdx) const |
| Returns the local storage term of a sub-control volume. More...
|
|
|
Implementation & | asImp_ () |
|
const Implementation & | asImp_ () const |
|
const Simulator & | simulator_ () const |
|
const Problem & | problem_ () const |
|
const Model & | model_ () const |
|
void | resize_ (const ElementContext &elemCtx) |
| Resize all internal attributes to the size of the element. More...
|
|
void | reset_ (const ElementContext &elemCtx) |
| Reset the all relevant internal attributes to 0. More...
|
|
void | evalPartialDerivative_ (ElementContext &elemCtx, int dofIdx, int pvIdx) |
| Compute the partial derivatives of a context's residual functions. More...
|
|
void | updateLocalLinearizer_ (const ElementContext &elemCtx, int primaryDofIdx, int pvIdx) |
| Updates the current local Jacobian matrix with the partial derivatives of all equations in regard to the primary variable 'pvIdx' at vertex 'dofIdx' . More...
|
|
template<class TypeTag>
class Ewoms::FvBaseFdLocalLinearizer< TypeTag >
Calculates the Jacobian of the local residual for finite volume spatial discretizations using a finite difference method.
The local Jacobian for a given context is defined as the derivatives of the residuals of all degrees of freedom featured by the stencil with regard to the primary variables of the stencil's "primary" degrees of freedom.
This class implements numeric differentiation using finite difference methods, i.e. forward or backward differences (2nd order), or central differences (3rd order). The method used is determined by the "NumericDifferenceMethod" property:
- If the value of this property is smaller than 0, backward differences are used, i.e.:
- If the value of this property is 0, central differences are used, i.e.:
- if the value of this property is larger than 0, forward differences are used, i.e.:
Here, is the residual function for all equations, is the value of a sub-control volume's primary variable at the evaluation point and is a small scalar value larger than 0.
Compute the partial derivatives of a context's residual functions.
This method can be overwritten by the implementation if a better scheme than numerical differentiation is available.
The default implementation of this method uses numeric differentiation, i.e. forward or backward differences (2nd order), or central differences (3rd order). The method used is determined by the "NumericDifferenceMethod" property:
- If the value of this property is smaller than 0, backward differences are used, i.e.:
- If the value of this property is 0, central differences are used, i.e.:
- if the value of this property is larger than 0, forward differences are used, i.e.:
Here, is the residual function for all equations, is the value of a sub-control volume's primary variable at the evaluation point and is a small value larger than 0.
- Parameters
-
elemCtx | The element context for which the local partial derivative ought to be calculated |
dofIdx | The sub-control volume index of the current finite element for which the partial derivative ought to be calculated |
pvIdx | The index of the primary variable at the dofIdx' sub-control volume of the current finite element for which the partial derivative ought to be calculated |
References Ewoms::FvBaseFdLocalLinearizer< TypeTag >::asImp_(), Ewoms::FvBaseFdLocalLinearizer< TypeTag >::derivResidual_, Ewoms::FvBaseFdLocalLinearizer< TypeTag >::derivStorage_, Ewoms::FvBaseFdLocalLinearizer< TypeTag >::localResidual_, Ewoms::FvBaseFdLocalLinearizer< TypeTag >::numericDifferenceMethod_(), Ewoms::FvBaseFdLocalLinearizer< TypeTag >::residual_, and Ewoms::FvBaseFdLocalLinearizer< TypeTag >::residualStorage_.
Returns the local Jacobian matrix of the residual of a sub-control volume.
- Parameters
-
domainScvIdx | The local index of the sub control volume which contains the independents |
rangeScvIdx | The local index of the sub control volume which contains the local residual |
References Ewoms::FvBaseFdLocalLinearizer< TypeTag >::jacobian_.
Compute an element's local Jacobian matrix and evaluate its residual.
The local Jacobian for a given context is defined as the derivatives of the residuals of all degrees of freedom featured by the stencil with regard to the primary variables of the stencil's "primary" degrees of freedom. Adding the local Jacobians for all elements in the grid will give the global Jacobian 'grad f(x)'.
- Parameters
-
element | The grid element for which the local residual and its local Jacobian should be calculated. |
References Ewoms::FvBaseFdLocalLinearizer< TypeTag >::internalElemContext_.
Returns the epsilon value which is added and removed from the current solution.
- Parameters
-
elemCtx | The element execution context for which the local residual and its gradient should be calculated. |
dofIdx | The local index of the element's vertex for which the local derivative ought to be calculated. |
pvIdx | The index of the primary variable which gets varied |
References Ewoms::FvBaseFdLocalLinearizer< TypeTag >::baseEpsilon().
The documentation for this class was generated from the following file:
|