]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add test
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 5 Jun 2017 14:45:21 +0000 (16:45 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 5 Jun 2017 14:53:02 +0000 (16:53 +0200)
tests/grid/real_to_unit_cell_affine_approx.cc [new file with mode: 0644]
tests/grid/real_to_unit_cell_affine_approx.output [new file with mode: 0644]

diff --git a/tests/grid/real_to_unit_cell_affine_approx.cc b/tests/grid/real_to_unit_cell_affine_approx.cc
new file mode 100644 (file)
index 0000000..36e62ed
--- /dev/null
@@ -0,0 +1,127 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Check TriaAccessor::real_to_unit_cell_affine_approximation
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_in.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/mapping_q.h>
+#include <fstream>
+#include <iomanip>
+
+
+template <int dim, int spacedim>
+void do_test(const Triangulation<dim,spacedim> &tria,
+             const Point<spacedim> &p)
+{
+  MappingQ<dim,spacedim> mapping(1);
+
+  for (typename Triangulation<dim,spacedim>::cell_iterator cell=tria.begin();
+       cell != tria.end(); ++cell)
+    {
+      deallog << "Point p=(" << p << ") on cell with id " << cell->id() << ": "
+              << cell->real_to_unit_cell_affine_approximation(p) << std::endl;
+      Point<dim> mapping_point;
+      try
+        {
+          mapping_point = mapping.transform_real_to_unit_cell(cell, p);
+          mapping_point -= cell->real_to_unit_cell_affine_approximation(p);
+          deallog << "Distance to mapping point: ";
+          for (unsigned int d=0; d<dim; ++d)
+            deallog << mapping_point[d] << " ";
+          deallog << std::endl;
+        }
+      catch (typename Mapping<dim,spacedim>::ExcTransformationFailed)
+        {
+          deallog << "No MappingQ transform possible for this cell and point." << std::endl;
+        }
+    }
+}
+
+
+
+template <int dim>
+void test1()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, -1, 1);
+  tria.refine_global(1);
+  Point<dim> p;
+  for (unsigned int d=0; d<dim; ++d)
+    p[d] = -0.2 + 0.3 * d;
+
+  do_test(tria, p);
+}
+
+
+
+template <int dim>
+void test2()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_ball(tria);
+  Point<dim> p;
+  for (unsigned int d=0; d<dim; ++d)
+    p[d] = -0.4 + 0.5 * d;
+
+  do_test(tria, p);
+}
+
+
+
+template <int dim, int spacedim>
+void test3()
+{
+  Triangulation<dim,spacedim> triangulation;
+  GridIn<dim, spacedim> grid_in;
+  grid_in.attach_triangulation (triangulation);
+  if (dim == 1)
+    {
+      std::ifstream fname(SOURCE_DIR "/../codim_one/grids/circle_1.inp");
+      grid_in.read_ucd (fname);
+    }
+  else
+    {
+      std::ifstream fname(SOURCE_DIR "/../codim_one/grids/sphere_0.inp");
+      grid_in.read_ucd (fname);
+    }
+
+  Point<spacedim> p;
+  for (unsigned int d=0; d<dim; ++d)
+    p[d] = -0.4 + 0.5 * d;
+
+  do_test(triangulation, p);
+}
+
+
+int main()
+{
+  initlog();
+
+  test1<1>();
+  test1<2>();
+  test1<3>();
+
+  test2<2>();
+  test2<3>();
+
+  test3<1,2>();
+  test3<2,3>();
+}
diff --git a/tests/grid/real_to_unit_cell_affine_approx.output b/tests/grid/real_to_unit_cell_affine_approx.output
new file mode 100644 (file)
index 0000000..b4685f8
--- /dev/null
@@ -0,0 +1,91 @@
+
+DEAL::Point p=(-0.200000) on cell with id 0_0:: 0.400000
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.200000) on cell with id 0_1:0: 0.800000
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.200000) on cell with id 0_1:1: -0.200000
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.200000 0.100000) on cell with id 0_0:: 0.400000 0.550000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000) on cell with id 0_1:0: 0.800000 1.10000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000) on cell with id 0_1:1: -0.200000 1.10000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000) on cell with id 0_1:2: 0.800000 0.100000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000) on cell with id 0_1:3: -0.200000 0.100000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_0:: 0.400000 0.550000 0.700000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:0: 0.800000 1.10000 1.40000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:1: -0.200000 1.10000 1.40000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:2: 0.800000 0.100000 1.40000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:3: -0.200000 0.100000 1.40000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:4: 0.800000 1.10000 0.400000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:5: -0.200000 1.10000 0.400000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:6: 0.800000 0.100000 0.400000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.200000 0.100000 0.400000) on cell with id 0_1:7: -0.200000 0.100000 0.400000
+DEAL::Distance to mapping point: 0.00000 0.00000 0.00000 
+DEAL::Point p=(-0.400000 0.100000) on cell with id 0_0:: 0.100000 1.94853
+DEAL::No MappingQ transform possible for this cell and point.
+DEAL::Point p=(-0.400000 0.100000) on cell with id 1_0:: 0.741421 0.600000
+DEAL::Distance to mapping point: -1.11022e-16 0.0250000 
+DEAL::Point p=(-0.400000 0.100000) on cell with id 2_0:: -0.182843 0.670711
+DEAL::Distance to mapping point: 2.77556e-17 0.00000 
+DEAL::Point p=(-0.400000 0.100000) on cell with id 3_0:: 0.600000 2.67279
+DEAL::No MappingQ transform possible for this cell and point.
+DEAL::Point p=(-0.400000 0.100000) on cell with id 4_0:: 1.46569 0.100000
+DEAL::Distance to mapping point: 2.22045e-16 -1.60000 
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 0_0:: -0.446410 0.736603 1.91962
+DEAL::Distance to mapping point: 5.55112e-17 1.11022e-16 -2.22045e-16 
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 1_0:: -0.00717968 0.626795 3.21658
+DEAL::No MappingQ transform possible for this cell and point.
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 2_0:: 0.626795 2.67017 1.26077
+DEAL::No MappingQ transform possible for this cell and point.
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 3_0:: -0.00717968 -0.0618802 0.626795
+DEAL::Distance to mapping point: 0.173846 -4.85723e-17 -0.0434616 
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 4_0:: 0.484530 0.626795 1.26077
+DEAL::Distance to mapping point: 0.00000 -0.00179492 -0.0107695 
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 5_0:: -0.00717968 1.85056 1.26077
+DEAL::No MappingQ transform possible for this cell and point.
+DEAL::Point p=(-0.400000 0.100000 0.600000) on cell with id 6_0:: 1.30415 -0.00717968 1.26077
+DEAL::Distance to mapping point: 2.22045e-16 -1.49282 2.23923 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 0_0:: -0.115537
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 1_0:: 0.119577
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 2_0:: 0.500000
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 3_0:: 0.880423
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 4_0:: 1.11554
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 5_0:: 1.11554
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 6_0:: 0.880423
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 7_0:: 0.500000
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 8_0:: 0.119577
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.00000) on cell with id 9_0:: -0.115537
+DEAL::Distance to mapping point: 0.00000 
+DEAL::Point p=(-0.400000 0.100000 0.00000) on cell with id 0_0:: 0.846410 0.586603
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.400000 0.100000 0.00000) on cell with id 1_0:: 0.586603 0.846410
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.400000 0.100000 0.00000) on cell with id 2_0:: 0.846410 0.500000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.400000 0.100000 0.00000) on cell with id 3_0:: 0.500000 0.846410
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.400000 0.100000 0.00000) on cell with id 4_0:: 0.586603 0.500000
+DEAL::Distance to mapping point: 0.00000 0.00000 
+DEAL::Point p=(-0.400000 0.100000 0.00000) on cell with id 5_0:: 0.500000 0.586603
+DEAL::Distance to mapping point: 0.00000 0.00000 

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.