]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
FEValues with mapping Q1, now uses direction_flag to return the
authorpauletti <pauletti@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 9 Dec 2010 17:06:19 +0000 (17:06 +0000)
committerpauletti <pauletti@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 9 Dec 2010 17:06:19 +0000 (17:06 +0000)
oriented normal for meshes of dim<spacedim

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

deal.II/include/deal.II/fe/fe_update_flags.h
deal.II/source/fe/fe_values.cc
deal.II/source/fe/mapping_q1.cc

index eabd2e4e9e7e32e1543100fd776a22ef2157e7be..e292ac9795941e31e566b54c2ab512e3a3f29c5b 100644 (file)
@@ -340,7 +340,8 @@ operator &= (UpdateFlags &f1, UpdateFlags f2)
  * to the previously visited cell. This information is used for reusing data
  * when calling the method FEValues::reinit() (like derivatives, which do
  * not change if one cell is just a translation of the previous). Currently,
- * this variable does only recognize a translation. However, this concept
+ * this variable does only recognize a translation and a inverted traslation
+ * (if dim<spacedim). However, this concept
  * makes it easy to add additional staties to be detected in
  * FEValues/FEFaceValues for making use of these similarities as well.
  */
@@ -350,6 +351,7 @@ namespace CellSimilarity
     {
       none, 
       translation,
+      inverted_translation,
       invalid_next_cell
     };
 }
index d09ea748bda784f4cc9f73cc7e55929e263cb022..c5f3a2a6a9e578b25cdbbe8f333a1a5b1f49d769 100644 (file)
@@ -3263,6 +3263,12 @@ FEValuesBase<dim,spacedim>::check_cell_similarity
                         :
                         CellSimilarity::none);
 
+  if ( (dim<spacedim) &&  (cell_similarity == CellSimilarity::translation) ) {
+    if (static_cast<const typename Triangulation<dim,spacedim>::cell_iterator &>
+       (*this->present_cell)->direction_flag()
+       != cell->direction_flag() )
+      cell_similarity =  CellSimilarity::inverted_translation;
+  }
                                   // TODO: here, one could implement other
                                   // checks for similarity, e.g. for
                                   // children of a parallelogram.
index a842c3ecc7c9b17e8ad4ed37f2705b40df5d21a8..5f13f7f0ad4d949637c112eb31f63cf11a6d1fa5 100644 (file)
@@ -774,49 +774,60 @@ MappingQ1<dim,spacedim>::fill_fe_values (
             ExcDimensionMismatch(normal_vectors.size(), n_q_points));
 
       if (cell_similarity != CellSimilarity::translation)
-       for (unsigned int point=0; point<n_q_points; ++point)
-         {
-
-           if (dim==spacedim)
-             JxW_values[point]
-               = determinant(data.contravariant[point])*weights[point];
-           else if ( (dim==1) && (spacedim==2) )
-             {
-               data.contravariant[point]=transpose(data.contravariant[point]);
+       for (unsigned int point=0; point<n_q_points; ++point) {
+
+         if (dim==spacedim)
+           JxW_values[point]
+             = determinant(data.contravariant[point])*weights[point];
+           
+         else {
+           if (cell_similarity == CellSimilarity::inverted_translation) {
+             // we only need to flip the normal
+             if(update_flags & update_normal_vectors)
+               normal_vectors[point] *= -1.;
+           }
+           else { 
+             if ( (dim==1) && (spacedim==2) ) {
+               data.contravariant[point]=transpose(data.contravariant[point]);
                JxW_values[point]
                  = data.contravariant[point][0].norm()*weights[point];
-               if(update_flags & update_normal_vectors)
-                 {
-                   normal_vectors[point][0]
-                     = -(data.contravariant[point][0][1]
-                         /
-                         data.contravariant[point][0].norm());
-                   normal_vectors[point][1]
-                     = (data.contravariant[point][0][0]
-                        /
-                        data.contravariant[point][0].norm());
-                 }
+               if(update_flags & update_normal_vectors) {
+                 normal_vectors[point][0]
+                   = -(data.contravariant[point][0][1]
+                       /
+                       data.contravariant[point][0].norm());
+                 normal_vectors[point][1]
+                   = (data.contravariant[point][0][0]
+                      /
+                      data.contravariant[point][0].norm());
+                 if (!cell->direction_flag())
+                   normal_vectors[point] *= -1.;
+               }
              }
-           else if ( (dim==2) && (spacedim==3) )
-             {
-               data.contravariant[point]=transpose(data.contravariant[point]);
-               cross_product(data.contravariant[point][2],
-                             data.contravariant[point][0],
-                             data.contravariant[point][1]);
-               JxW_values[point]
-                 = data.contravariant[point][2].norm()*weights[point];
-                                         //the cell normal vector
-                                         //(normal to the surface)
-                                          //is stored in the 3d
-                                         //subtensor of the contravariant tensor
-               data.contravariant[point][2] /= data.contravariant[point][2].norm();
-               if(update_flags & update_normal_vectors)
-                 normal_vectors[point]=data.contravariant[point][2];
+             else { 
+               if ( (dim==2) && (spacedim==3) ) {
+                 data.contravariant[point]=transpose(data.contravariant[point]);
+                 cross_product(data.contravariant[point][2],
+                               data.contravariant[point][0],
+                               data.contravariant[point][1]);
+                 JxW_values[point]
+                   = data.contravariant[point][2].norm()*weights[point];
+                 //the cell normal vector
+                 //(normal to the surface)
+                 //is stored in the 3d
+                 //subtensor of the contravariant tensor
+                 data.contravariant[point][2] /= data.contravariant[point][2].norm();
+                 if(update_flags & update_normal_vectors){
+                   normal_vectors[point]=data.contravariant[point][2];
+                   if (!cell->direction_flag())
+                     normal_vectors[point] *= -1.;
+                 }
+               }
              }
-
+           }
          }
+       }
     }
-
                                   // copy values from InternalData to vector
                                   // given by reference
   if (update_flags & update_jacobians)

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.