fix tangentil jump

git-svn-id: https://svn.dealii.org/trunk@31899 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/include/deal.II/integrators/laplace.h b/deal.II/include/deal.II/integrators/laplace.h
index 61a58afbb89950196eb2e681f44e560cd6456ad8..0148bbcaa5fd27d0f2d7fb72e9eba876ae557e9f 100644 (file)
--- a/deal.II/include/deal.II/integrators/laplace.h
+++ b/deal.II/include/deal.II/integrators/laplace.h
@@ -405,28 +405,40 @@ namespace LocalIntegrators
                  double v1dotn = 0.;
                  double u2dotn = 0.;
                  double v2dotn = 0.;
-                 
+
+                 double ngradu1n = 0.;
+                 double ngradv1n = 0.;
+                 double ngradu2n = 0.;
+                 double ngradv2n = 0.;
+
                  for (unsigned int d=0; d<dim; ++d)
                    {
                      u1dotn += n(d)*fe1.shape_value_component(j,k,d);
                      v1dotn += n(d)*fe1.shape_value_component(i,k,d);
                      u2dotn += n(d)*fe2.shape_value_component(j,k,d);
                      v2dotn += n(d)*fe2.shape_value_component(i,k,d);
+
+                     ngradu1n += n*fe1.shape_grad_component(j,k,d)*n(d);
+                     ngradv1n += n*fe1.shape_grad_component(i,k,d)*n(d);
+                     ngradu2n += n*fe2.shape_grad_component(j,k,d)*n(d);
+                     ngradv2n += n*fe2.shape_grad_component(i,k,d)*n(d);
                    }
                          
                  for (unsigned int d=0; d<fe1.get_fe().n_components(); ++d)
                    {
+                      const double vi = fe1.shape_value_component(i,k,d)-v1dotn*n(d);
+                     const double dnvi = n * fe1.shape_grad_component(i,k,d)-ngradv1n*n(d);
                      
-                                                      // multiply by 
-                      const double vi = fe1.shape_value_component(i,k,d)*(1-v1dotn);
-                      const double dnvi = n * fe1.shape_grad_component(i,k,d)*(1-v1dotn);
-                      const double ve = fe2.shape_value_component(i,k,d)*(1-v2dotn);
-                      const double dnve = n * fe2.shape_grad_component(i,k,d)*(1-v2dotn);
-                      const double ui = fe1.shape_value_component(j,k,d)*(1-u1dotn);
-                      const double dnui = n * fe1.shape_grad_component(j,k,d)*(1-u1dotn);
-                      const double ue = fe2.shape_value_component(j,k,d)*(1-u2dotn);
-                      const double dnue = n * fe2.shape_grad_component(j,k,d)*(1-u2dotn);
-                      M11(i,j) += dx*(-.5*nui*dnvi*ui-.5*nui*dnui*vi+nu*penalty*ui*vi);
+                      const double ve = fe2.shape_value_component(i,k,d)-v2dotn*n(d);
+                      const double dnve = n * fe2.shape_grad_component(i,k,d)-ngradv2n*n(d);
+                     
+                      const double ui = fe1.shape_value_component(j,k,d)-u1dotn*n(d);
+                      const double dnui = n * fe1.shape_grad_component(j,k,d)-ngradu1n*n(d);
+                     
+                      const double ue = fe2.shape_value_component(j,k,d)-u2dotn*n(d);
+                      const double dnue = n * fe2.shape_grad_component(j,k,d)-ngradu2n*n(d);
+
+                     M11(i,j) += dx*(-.5*nui*dnvi*ui-.5*nui*dnui*vi+nu*penalty*ui*vi);
                       M12(i,j) += dx*( .5*nui*dnvi*ue-.5*nue*dnue*vi-nu*penalty*vi*ue);
                       M21(i,j) += dx*(-.5*nue*dnve*ui+.5*nui*dnui*ve-nu*penalty*ui*ve);
                       M22(i,j) += dx*( .5*nue*dnve*ue+.5*nue*dnue*ve+nu*penalty*ue*ve);