]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New test: quadrature points reordering if the face is rotated against the standard...
authorTobias Leicht <tobias.leicht@dlr.de>
Fri, 26 Jan 2007 08:48:02 +0000 (08:48 +0000)
committerTobias Leicht <tobias.leicht@dlr.de>
Fri, 26 Jan 2007 08:48:02 +0000 (08:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@14381 0785d39b-7218-0410-832d-ea1e28bc413d

tests/bits/mesh_3d.h
tests/fail/mesh_3d_20.cc [new file with mode: 0644]

index 7fd097686b39ad8e23e19457f01be682227c8d9f..ef61093003f4819e2f78b1f7f95ea3ed47f10862 100644 (file)
 //
 //----------------------------  mesh_3d.h  ---------------------------
 
-// generate two cubes that are attached to each other in a way so that
-// the edges are all ok, but the normals of the common face don't
-// match up for the standard orientation of the normals. we thus have
-// to store the face orientation in each cell
-
 #include "../tests.h"
 #include <grid/tria.h>
 #include <grid/tria_accessor.h>
 #include <grid/grid_reordering.h>
 
 
+// generate two cubes that are attached to each other in a way so that
+// the edges are all ok, but the normals of the common face don't
+// match up for the standard orientation of the normals. we thus have
+// to store the face orientation in each cell
 
 void create_two_cubes (Triangulation<3> &coarse_grid)
 {
@@ -59,8 +58,52 @@ void create_two_cubes (Triangulation<3> &coarse_grid)
   GridReordering<3>::reorder_cells (cells);
   coarse_grid.create_triangulation_compatibility (vertices, cells, SubCellData());
 }
+
+
+
+// generate two cubes that are attached to each other in a way so that
+// the edges are not all ok and the common face is rotated. we thus have
+// to store the face rotation (and face flip) in each cell
+
+void create_two_cubes_rotation (Triangulation<3> &coarse_grid,
+                               const unsigned int n_rotations)
+{
+  Assert(n_rotations<4, ExcNotImplemented());
   
+  const Point<3> points[6] = { Point<3>(0,0,0),
+                               Point<3>(1,0,0),
+                               Point<3>(1,1,0),
+                               Point<3>(0,1,0),
+                               Point<3>(2,0,0),
+                               Point<3>(2,1,0)};
+  std::vector<Point<3> > vertices;
+  for (unsigned int i=0; i<6; ++i)
+    vertices.push_back (points[i]);
+  for (unsigned int i=0; i<6; ++i)
+    vertices.push_back (points[i] + Point<3>(0,0,-1));
 
+  std::vector<CellData<3> > cells(2);
+    
+  const unsigned int connectivity[5][8]
+                                    // first row: left cube
+                                    // second row: right cube, standard_orientation
+                                    // third row: right cube, rotated once
+                                    // forth row: right cube, rotated twice
+                                    // fifth row: right cube, rotated three times
+    = { { 0,1,2,3,6,7,8,9 },
+        { 1,4,5,2,7,10,11,8},
+        { 2,5,11,8,1,4,10,7},
+        { 8,11,10,7,2,5,4,1},
+        { 7,10,4,1,8,11,5,2} };
+  for (unsigned int j=0; j<8; ++j)
+    {
+      cells[0].vertices[j]   = connectivity[0][j];
+      cells[1].vertices[j]   = connectivity[1+n_rotations][j];
+    }
+                                   // finally generate a triangulation
+                                   // out of this
+  coarse_grid.create_triangulation_compatibility (vertices, cells, SubCellData());
+}
 
 
 
diff --git a/tests/fail/mesh_3d_20.cc b/tests/fail/mesh_3d_20.cc
new file mode 100644 (file)
index 0000000..6296c91
--- /dev/null
@@ -0,0 +1,147 @@
+//----------------------------  mesh_3d_20.cc  ---------------------------
+//    $Id: mesh_3d_7.cc 11749 2005-11-09 19:11:20Z wolf $
+//    Version: $Name$ 
+//
+//    Copyright (C) 2003, 2004, 2005 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  mesh_3d_20.cc  ---------------------------
+
+
+// check that face_rotation and face_flip flags work by looping over all cells
+// and check on all faces that quadrature points match up. based on mesh_3d_7
+
+#include "../tests.h"
+#include "../bits/mesh_3d.h"
+
+#include <base/logstream.h>
+#include <base/quadrature_lib.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_reordering.h>
+#include <grid/grid_generator.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+#include <numerics/data_out.h>
+#include <lac/vector.h>
+
+#include <fstream>
+#include <iostream>
+
+void check_this (Triangulation<3> &tria)
+{
+  QTrapez<2> quadrature;
+  FE_Q<3> fe(1);
+  FEFaceValues<3> fe_face_values1 (fe, quadrature,
+                                   update_q_points | update_JxW_values);
+  FEFaceValues<3> fe_face_values2 (fe, quadrature,
+                                   update_q_points | update_JxW_values);
+  
+  DoFHandler<3> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+
+  unsigned int global_datum = 0;
+
+                                   // look at all faces, not only
+                                   // active ones
+  for (DoFHandler<3>::cell_iterator cell=dof_handler.begin();
+       cell != dof_handler.end(); ++cell)
+    for (unsigned int f=0; f<GeometryInfo<3>::faces_per_cell; ++f)
+      if (!cell->at_boundary(f))
+        {
+          const unsigned int nn
+            = cell->neighbor_of_neighbor (f);
+          fe_face_values1.reinit (cell, f);
+          fe_face_values2.reinit (cell->neighbor(f), nn);
+
+         for (unsigned int q=0; q<quadrature.n_quadrature_points; ++q)
+            {
+                                               // in order to reduce
+                                               // output file size,
+                                               // only write every
+                                               // 289th datum. if the
+                                               // values differ
+                                               // anyway, then the
+                                               // assertion below will
+                                               // catch this, and if
+                                               // we compute _all_
+                                               // values wrongly, then
+                                               // outputting some will
+                                               // be ok, I guess
+              if (global_datum++ % 17*17 == 0)
+               deallog << "Cell " << cell << ", face " << f
+                       << std::endl
+                       << "  " << fe_face_values1.quadrature_point(q)
+                       << ", " << fe_face_values1.JxW(q)
+                       << std::endl;
+              
+               Assert ((fe_face_values1.quadrature_point(q)-
+                        fe_face_values2.quadrature_point(q)).square()
+                       < 1e-20,
+                       ExcInternalError());
+
+               Assert (std::fabs(fe_face_values1.JxW(q)-
+                                 fe_face_values2.JxW(q)) < 1e-15,
+                       ExcInternalError());
+            }
+        }
+}
+
+
+void check (Triangulation<3> &tria)
+{
+  deallog << "Initial check" << std::endl;
+  check_this (tria);
+  
+   for (unsigned int r=0; r<3; ++r)
+     {
+       tria.refine_global (1);
+       deallog << "Check " << r << std::endl;
+       check_this (tria);
+     }
+
+   coarsen_global (tria);
+   deallog << "Check " << 1 << std::endl;
+   check_this (tria);
+
+   tria.refine_global (1);
+   deallog << "Check " << 2 << std::endl;
+   check_this (tria);
+}
+
+
+int main () 
+{
+  std::ofstream logfile("mesh_3d_20/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  {
+    Triangulation<3> coarse_grid;
+    create_two_cubes_rotation (coarse_grid,1);
+    check (coarse_grid);
+  }
+
+  {
+    Triangulation<3> coarse_grid;
+    create_two_cubes_rotation (coarse_grid,2);
+    check (coarse_grid);
+  }
+
+  {
+    Triangulation<3> coarse_grid;
+    create_two_cubes_rotation (coarse_grid,3);
+    check (coarse_grid);
+  }
+  
+}
+
+  
+  

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.