]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Two tests for generating a parallelepiped grid.
authorToby D. Young <tyoung@ippt.pan.pl>
Sun, 13 Jan 2013 00:10:53 +0000 (00:10 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Sun, 13 Jan 2013 00:10:53 +0000 (00:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@28034 0785d39b-7218-0410-832d-ea1e28bc413d

tests/deal.II/grid_parallelepiped.cc [new file with mode: 0644]
tests/deal.II/grid_parallelepiped/cmp/generic [new file with mode: 0644]
tests/deal.II/grid_parallelepiped_01.cc [new file with mode: 0644]
tests/deal.II/grid_parallelepiped_01/cmp/generic [new file with mode: 0644]

diff --git a/tests/deal.II/grid_parallelepiped.cc b/tests/deal.II/grid_parallelepiped.cc
new file mode 100644 (file)
index 0000000..ee23e49
--- /dev/null
@@ -0,0 +1,83 @@
+//----------------------------------------------------------------------
+//    $Id: grid_parallelepiped.cc 2013-01-07 young $
+//    Version: $Name$ 
+//
+//    Copyright (C) 2013 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.
+//
+//----------------------------------------------------------------------
+
+
+#include "../tests.h"
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/logstream.h>
+#include <cmath>
+#include <cstdlib>
+
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+
+// Output
+std::ofstream logfile ("grid_parallelepiped/output");
+
+// The simplest test case is to create a parallelepiped grid, output
+// the result, and hope for the best.
+template<int dim>
+void check_parallelepiped (bool colorize, bool log)
+{
+  Tensor<2, dim> corners = Tensor<2, dim> ();
+
+  // build corners for this particular dim:
+  switch (dim)
+    {
+    case 1:
+      corners[0] = Point<dim> (0.5);
+      break;
+
+    case 2:
+      corners[0] = Point<dim> (0.0, 0.5);
+      corners[1] = Point<dim> (0.5, 0.0);
+      break;
+      
+    case 3:
+      corners[0] = Point<dim> (0.0, 0.5, 0.5);
+      corners[1] = Point<dim> (0.5, 0.0, 0.5);
+      corners[2] = Point<dim> (0.5, 0.5, 0.0);
+      break;
+
+    default:
+      Assert (false, ExcInternalError ());
+    }
+  
+  Triangulation<dim> triangulation;
+  GridGenerator::parallelepiped (triangulation, corners, colorize);
+  
+  GridOut grid_out;
+
+  if (log)
+    grid_out.write_gnuplot (triangulation, logfile);
+
+  else
+    grid_out.write_gnuplot (triangulation, std::cout);
+
+  triangulation.clear ();
+}
+
+int main ()
+{
+  // Check parallelepiped 
+  check_parallelepiped<1> (false, true);
+  check_parallelepiped<2> (false, true);
+  check_parallelepiped<3> (true,  true);
+}
diff --git a/tests/deal.II/grid_parallelepiped/cmp/generic b/tests/deal.II/grid_parallelepiped/cmp/generic
new file mode 100644 (file)
index 0000000..8034c6e
--- /dev/null
@@ -0,0 +1,34 @@
+0 0 0
+0.5 0 0
+
+0 0 0 0
+0 0.5 0 0
+0.5 0.5 0 0
+0.5 0 0 0
+0 0 0 0
+
+
+0 0 0 0 0
+0 0.5 0.5 0 0
+0.5 1 0.5 0 0
+0.5 0.5 0 0 0
+0 0 0 0 0
+
+0.5 0 0.5 0 0
+0.5 0.5 1 0 0
+1 1 1 0 0
+1 0.5 0.5 0 0
+0.5 0 0.5 0 0
+
+0 0 0 0 0
+0.5 0 0.5 0 0
+
+0 0.5 0.5 0 0
+0.5 0.5 1 0 0
+
+0.5 1 0.5 0 0
+1 1 1 0 0
+
+0.5 0.5 0 0 0
+1 0.5 0.5 0 0
+
diff --git a/tests/deal.II/grid_parallelepiped_01.cc b/tests/deal.II/grid_parallelepiped_01.cc
new file mode 100644 (file)
index 0000000..faaf2b6
--- /dev/null
@@ -0,0 +1,96 @@
+//----------------------------------------------------------------------
+//    $Id: grid_parallelepiped.cc 2013-01-12 young $
+//    Version: $Name$ 
+//
+//    Copyright (C) 2013 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.
+//
+//----------------------------------------------------------------------
+
+
+#include "../tests.h"
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+#include <deal.II/base/logstream.h>
+#include <cmath>
+#include <cstdlib>
+
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+
+// Output
+std::ofstream logfile ("grid_parallelepiped_01/output");
+
+// As sketched in the deal.II docs, the parallelepiped class is just a
+// hyper_rectangle in 1d and a parallelogram in 2d. That can checked
+// by simple comparison to a reference triangulation.
+
+// Here is the implementation in 1d:
+void check_1d_parallelepiped_by_comparison (bool log)
+{
+  // build corners for this particular dim:
+  Tensor<2, 1> corners = Tensor<2, 1> ();
+  corners[0] = Point<1> (0.5);
+  
+  Triangulation<1> triangulation_parallelepiped;
+  GridGenerator::parallelepiped (triangulation_parallelepiped, corners, false);
+  
+  Triangulation<1> triangulation_cube;
+  GridGenerator::hyper_cube (triangulation_cube, 0., 0.5);
+
+  if (log)
+    {
+      logfile << "\ncheck 1d parallelepiped (hyper_cube): ";
+      if (GridTools::have_same_coarse_mesh (triangulation_parallelepiped,
+                                           triangulation_cube))
+       logfile << "OK" << std::endl;
+      else
+       logfile << "not OK... coarse grid are different but they should be the same!" 
+               << std::endl;
+    }
+}
+
+// Here is the implementation in 2d:
+void check_2d_parallelepiped_by_comparison (bool log)
+{
+  // build corners for this particular dim that are known to give the
+  // same output order as parallelogram:
+  Tensor<2, 2> corners = Tensor<2, 2> ();
+  corners[0] = Point<2> (0.0, 0.5);
+  corners[1] = Point<2> (0.5, 0.0);
+  
+  Triangulation<2> triangulation_parallelepiped;
+  GridGenerator::parallelepiped (triangulation_parallelepiped, corners, false);
+  
+  Triangulation<2> triangulation_parallelogram;
+  GridGenerator::parallelogram (triangulation_parallelogram, corners, false);
+
+  if (log)
+    {
+      logfile << "\ncheck 2d parallelepiped (parallelogram): ";
+      if (GridTools::have_same_coarse_mesh (triangulation_parallelepiped,
+                                           triangulation_parallelogram))
+       logfile << "OK" << std::endl;
+      else
+       logfile << "not OK... coarse grid are different but they should be the same!" 
+               << std::endl;
+    }
+}
+
+int main ()
+{
+  // Check parallelepiped 
+  check_1d_parallelepiped_by_comparison (true);
+  check_2d_parallelepiped_by_comparison (true);
+}
diff --git a/tests/deal.II/grid_parallelepiped_01/cmp/generic b/tests/deal.II/grid_parallelepiped_01/cmp/generic
new file mode 100644 (file)
index 0000000..8c27f27
--- /dev/null
@@ -0,0 +1,4 @@
+
+check 1d parallelepiped (hyper_cube): OK
+
+check 2d parallelepiped (parallelogram): OK

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.