]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added test with MappingQEulerian.
authorLuca Heltai <luca.heltai@sissa.it>
Tue, 19 Sep 2017 13:22:45 +0000 (15:22 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 27 Sep 2017 13:07:33 +0000 (15:07 +0200)
tests/grid/grid_tools_extract_used_vertices_03.cc [new file with mode: 0644]
tests/grid/grid_tools_extract_used_vertices_03.output [new file with mode: 0644]

diff --git a/tests/grid/grid_tools_extract_used_vertices_03.cc b/tests/grid/grid_tools_extract_used_vertices_03.cc
new file mode 100644 (file)
index 0000000..f6d4ada
--- /dev/null
@@ -0,0 +1,72 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 2015 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 extract used_vertices and find_closest_vertex when a mapping
+// is present
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/mapping_q_eulerian.h>
+#include <deal.II/dofs/dof_handler.h>
+
+
+template <int dim, int spacedim>
+void test (const Point<spacedim> &p)
+{
+  deallog << "dim: " << dim << ", spacedim: "
+          << spacedim << std::endl;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(1);
+
+  // Vector fe
+  FESystem<dim,spacedim> fe(FE_Q<dim, spacedim>(1)^spacedim);
+  DoFHandler<dim,spacedim> dh(tria);
+  dh.distribute_dofs(fe);
+
+  // Displacement vector.
+  Vector<double> disp(dh.n_dofs());
+  for (unsigned int i=0; i<dh.n_dofs(); ++i)
+    disp[i] = .5;
+
+  // mapping
+  MappingQEulerian<dim, Vector<double>, spacedim> mapping(1, dh, disp);
+
+  auto m = GridTools::extract_used_vertices(tria, mapping);
+
+  for (auto &e: m)
+    deallog << "Vertex: " << e.first << ": " << e.second << std::endl;
+
+  auto i = GridTools::find_closest_vertex(m,p);
+  deallog << "Closest vertex to " << p
+          << ", v["<< i << "] :"
+          << m[i] << std::endl;
+};
+
+
+int main ()
+{
+  initlog();
+  test<2,2> (Point<2>(.2,.2));
+  test<2,2> (Point<2>(.6,.9));
+  return 0;
+}
+
diff --git a/tests/grid/grid_tools_extract_used_vertices_03.output b/tests/grid/grid_tools_extract_used_vertices_03.output
new file mode 100644 (file)
index 0000000..0331bd7
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL::dim: 2, spacedim: 2
+DEAL::Vertex: 0: 0.500000 0.500000
+DEAL::Vertex: 1: 1.50000 0.500000
+DEAL::Vertex: 2: 0.500000 1.50000
+DEAL::Vertex: 3: 1.50000 1.50000
+DEAL::Vertex: 4: 1.00000 0.500000
+DEAL::Vertex: 5: 0.500000 1.00000
+DEAL::Vertex: 6: 1.50000 1.00000
+DEAL::Vertex: 7: 1.00000 1.50000
+DEAL::Vertex: 8: 1.00000 1.00000
+DEAL::Closest vertex to 0.200000 0.200000, v[0] :0.500000 0.500000
+DEAL::dim: 2, spacedim: 2
+DEAL::Vertex: 0: 0.500000 0.500000
+DEAL::Vertex: 1: 1.50000 0.500000
+DEAL::Vertex: 2: 0.500000 1.50000
+DEAL::Vertex: 3: 1.50000 1.50000
+DEAL::Vertex: 4: 1.00000 0.500000
+DEAL::Vertex: 5: 0.500000 1.00000
+DEAL::Vertex: 6: 1.50000 1.00000
+DEAL::Vertex: 7: 1.00000 1.50000
+DEAL::Vertex: 8: 1.00000 1.00000
+DEAL::Closest vertex to 0.600000 0.900000, v[5] :0.500000 1.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.