Reference documentation for deal.II version 9.1.1
|
Local integrators related to L2-inner products. More...
Functions | |
template<int dim> | |
void | mass_matrix (FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.) |
template<int dim> | |
void | weighted_mass_matrix (FullMatrix< double > &M, const FEValuesBase< dim > &fe, const std::vector< double > &weights) |
template<int dim, typename number > | |
void | L2 (Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.) |
template<int dim, typename number > | |
void | L2 (Vector< number > &result, const FEValuesBase< dim > &fe, const ArrayView< const std::vector< double >> &input, const double factor=1.) |
template<int dim> | |
void | jump_matrix (FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const double factor1=1., const double factor2=1.) |
Local integrators related to L2-inner products.
void LocalIntegrators::L2::mass_matrix | ( | FullMatrix< double > & | M, |
const FEValuesBase< dim > & | fe, | ||
const double | factor = 1. |
||
) |
The mass matrix for scalar or vector values finite elements.
\[ \int_Z uv\,dx \quad \text{or} \quad \int_Z \mathbf u\cdot \mathbf v\,dx \]
Likewise, this term can be used on faces, where it computes the integrals
\[ \int_F uv\,ds \quad \text{or} \quad \int_F \mathbf u\cdot \mathbf v\,ds \]
M | The mass matrix obtained as result. |
fe | The FEValues object describing the local trial function space. update_values and update_JxW_values must be set. |
factor | A constant that multiplies the mass matrix. |
void LocalIntegrators::L2::weighted_mass_matrix | ( | FullMatrix< double > & | M, |
const FEValuesBase< dim > & | fe, | ||
const std::vector< double > & | weights | ||
) |
The weighted mass matrix for scalar or vector values finite elements.
\[ \int_Z \omega(x) uv\,dx \quad \text{or} \quad \int_Z \omega(x) \mathbf u\cdot \mathbf v\,dx \]
Likewise, this term can be used on faces, where it computes the integrals
\[ \int_F \omega(x) uv\,ds \quad \text{or} \quad \int_F \omega(x) \mathbf u\cdot \mathbf v\,ds \]
M | The weighted mass matrix obtained as result. |
fe | The FEValues object describing the local trial function space. update_values and update_JxW_values must be set. |
weights | The weights, \(\omega(x)\), evaluated at the quadrature points in the finite element (size must be equal to the number of quadrature points in the element). |
void LocalIntegrators::L2::L2 | ( | Vector< number > & | result, |
const FEValuesBase< dim > & | fe, | ||
const std::vector< double > & | input, | ||
const double | factor = 1. |
||
) |
L2-inner product for scalar functions.
\[ \int_Z fv\,dx \quad \text{or} \quad \int_F fv\,ds \]
result | The vector obtained as result. |
fe | The FEValues object describing the local trial function space. update_values and update_JxW_values must be set. |
input | The representation of \(f\) evaluated at the quadrature points in the finite element (size must be equal to the number of quadrature points in the element). |
factor | A constant that multiplies the result. |
void LocalIntegrators::L2::L2 | ( | Vector< number > & | result, |
const FEValuesBase< dim > & | fe, | ||
const ArrayView< const std::vector< double >> & | input, | ||
const double | factor = 1. |
||
) |
L2-inner product for a slice of a vector valued right hand side.
\[ \int_Z \mathbf f\cdot \mathbf v\,dx \quad \text{or} \quad \int_F \mathbf f\cdot \mathbf v\,ds \]
result | The vector obtained as result. |
fe | The FEValues object describing the local trial function space. update_values and update_JxW_values must be set. |
input | The vector valued representation of \(\mathbf f\) evaluated at the quadrature points in the finite element (size of each component must be equal to the number of quadrature points in the element). |
factor | A constant that multiplies the result. |
void LocalIntegrators::L2::jump_matrix | ( | FullMatrix< double > & | M11, |
FullMatrix< double > & | M12, | ||
FullMatrix< double > & | M21, | ||
FullMatrix< double > & | M22, | ||
const FEValuesBase< dim > & | fe1, | ||
const FEValuesBase< dim > & | fe2, | ||
const double | factor1 = 1. , |
||
const double | factor2 = 1. |
||
) |
The jump matrix between two cells for scalar or vector values finite elements. Note that the factor \(\gamma\) can be used to implement weighted jumps.
\[ \int_F [\gamma u][\gamma v]\,ds \quad \text{or} \int_F [\gamma \mathbf u]\cdot [\gamma \mathbf v]\,ds \]
Using appropriate weights, this term can be used to penalize violation of conformity in H1.
Note that for the parameters that follow, the external matrix refers to the flux between cells, while the internal matrix refers to entries coupling inside the cell.
M11 | The internal matrix for the first cell obtained as result. |
M12 | The external matrix for the first cell obtained as result. |
M12 | The external matrix for the second cell obtained as result. |
M22 | The internal matrix for the second cell obtained as result. |
fe1 | The FEValues object describing the local trial function space for the first cell. update_values and update_JxW_values must be set. |
fe2 | The FEValues object describing the local trial function space for the second cell. update_values and update_JxW_values must be set. |
factor1 | A constant that multiplies the shape functions for the first cell. |
factor2 | A constant that multiplies the shape functions for the second cell. |