]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add some new tests.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 29 Sep 2003 14:34:26 +0000 (14:34 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 29 Sep 2003 14:34:26 +0000 (14:34 +0000)
git-svn-id: https://svn.dealii.org/trunk@8047 0785d39b-7218-0410-832d-ea1e28bc413d

tests/bits/Makefile
tests/bits/mesh_3d_10.cc [new file with mode: 0644]
tests/bits/mesh_3d_8.cc [new file with mode: 0644]
tests/bits/mesh_3d_9.cc [new file with mode: 0644]

index d0d52be1fa1c75a40dab50ee77657f83ca4e0038..d7f143d1fd97881f79ea4c9db3b81ff24f0e0c03 100644 (file)
@@ -112,6 +112,7 @@ 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)
+mesh_3d_10.exe          : mesh_3d_10.g.$(OBJEXT)           $(libraries)
 
 
 tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
@@ -142,7 +143,7 @@ tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
        sparse_lu_decomposition_1 block_sparse_matrix_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
+       mesh_3d_6 mesh_3d_7 mesh_3d_8 mesh_3d_9 mesh_3d_10
 
 ############################################################
 
diff --git a/tests/bits/mesh_3d_10.cc b/tests/bits/mesh_3d_10.cc
new file mode 100644 (file)
index 0000000..df60d67
--- /dev/null
@@ -0,0 +1,105 @@
+//----------------------------  mesh_3d_10.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_10.cc  ---------------------------
+
+
+// check that face orientation flags work. at one point in time we
+// mixed up the generation of edges associated with faces for the
+// non-standard case when refining the mesh. like in mesh_3d_10, but
+// make sure that each cell has 8 distinct vertices if looking at the
+// vertices of the faces
+
+#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>
+#include <set>
+
+void check_this (Triangulation<3> &tria)
+{
+  for (Triangulation<3>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      std::set<unsigned int> vertices;
+      for (unsigned int f=0; f<GeometryInfo<3>::faces_per_cell; ++f)
+        {
+          vertices.insert (cell->face(f)->vertex_index(0));
+          vertices.insert (cell->face(f)->vertex_index(1));
+          vertices.insert (cell->face(f)->vertex_index(2));
+          vertices.insert (cell->face(f)->vertex_index(3));
+        }
+      Assert (vertices.size() == 8, ExcInternalError());
+    }
+  deallog << "    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_10.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_8.cc b/tests/bits/mesh_3d_8.cc
new file mode 100644 (file)
index 0000000..9b7be4a
--- /dev/null
@@ -0,0 +1,104 @@
+//----------------------------  mesh_3d_8.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_8.cc  ---------------------------
+
+
+// check that face orientation flags work. at one point in time we
+// mixed up the generation of edges associated with faces for the
+// non-standard case when refining the mesh. thus make sure that each
+// edge only appears once in the list of edges of every cell.
+
+#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>
+#include <set>
+
+void check_this (Triangulation<3> &tria)
+{
+  for (Triangulation<3>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      deallog << "  Cell=" << cell << std::endl;
+      
+      std::set<Triangulation<3>::line_iterator> lines;
+      for (unsigned int l=0; l<GeometryInfo<3>::lines_per_cell; ++l)
+        {
+          Assert (lines.find (cell->line(l)) == lines.end(),
+                  ExcInternalError());
+          lines.insert (cell->line(l));
+          deallog << "    line(" << l << ")=" << cell->line(l) << 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_8.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_9.cc b/tests/bits/mesh_3d_9.cc
new file mode 100644 (file)
index 0000000..e8cf070
--- /dev/null
@@ -0,0 +1,103 @@
+//----------------------------  mesh_3d_9.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_9.cc  ---------------------------
+
+
+// check that face orientation flags work. at one point in time we
+// mixed up the generation of edges associated with faces for the
+// non-standard case when refining the mesh. like in mesh_3d_8, but
+// make sure that each cell has 8 distinct vertices if looking at the
+// vertices of the edges
+
+#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>
+#include <set>
+
+void check_this (Triangulation<3> &tria)
+{
+  for (Triangulation<3>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      std::set<unsigned int> vertices;
+      for (unsigned int l=0; l<GeometryInfo<3>::lines_per_cell; ++l)
+        {
+          vertices.insert (cell->line(l)->vertex_index(0));
+          vertices.insert (cell->line(l)->vertex_index(1));
+        }
+      Assert (vertices.size() == 8, ExcInternalError());
+    }
+  deallog << "    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_9.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.