From: Guido Kanschat Date: Sun, 2 Feb 2014 22:41:39 +0000 (+0000) Subject: add weighted mass matrix X-Git-Tag: v8.2.0-rc1~902 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5d150873656a423aa84526a7f99e89a824cbe3ef;p=dealii.git add weighted mass matrix git-svn-id: https://svn.dealii.org/trunk@32381 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/integrators/l2.h b/deal.II/include/deal.II/integrators/l2.h index 2313997c04..c3004be3e0 100644 --- a/deal.II/include/deal.II/integrators/l2.h +++ b/deal.II/include/deal.II/integrators/l2.h @@ -75,6 +75,48 @@ namespace LocalIntegrators } } + /** + * The weighted mass matrix for scalar or vector values finite elements. + * \f[ + * \int_Z \omega(x) uv\,dx \quad \text{or} \quad \int_Z \omega(x) \mathbf u\cdot \mathbf v\,dx + * \f] + * + * Likewise, this term can be used on faces, where it computes the integrals + * \f[ + * \int_F \omega(x) uv\,ds \quad \text{or} \quad \int_F \omega(x) \mathbf u\cdot \mathbf v\,ds + * \f] + * + * The size of the vector weights must be equal to the + * number of quadrature points in the finite element. + * + * @author Guido Kanschat + * @date 2014 + */ + template + void weighted_mass_matrix ( + FullMatrix &M, + const FEValuesBase &fe, + const std::vector weights) + { + const unsigned int n_dofs = fe.dofs_per_cell; + const unsigned int n_components = fe.get_fe().n_components(); + AssertDimension(M.m(), n_dofs); + AssertDimension(M.n(), n_dofs); + AssertDimension(weights.size(), fe.n_quadrature_points); + + for (unsigned int k=0; kL2-inner product for scalar functions. *