From b71b0f7807a187e9d66f0a7a75ed1dfadedb0cf9 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Wed, 8 Aug 2012 15:47:32 +0000 Subject: [PATCH] add cell residual git-svn-id: https://svn.dealii.org/trunk@25776 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/integrators/laplace.h | 70 +++++++++++++++++++ 1 file changed, 70 insertions(+) diff --git a/deal.II/include/deal.II/integrators/laplace.h b/deal.II/include/deal.II/integrators/laplace.h index bebeb55341..08e5f1c929 100644 --- a/deal.II/include/deal.II/integrators/laplace.h +++ b/deal.II/include/deal.II/integrators/laplace.h @@ -69,6 +69,71 @@ namespace LocalIntegrators } } } + +/** + * Laplacian residual operator in weak form + * + * \f[ + * -(\nabla u, \nabla v) + * \f] + */ + template + inline void + cell_residual ( + Vector& result, + const FEValuesBase& fe, + const std::vector >& input, + double factor = 1.) + { + const unsigned int nq = fe.n_quadrature_points; + const unsigned int n_dofs = fe.dofs_per_cell; + Assert(input.size() == nq, ExcDimensionMismatch(input.size(), nq)); + Assert(result.size() == n_dofs, ExcDimensionMismatch(result.size(), n_dofs)); + + for (unsigned k=0;k + inline void + cell_residual ( + Vector& result, + const FEValuesBase& fe, + const VectorSlice > > >& input, + double factor = 1.) + { + const unsigned int nq = fe.n_quadrature_points; + const unsigned int n_dofs = fe.dofs_per_cell; + const unsigned int n_comp = fe.get_fe().n_components(); + + AssertVectorVectorDimension(input, n_comp, fe.n_quadrature_points); + Assert(result.size() == n_dofs, ExcDimensionMismatch(result.size(), n_dofs)); + + for (unsigned k=0;kF the matrix * @f[ @@ -174,6 +239,10 @@ namespace LocalIntegrators * Dinput, respectively. g is the inhomogeneous * boundary value in the argument data. $\gamma$ is the usual * penalty parameter. + * + * @ingroup Integrators + * @author Guido Kanschat + * @date 2008, 2009, 2010 */ template void nitsche_residual ( @@ -291,6 +360,7 @@ namespace LocalIntegrators * hi is the length of Zi * orthogonal to the current face. * + * @ingroup Integrators * @author Guido Kanschat * @date 2010 */ -- 2.39.5