From 3b415bbc6709671e308123d5847815d5d7139f93 Mon Sep 17 00:00:00 2001 From: kanschat Date: Sun, 28 Nov 2010 23:19:04 +0000 Subject: [PATCH] begin operators for LDG git-svn-id: https://svn.dealii.org/trunk@22869 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/integrators/differential.h | 4 +- deal.II/include/deal.II/integrators/l2.h | 66 +++++++++++++++++-- 2 files changed, 64 insertions(+), 6 deletions(-) diff --git a/deal.II/include/deal.II/integrators/differential.h b/deal.II/include/deal.II/integrators/differential.h index 00af25a290..e9eac34ed8 100644 --- a/deal.II/include/deal.II/integrators/differential.h +++ b/deal.II/include/deal.II/integrators/differential.h @@ -39,7 +39,7 @@ namespace LocalIntegrators * Cell matrix for divergence. The derivative is on the trial function. * * \f[ - * \int_Z v\nabla \cdot u \,dx + * \int_Z v\nabla \cdot \mathbf u \,dx * \f] */ template @@ -76,7 +76,7 @@ namespace LocalIntegrators * Cell matrix for divergence. The derivative is on the test function. * * \f[ - * \int_Z \nabla u \cdot v\,dx + * \int_Z \nabla u \cdot \mathbf v\,dx * \f] */ template diff --git a/deal.II/include/deal.II/integrators/l2.h b/deal.II/include/deal.II/integrators/l2.h index 6b3bc2fa3f..4deead2654 100644 --- a/deal.II/include/deal.II/integrators/l2.h +++ b/deal.II/include/deal.II/integrators/l2.h @@ -36,10 +36,10 @@ namespace LocalIntegrators namespace L2 { /** - * The mass matrix. + * The mass matrix for scalar or vector values finite elements. * * \f[ - * (a u,v) + * \int_Z uv\,dx \quad \text{or} \int_Z \mathbf u\cdot \mathbf v\,dx * \f] */ template @@ -68,7 +68,7 @@ namespace LocalIntegrators * L2-inner product for scalar functions. * * \f[ - * (f,v) + * \int_Z fv\,dx * \f] */ template @@ -93,7 +93,7 @@ namespace LocalIntegrators * right hand side. * * \f[ - * \int_Z f\cdot v\,dx + * \int_Z \mathbf f\cdot \mathbf v\,dx * \f] */ template @@ -116,6 +116,64 @@ namespace LocalIntegrators result(i) += fe.JxW(k) * factor * fe.shape_value_component(i,k,d) * input[d][k]; } } + +/** + * 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. + * + * \f[ + * \int_F [\gamma u][\gamma v]\,dx + * \quad \text{or} + * \int_F [\gamma \mathbf u]\cdot [\gamma \mathbf v]\,dx + * \f] + */ + template + void jump_matrix ( + FullMatrix& M11, + FullMatrix& M12, + FullMatrix& M21, + FullMatrix& M22, + const FEValuesBase& fe1, + const FEValuesBase& fe2, + const double factor1 = 1., + const double factor2 = 1.) + { + const unsigned int n1_dofs = fe1.dofs_per_cell; + const unsigned int n2_dofs = fe2.dofs_per_cell; + const unsigned int n_components = fe1.get_fe().n_components(); + + Assert(n1_dofs == n2_dofs, ExcNotImplemented()); + AssertDimension(n_components, fe2.get_fe().n_components()); + AssertDimension(M11.m(), n1_dofs); + AssertDimension(M12.m(), n1_dofs); + AssertDimension(M21.m(), n2_dofs); + AssertDimension(M22.m(), n2_dofs); + AssertDimension(M11.n(), n1_dofs); + AssertDimension(M12.n(), n2_dofs); + AssertDimension(M21.n(), n1_dofs); + AssertDimension(M22.n(), n2_dofs); + + for (unsigned k=0;k