From: kanschat Date: Tue, 24 Sep 2013 09:04:36 +0000 (+0000) Subject: add homogeneous boundary operator X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=16608302cec61df8ea360ff40646b57da363d0d1;p=dealii-svn.git add homogeneous boundary operator git-svn-id: https://svn.dealii.org/trunk@30907 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/integrators/elasticity.h b/deal.II/include/deal.II/integrators/elasticity.h index d772cf695c..6623ba9390 100644 --- a/deal.II/include/deal.II/integrators/elasticity.h +++ b/deal.II/include/deal.II/integrators/elasticity.h @@ -214,6 +214,60 @@ namespace LocalIntegrators } } + /** + * Homogeneous weak boundary condition for the elasticity operator by Nitsche, + * namely on the face F the vector + * @f[ + * \int_F \Bigl(\gamma u \cdot v - n^T \epsilon(u) v - u \epsilon(v) n^T\Bigr)\;ds. + * @f] + * + * Here, u is the finite element function whose values and + * gradient are given in the arguments input and + * Dinput, respectively. $n$ is the outer + * normal vector and $\gamma$ is the usual penalty parameter. + * + * @author Guido Kanschat + * @date 2013 + */ + template + void nitsche_residual_homogeneous ( + Vector &result, + const FEValuesBase &fe, + const VectorSlice > > &input, + const VectorSlice > > > &Dinput, + double penalty, + double factor = 1.) + { + const unsigned int n_dofs = fe.dofs_per_cell; + AssertVectorVectorDimension(input, dim, fe.n_quadrature_points); + AssertVectorVectorDimension(Dinput, dim, fe.n_quadrature_points); + + for (unsigned int k=0; k &n = fe.normal_vector(k); + for (unsigned int i=0; i