]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Lots of new tests about 3d meshes in which faces are not orientable.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 28 Sep 2003 14:27:38 +0000 (14:27 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 28 Sep 2003 14:27:38 +0000 (14:27 +0000)
git-svn-id: https://svn.dealii.org/trunk@8042 0785d39b-7218-0410-832d-ea1e28bc413d

tests/bits/Makefile
tests/bits/coarsening_3d.cc
tests/bits/mesh_3d.h [new file with mode: 0644]
tests/bits/mesh_3d_1.cc
tests/bits/mesh_3d_2.cc [new file with mode: 0644]
tests/bits/mesh_3d_3.cc [new file with mode: 0644]
tests/bits/mesh_3d_4.cc [new file with mode: 0644]
tests/bits/mesh_3d_5.cc [new file with mode: 0644]
tests/bits/mesh_3d_6.cc [new file with mode: 0644]
tests/bits/mesh_3d_7.cc [new file with mode: 0644]

index b9e69344a59f121ab54b6991611f13c752646637..d0d52be1fa1c75a40dab50ee77657f83ca4e0038 100644 (file)
@@ -104,6 +104,14 @@ hyper_ball_3d.exe       : hyper_ball_3d.g.$(OBJEXT)        $(libraries)
 cylinder.exe            : cylinder.g.$(OBJEXT)             $(libraries)
 coarsening_3d.exe       : coarsening_3d.g.$(OBJEXT)        $(libraries)
 mesh_3d_1.exe           : mesh_3d_1.g.$(OBJEXT)            $(libraries)
+mesh_3d_2.exe           : mesh_3d_2.g.$(OBJEXT)            $(libraries)
+mesh_3d_3.exe           : mesh_3d_3.g.$(OBJEXT)            $(libraries)
+mesh_3d_4.exe           : mesh_3d_4.g.$(OBJEXT)            $(libraries)
+mesh_3d_5.exe           : mesh_3d_5.g.$(OBJEXT)            $(libraries)
+mesh_3d_6.exe           : mesh_3d_6.g.$(OBJEXT)            $(libraries)
+mesh_3d_7.exe           : mesh_3d_7.g.$(OBJEXT)            $(libraries)
+mesh_3d_8.exe           : mesh_3d_8.g.$(OBJEXT)            $(libraries)
+mesh_3d_9.exe           : mesh_3d_9.g.$(OBJEXT)            $(libraries)
 
 
 tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
@@ -132,7 +140,9 @@ tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
        parameter_handler_3 parameter_handler_4 \
        parameter_handler_5 \
        sparse_lu_decomposition_1 block_sparse_matrix_1 \
-       hyper_ball_3d cylinder coarsening_3d mesh_3d_1
+       hyper_ball_3d cylinder coarsening_3d \
+       mesh_3d_1 mesh_3d_2 mesh_3d_3 mesh_3d_4 mesh_3d_5 \
+       mesh_3d_6 mesh_3d_7
 
 ############################################################
 
index e45a65440de4bd389b341652172ed738df3abc93..798367920c9a66005c5a0aa0cc9faad9186c3cf2 100644 (file)
@@ -15,6 +15,7 @@
 // this test failed with an internal error somewhere in the coarsening
 // functions
 
+#include "mesh_3d.h"
 
 #include <base/logstream.h>
 #include <grid/tria.h>
 #include <fstream>
 
 
-
-void create_coarse_grid (Triangulation<3> &coarse_grid)
-{
-  std::vector<Point<3> >    vertices;
-  std::vector<CellData<3> > cells;
-  SubCellData               sub_cell_data;
-  
-  const Point<3> outer_points[8] = { Point<3>(-1,0,0),
-                                     Point<3>(-1,-1,0),
-                                     Point<3>(0,-1,0),
-                                     Point<3>(+1,-1,0),
-                                     Point<3>(+1,0,0),
-                                     Point<3>(+1,+1,0),
-                                     Point<3>(0,+1,0),
-                                     Point<3>(-1,+1,0) };
-
-                                   // first the point in the middle
-                                   // and the rest of those on the
-                                   // upper surface
-  vertices.push_back (Point<3>(0,0,0));
-  for (unsigned int i=0; i<7; ++i)
-    vertices.push_back (outer_points[i]);
-
-                                   // then points on lower surface
-  vertices.push_back (Point<3>(0,0,-1));
-  for (unsigned int i=0; i<7; ++i)
-    vertices.push_back (outer_points[i]
-                        +
-                        Point<3>(0,0,-1));
-
-  const unsigned int n_vertices_per_surface = 8;
-  Assert (vertices.size() == n_vertices_per_surface*2,
-          ExcInternalError());
-    
-  const unsigned int connectivity[3][4]
-    = { { 1, 2, 3, 0 },
-        { 3, 4, 5, 0 },
-        { 0, 5, 6, 7 } };
-  for (unsigned int i=0; i<3; ++i)
-    {
-      CellData<3> cell;
-      for (unsigned int j=0; j<4; ++j)
-        {
-          cell.vertices[j]   = connectivity[i][j];
-          cell.vertices[j+4] = connectivity[i][j]+n_vertices_per_surface;
-        }
-      cells.push_back (cell);
-    }
-
-                                   // finally generate a triangulation
-                                   // out of this
-  GridReordering<3>::reorder_cells (cells);
-  coarse_grid.create_triangulation (vertices, cells, sub_cell_data);
-}
-
-
 int main () 
 {
   std::ofstream logfile("coarsening_3d.output");
@@ -89,7 +34,7 @@ int main ()
   deallog.depth_console(0);
 
   Triangulation<3> coarse_grid;
-  create_coarse_grid (coarse_grid);
+  create_L_shape (coarse_grid);
 
                                    // refine once, then unrefine again
   coarse_grid.refine_global (1);
diff --git a/tests/bits/mesh_3d.h b/tests/bits/mesh_3d.h
new file mode 100644 (file)
index 0000000..f92529d
--- /dev/null
@@ -0,0 +1,126 @@
+//----------------------------  mesh_3d.h  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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.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 <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_reordering.h>
+
+
+
+void create_two_cubes (Triangulation<3> &coarse_grid)
+{
+  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;
+    
+  const unsigned int connectivity[2][4]
+    = { { 0,1,2,3 },
+        { 4,5,2,1 } };
+  for (unsigned int i=0; i<2; ++i)
+    {
+      CellData<3> cell;
+      for (unsigned int j=0; j<4; ++j)
+        {
+          cell.vertices[j]   = connectivity[i][j];
+          cell.vertices[j+4] = connectivity[i][j]+6;
+        }
+      cells.push_back (cell);
+    }
+
+                                   // finally generate a triangulation
+                                   // out of this
+  GridReordering<3>::reorder_cells (cells);
+  coarse_grid.create_triangulation (vertices, cells, SubCellData());
+}
+  
+
+
+
+
+void create_L_shape (Triangulation<3> &coarse_grid)
+{
+  std::vector<Point<3> >    vertices;
+  std::vector<CellData<3> > cells;
+  
+  const Point<3> outer_points[8] = { Point<3>(-1,0,0),
+                                     Point<3>(-1,-1,0),
+                                     Point<3>(0,-1,0),
+                                     Point<3>(+1,-1,0),
+                                     Point<3>(+1,0,0),
+                                     Point<3>(+1,+1,0),
+                                     Point<3>(0,+1,0),
+                                     Point<3>(-1,+1,0) };
+
+                                   // first the point in the middle
+                                   // and the rest of those on the
+                                   // upper surface
+  vertices.push_back (Point<3>(0,0,0));
+  for (unsigned int i=0; i<7; ++i)
+    vertices.push_back (outer_points[i]);
+
+                                   // then points on lower surface
+  vertices.push_back (Point<3>(0,0,-1));
+  for (unsigned int i=0; i<7; ++i)
+    vertices.push_back (outer_points[i]
+                        +
+                        Point<3>(0,0,-1));
+
+  const unsigned int n_vertices_per_surface = 8;
+  Assert (vertices.size() == n_vertices_per_surface*2,
+          ExcInternalError());
+    
+  const unsigned int connectivity[3][4]
+    = { { 1, 2, 3, 0 },
+        { 3, 4, 5, 0 },
+        { 0, 5, 6, 7 } };
+  for (unsigned int i=0; i<3; ++i)
+    {
+      CellData<3> cell;
+      for (unsigned int j=0; j<4; ++j)
+        {
+          cell.vertices[j]   = connectivity[i][j];
+          cell.vertices[j+4] = connectivity[i][j]+n_vertices_per_surface;
+        }
+      cells.push_back (cell);
+    }
+
+                                   // finally generate a triangulation
+                                   // out of this
+  GridReordering<3>::reorder_cells (cells);
+  coarse_grid.create_triangulation (vertices, cells, SubCellData());
+}
+
+
+void coarsen_global (Triangulation<3> &grid)
+{
+  for (Triangulation<3>::active_cell_iterator c=grid.begin_active();
+       c != grid.end(); ++c)
+    c->set_coarsen_flag ();
+  grid.execute_coarsening_and_refinement ();
+}
index 4c0b2771c1be69973bbce0f492e63e55d5bbdade..92d5884cedecb17a3e4051b0fa3e1be8654cc151 100644 (file)
 // coarsening_3d), but we can check this fact much earlier already (done
 // here)
 
+#include "mesh_3d.h"
 
 #include <base/logstream.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 <fstream>
 
 
 
-void create_coarse_grid (Triangulation<3> &coarse_grid)
-{
-  std::vector<Point<3> >    vertices;
-  std::vector<CellData<3> > cells;
-  SubCellData               sub_cell_data;
-  
-  const Point<3> outer_points[8] = { Point<3>(-1,0,0),
-                                     Point<3>(-1,-1,0),
-                                     Point<3>(0,-1,0),
-                                     Point<3>(+1,-1,0),
-                                     Point<3>(+1,0,0),
-                                     Point<3>(+1,+1,0),
-                                     Point<3>(0,+1,0),
-                                     Point<3>(-1,+1,0) };
-
-                                   // first the point in the middle
-                                   // and the rest of those on the
-                                   // upper surface
-  vertices.push_back (Point<3>(0,0,0));
-  for (unsigned int i=0; i<7; ++i)
-    vertices.push_back (outer_points[i]);
-
-                                   // then points on lower surface
-  vertices.push_back (Point<3>(0,0,-1));
-  for (unsigned int i=0; i<7; ++i)
-    vertices.push_back (outer_points[i]
-                        +
-                        Point<3>(0,0,-1));
-
-  const unsigned int n_vertices_per_surface = 8;
-  Assert (vertices.size() == n_vertices_per_surface*2,
-          ExcInternalError());
-    
-  const unsigned int connectivity[3][4]
-    = { { 1, 2, 3, 0 },
-        { 3, 4, 5, 0 },
-        { 0, 5, 6, 7 } };
-  for (unsigned int i=0; i<3; ++i)
-    {
-      CellData<3> cell;
-      for (unsigned int j=0; j<4; ++j)
-        {
-          cell.vertices[j]   = connectivity[i][j];
-          cell.vertices[j+4] = connectivity[i][j]+n_vertices_per_surface;
-        }
-      cells.push_back (cell);
-    }
-
-                                   // finally generate a triangulation
-                                   // out of this
-  GridReordering<3>::reorder_cells (cells);
-  coarse_grid.create_triangulation (vertices, cells, sub_cell_data);
-}
-
-
 int main () 
 {
   std::ofstream logfile("mesh_3d_1.output");
@@ -93,7 +38,7 @@ int main ()
   deallog.depth_console(0);
 
   Triangulation<3> coarse_grid;
-  create_coarse_grid (coarse_grid);
+  create_L_shape (coarse_grid);
 
                                    // output all lines and faces
   for (Triangulation<3>::active_cell_iterator cell=coarse_grid.begin_active();
diff --git a/tests/bits/mesh_3d_2.cc b/tests/bits/mesh_3d_2.cc
new file mode 100644 (file)
index 0000000..5817ae2
--- /dev/null
@@ -0,0 +1,75 @@
+//----------------------------  mesh_3d_2.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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_2.cc  ---------------------------
+
+
+// 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
+//
+// this test just checks that such a mesh can be created. further
+// tests check that we are still producing consistent states with this
+
+#include "mesh_3d.h"
+
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_reordering.h>
+
+#include <fstream>
+
+
+int main () 
+{
+  std::ofstream logfile("mesh_3d_2.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  Triangulation<3> coarse_grid;
+  create_two_cubes (coarse_grid);
+
+                                   // output all lines and faces
+  for (Triangulation<3>::active_cell_iterator cell=coarse_grid.begin_active();
+       cell != coarse_grid.end(); ++cell)
+    {
+      deallog << "Cell = " << cell << std::endl;
+      for (unsigned int i=0; i<GeometryInfo<3>::lines_per_cell; ++i)
+        deallog << "    Line = " << cell->line(i)
+                << " : " << cell->line(i)->vertex_index(0)
+                << " -> " << cell->line(i)->vertex_index(1)
+                << std::endl;
+
+      for (unsigned int i=0; i<GeometryInfo<3>::quads_per_cell; ++i)
+        deallog << "    Quad = " << cell->quad(i)
+                << " : " << cell->quad(i)->vertex_index(0)
+                << " -> " << cell->quad(i)->vertex_index(1)
+                << " -> " << cell->quad(i)->vertex_index(2)
+                << " -> " << cell->quad(i)->vertex_index(3)
+                << std::endl
+                << "           orientation = "
+                << (cell->get_face_orientation(i) ? "true" : "false")
+                << std::endl;
+    }
+
+                                   // we know that from the second
+                                   // cell, the common face must have
+                                   // wrong orientation. check this
+  Assert ((++coarse_grid.begin_active())->get_face_orientation(4)
+          == false,
+          ExcInternalError());
+}
+
+  
+  
diff --git a/tests/bits/mesh_3d_3.cc b/tests/bits/mesh_3d_3.cc
new file mode 100644 (file)
index 0000000..9d19667
--- /dev/null
@@ -0,0 +1,67 @@
+//----------------------------  mesh_3d_3.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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_3.cc  ---------------------------
+
+
+// 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
+//
+// for this grid, check that vertex numbers still match up
+
+#include "mesh_3d.h"
+
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_reordering.h>
+
+#include <fstream>
+
+
+int main () 
+{
+  std::ofstream logfile("mesh_3d_3.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  Triangulation<3> coarse_grid;
+  create_two_cubes (coarse_grid);
+
+  const Triangulation<3>::active_cell_iterator
+    cells[2] = { coarse_grid.begin_active(),
+                 ++coarse_grid.begin_active() };
+
+                                   // output all vertices
+  for (unsigned int c=0; c<2; ++c)
+    for (unsigned int v=0; v<GeometryInfo<3>::vertices_per_cell; ++v)
+      deallog << "Cell " << c << ", vertex " << v
+              << ": " << cells[c]->vertex_index(v)
+              << "  @  " << cells[c]->vertex(v)
+              << std::endl;
+
+                                   // make sure by hand that certain
+                                   // vertices match up
+  Assert (cells[0]->vertex(1) == cells[1]->vertex(3),
+          ExcInternalError());
+  Assert (cells[0]->vertex(2) == cells[1]->vertex(2),
+          ExcInternalError());
+  Assert (cells[0]->vertex(5) == cells[1]->vertex(7),
+          ExcInternalError());
+  Assert (cells[0]->vertex(6) == cells[1]->vertex(6),
+          ExcInternalError());
+}
+
+  
+  
diff --git a/tests/bits/mesh_3d_4.cc b/tests/bits/mesh_3d_4.cc
new file mode 100644 (file)
index 0000000..f0f3ca5
--- /dev/null
@@ -0,0 +1,104 @@
+//----------------------------  mesh_3d_4.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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_4.cc  ---------------------------
+
+
+// check that face orientation flags are properly inherited by
+// counting them
+
+#include "mesh_3d.h"
+
+#include <base/logstream.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 <fstream>
+
+
+unsigned int count_wrong_faces (const Triangulation<3> &tria)
+{
+  unsigned int count = 0;
+  
+                                   // count faces with "wrong"
+                                   // orientation
+  for (Triangulation<3>::active_cell_iterator cell=tria.begin_active();
+       cell != tria.end(); ++cell)
+    for (unsigned int f=0; f<GeometryInfo<3>::faces_per_cell; ++f)
+      if (cell->get_face_orientation(f) == false)
+        ++count;
+  return count;
+}
+
+  
+
+void check (Triangulation<3> &tria)
+{
+  const unsigned int initial_count = count_wrong_faces (tria);
+  for (unsigned int r=0; r<3; ++r)
+    {
+      tria.refine_global (1);
+      const unsigned int count = count_wrong_faces (tria);
+      deallog << "'Wrong' faces = " << count << std::endl;
+      Assert (count == initial_count * (4<<(2*r)),
+              ExcInternalError());
+    }
+  
+  {
+    coarsen_global (tria);
+    const unsigned int count = count_wrong_faces (tria);
+    deallog << "'Wrong' faces = " << count << std::endl;
+    Assert (count == initial_count * (4<<(2*1)),
+            ExcInternalError());
+  }
+  
+  {
+    tria.refine_global (1);
+    const unsigned int count = count_wrong_faces (tria);
+    deallog << "'Wrong' faces = " << count << std::endl;
+    Assert (count == initial_count * (4<<(2*2)),
+            ExcInternalError());
+  }
+}
+
+  
+
+int main () 
+{
+  std::ofstream logfile("mesh_3d_4.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  {  
+    Triangulation<3> coarse_grid;
+    create_two_cubes (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    create_L_shape (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    GridGenerator::hyper_ball (coarse_grid);
+    check (coarse_grid);
+  }
+  
+}
+
+  
+  
diff --git a/tests/bits/mesh_3d_5.cc b/tests/bits/mesh_3d_5.cc
new file mode 100644 (file)
index 0000000..419ed55
--- /dev/null
@@ -0,0 +1,101 @@
+//----------------------------  mesh_3d_5.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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_5.cc  ---------------------------
+
+
+// check that face orientation flags are properly inherited by
+// looking at them and checking their children
+
+#include "mesh_3d.h"
+
+#include <base/logstream.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 <fstream>
+
+
+void check_this (Triangulation<3> &tria)
+{
+                                   // look at all faces, not only
+                                   // active ones
+  for (Triangulation<3>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    for (unsigned int f=0; f<GeometryInfo<3>::faces_per_cell; ++f)
+      if (cell->has_children())
+        for (unsigned int c=0; c<GeometryInfo<3>::subfaces_per_face; ++c)
+          {
+            Assert (cell->get_face_orientation(f) ==
+                    cell->child(GeometryInfo<3>::child_cell_on_face(f,c))
+                    ->get_face_orientation(f),
+                    ExcInternalError());
+            deallog << "Cell << " << cell
+                    << ", face " << f
+                    << " subface " << c << " is ok."
+                    << std::endl;
+          }
+}
+
+
+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_5.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  {  
+    Triangulation<3> coarse_grid;
+    create_two_cubes (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    create_L_shape (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    GridGenerator::hyper_ball (coarse_grid);
+    check (coarse_grid);
+  }
+  
+}
+
+  
+  
diff --git a/tests/bits/mesh_3d_6.cc b/tests/bits/mesh_3d_6.cc
new file mode 100644 (file)
index 0000000..b912f09
--- /dev/null
@@ -0,0 +1,118 @@
+//----------------------------  mesh_3d_6.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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_6.cc  ---------------------------
+
+
+// check that face orientation flags work by looping over all cells
+// and check on all faces that if we look from both sides that normal
+// vectors point in opposite directions
+
+#include "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 <fstream>
+
+
+void check_this (Triangulation<3> &tria)
+{
+  QMidpoint<2> q;
+  FE_Q<3> fe(1);
+  FEFaceValues<3> fe_face_values1 (fe, q, update_normal_vectors);
+  FEFaceValues<3> fe_face_values2 (fe, q, update_normal_vectors);
+
+  DoFHandler<3> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+  
+                                   // 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);
+
+          deallog << "Cell " << cell << ", face " << f
+                  << " << n=" << fe_face_values1.normal_vector(0)
+                  << " << nx=" << fe_face_values2.normal_vector(0)
+                  << std::endl;
+
+          Assert (fe_face_values1.normal_vector(0) ==
+                  - fe_face_values2.normal_vector(0),
+                  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_6.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  {  
+    Triangulation<3> coarse_grid;
+    create_two_cubes (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    create_L_shape (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    GridGenerator::hyper_ball (coarse_grid);
+    check (coarse_grid);
+  }
+  
+}
+
+  
+  
diff --git a/tests/bits/mesh_3d_7.cc b/tests/bits/mesh_3d_7.cc
new file mode 100644 (file)
index 0000000..4e35cf2
--- /dev/null
@@ -0,0 +1,128 @@
+//----------------------------  mesh_3d_7.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2003 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_7.cc  ---------------------------
+
+
+// check that face orientation flags work by looping over all cells
+// and check on all faces that quadrature points match up
+
+#include "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 <fstream>
+
+
+void check_this (Triangulation<3> &tria)
+{
+  QGauss3<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);
+
+                                   // 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);
+
+          deallog << "Cell " << cell << ", face " << f
+                  << std::endl;
+
+          for (unsigned int q=0; q<quadrature.n_quadrature_points; ++q)
+            {
+              deallog << "  " << 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),
+                      ExcInternalError());
+
+              Assert (fe_face_values1.JxW(q) ==
+                      fe_face_values2.JxW(q),
+                      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_7.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  {  
+    Triangulation<3> coarse_grid;
+    create_two_cubes (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    create_L_shape (coarse_grid);
+    check (coarse_grid);
+  }
+  
+  {  
+    Triangulation<3> coarse_grid;
+    GridGenerator::hyper_ball (coarse_grid);
+    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.