]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add scalar version of Nitsche boundary residual
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 14 Nov 2010 22:29:29 +0000 (22:29 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Sun, 14 Nov 2010 22:29:29 +0000 (22:29 +0000)
git-svn-id: https://svn.dealii.org/trunk@22723 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/integrators/laplace.h

index 9d213be76d4204ccf5565476251806f2ffbcbc32..47ed06e1dca0f887d5d60220d59476fee9c186f7 100644 (file)
@@ -109,7 +109,8 @@ namespace LocalIntegrators
     }
 
 /**
- * Weak boundary condition for the Laplace operator by Nitsche, namely on the face <i>F</i>
+ * Weak boundary condition for the Laplace operator by Nitsche, vector
+ * valued version, namely on the face <i>F</i>
  * the vector
  * @f[
  * \int_F \Bigl(\gamma (u-g) v - \partial_n u v - (u-g) \partial_n v\Bigr)\;ds.
@@ -156,6 +157,51 @@ namespace LocalIntegrators
          }
       }
 
+/**
+ * Weak boundary condition for the Laplace operator by Nitsche, scalar
+ * version, namely on the face <i>F</i> the vector
+ * @f[
+ * \int_F \Bigl(\gamma (u-g) v - \partial_n u v - (u-g) \partial_n v\Bigr)\;ds.
+ * @f]
+ *
+ * Here, <i>u</i> is the finite element function whose values and
+ * gradient are given in the arguments <tt>input</tt> and
+ * <tt>Dinput</tt>, respectively. <i>g</i> is the inhomogeneous
+ * boundary value in the argument <tt>data</tt>. $\gamma$ is the usual
+ * penalty parameter.
+ */
+      template <int dim>
+      void nitsche_residual (
+       Vector<double>& result,
+       const FEValuesBase<dim>& fe,
+       const std::vector<double>& input,
+       const std::vector<Tensor<1,dim> >& Dinput,
+       const std::vector<double>& data,
+       double penalty,
+       double factor = 1.)
+      {
+       const unsigned int n_dofs = fe.dofs_per_cell;
+       AssertDimension(input.size(), fe.n_quadrature_points);
+       AssertDimension(Dinput.size(), fe.n_quadrature_points);
+       AssertDimension(data.size(), fe.n_quadrature_points);
+
+       for (unsigned k=0;k<fe.n_quadrature_points;++k)
+         {
+           const double dx = factor * fe.JxW(k);
+           const Point<dim>& n = fe.normal_vector(k);
+           for (unsigned i=0;i<n_dofs;++i)
+             {
+               const double dnv = fe.shape_grad(i,k) * n;
+               const double dnu = Dinput[k] * n;
+               const double v= fe.shape_value(i,k);
+               const double u= input[k];
+               const double g= data[k];
+
+               result(i) += dx*(2.*penalty*(u-g)*v - dnv*(u-g) - dnu*v);
+             }
+         }
+      }
+
 /**
  * Flux for the interior penalty method for the Laplacian, namely on
  * the face <i>F</i> the matrices associated with the bilinear form

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.