]> https://gitweb.dealii.org/ - dealii.git/commitdiff
VectorTools::get_position_vector(): add mapping as argument 11880/head
authorPeter Munch <peterrmuench@gmail.com>
Wed, 10 Mar 2021 21:02:19 +0000 (22:02 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Thu, 11 Mar 2021 10:36:37 +0000 (11:36 +0100)
doc/news/changes/minor/20210310Munch [new file with mode: 0644]
include/deal.II/numerics/vector_tools_interpolate.h
include/deal.II/numerics/vector_tools_interpolate.templates.h
tests/vector_tools/get_position_vector_01.cc [new file with mode: 0644]
tests/vector_tools/get_position_vector_01.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20210310Munch b/doc/news/changes/minor/20210310Munch
new file mode 100644 (file)
index 0000000..9d08683
--- /dev/null
@@ -0,0 +1,5 @@
+Improved: The function VectorTools::get_position_vector() can now also take 
+a Mapping object as input parameter.
+<br>
+(Peter Munch, 2021/03/10)
+
index b13b173a5024a9e07a3794afb72981cc4d4290c9..8c5b60ee585d63903543758df417a68893191f7f 100644 (file)
@@ -284,6 +284,21 @@ namespace VectorTools
                       VectorType &                     vector,
                       const ComponentMask &            mask = ComponentMask());
 
+  /**
+   * Like the above function but also taking @p mapping as argument.
+   * This will introduce an additional approximation between the true geometry
+   * specified by the manifold if the degree of the mapping is lower than the
+   * degree of the finite finite element in the DoFHandler @p dh, but more
+   * importantly it allows to fill location vectors for mappings that do not
+   * preserve vertex locations (like Eulerian mappings).
+   */
+  template <int dim, int spacedim, typename VectorType>
+  void
+  get_position_vector(const Mapping<dim, spacedim> &   mapping,
+                      const DoFHandler<dim, spacedim> &dh,
+                      VectorType &                     vector,
+                      const ComponentMask &            mask = ComponentMask());
+
   //@}
 } // namespace VectorTools
 
index dfde4a80c6b41759abcea1ec95c5471f77b93d9a..c7e73c1f14dc80d0f928554f185796af41e9f72b 100644 (file)
@@ -618,11 +618,30 @@ namespace VectorTools
   }
 
 
+
   template <int dim, int spacedim, typename VectorType>
   void
   get_position_vector(const DoFHandler<dim, spacedim> &dh,
                       VectorType &                     vector,
                       const ComponentMask &            mask)
+  {
+    const FiniteElement<dim, spacedim> &fe = dh.get_fe();
+    get_position_vector(
+      *fe.reference_cell().template get_default_mapping<dim, spacedim>(
+        fe.degree),
+      dh,
+      vector,
+      mask);
+  }
+
+
+
+  template <int dim, int spacedim, typename VectorType>
+  void
+  get_position_vector(const Mapping<dim, spacedim> &   map_q,
+                      const DoFHandler<dim, spacedim> &dh,
+                      VectorType &                     vector,
+                      const ComponentMask &            mask)
   {
     AssertDimension(vector.size(), dh.n_dofs());
     const FiniteElement<dim, spacedim> &fe = dh.get_fe();
@@ -652,13 +671,7 @@ namespace VectorTools
       {
         const Quadrature<dim> quad(fe.get_unit_support_points());
 
-        const auto map_q =
-          fe.reference_cell().template get_default_mapping<dim, spacedim>(
-            fe.degree);
-        FEValues<dim, spacedim>              fe_v(*map_q,
-                                     fe,
-                                     quad,
-                                     update_quadrature_points);
+        FEValues<dim, spacedim> fe_v(map_q, fe, quad, update_quadrature_points);
         std::vector<types::global_dof_index> dofs(fe.n_dofs_per_cell());
 
         AssertDimension(fe.n_dofs_per_cell(),
@@ -719,7 +732,7 @@ namespace VectorTools
         dhq.distribute_dofs(feq);
         Vector<double>      eulerq(dhq.n_dofs());
         const ComponentMask maskq(spacedim, true);
-        get_position_vector(dhq, eulerq);
+        get_position_vector(map_q, dhq, eulerq);
 
         FullMatrix<double>             transfer(fe.n_dofs_per_cell(),
                                     feq.n_dofs_per_cell());
diff --git a/tests/vector_tools/get_position_vector_01.cc b/tests/vector_tools/get_position_vector_01.cc
new file mode 100644 (file)
index 0000000..55a02ef
--- /dev/null
@@ -0,0 +1,124 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Test VectorTools::get_position_vector() with Mapping as argument.
+
+#include <deal.II/fe/fe_nothing.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping_fe_field.h>
+#include <deal.II/fe/mapping_q_generic.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+template <int dim>
+void
+test()
+{
+  // surface mesh
+  const unsigned int spacedim      = dim + 1;
+  const unsigned int fe_degree     = 3;
+  const unsigned int n_refinements = 1;
+
+  Triangulation<dim, spacedim> tria;
+  GridGenerator::hyper_sphere(tria, Point<spacedim>(), 0.5);
+  tria.refine_global(n_refinements);
+
+  FE_Q<dim, spacedim>       fe(fe_degree);
+  DoFHandler<dim, spacedim> dof_handler(tria);
+  dof_handler.distribute_dofs(fe);
+
+  FESystem<dim, spacedim>   fe_dim(fe, spacedim);
+  DoFHandler<dim, spacedim> dof_handler_dim(tria);
+  dof_handler_dim.distribute_dofs(fe_dim);
+
+  // set up base high-order mapping
+  Vector<double> euler_vector_base(dof_handler_dim.n_dofs());
+  VectorTools::get_position_vector(MappingQGeneric<dim, spacedim>(4),
+                                   dof_handler_dim,
+                                   euler_vector_base);
+  MappingFEField<dim, spacedim> mapping_base(dof_handler_dim,
+                                             euler_vector_base);
+
+  // clear manifold
+  tria.reset_all_manifolds();
+
+  // output mesh with with MappingQGeneric(degree=4)
+  {
+    DataOutBase::VtkFlags flags;
+
+    DataOut<dim, DoFHandler<dim, spacedim>> data_out;
+    data_out.set_flags(flags);
+    data_out.attach_dof_handler(dof_handler);
+
+    data_out.build_patches(
+      MappingQGeneric<dim, spacedim>(4),
+      fe_degree + 1,
+      DataOut<dim,
+              DoFHandler<dim, spacedim>>::CurvedCellRegion::curved_inner_cells);
+
+#if false
+    std::ofstream output("test.0.vtk");
+    data_out.write_vtk(output);
+#else
+    data_out.write_vtk(deallog.get_file_stream());
+#endif
+  }
+
+  // output mesh with with MappingFEField based on a vector constructed
+  // with the triangulation without manifolds and the high-order base mapping
+  {
+    Vector<double> euler_vector(dof_handler_dim.n_dofs());
+    VectorTools::get_position_vector(mapping_base,
+                                     dof_handler_dim,
+                                     euler_vector);
+    MappingFEField<dim, spacedim> mapping(dof_handler_dim, euler_vector);
+    DataOutBase::VtkFlags         flags;
+
+    DataOut<dim, DoFHandler<dim, spacedim>> data_out;
+    data_out.set_flags(flags);
+    data_out.attach_dof_handler(dof_handler);
+
+    data_out.build_patches(
+      mapping,
+      fe_degree + 1,
+      DataOut<dim,
+              DoFHandler<dim, spacedim>>::CurvedCellRegion::curved_inner_cells);
+
+#if false
+    std::ofstream output("test.1.vtk");
+    data_out.write_vtk(output);
+#else
+    data_out.write_vtk(deallog.get_file_stream());
+#endif
+  }
+}
+
+
+
+int
+main(int argc, char **argv)
+{
+  initlog();
+
+  test<1>();
+}
diff --git a/tests/vector_tools/get_position_vector_01.output b/tests/vector_tools/get_position_vector_01.output
new file mode 100644 (file)
index 0000000..8ae4910
--- /dev/null
@@ -0,0 +1,169 @@
+
+# vtk DataFile Version 3.0
+#This file was generated 
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 40 double
+0.353553 -0.353553 0
+0.265165 -0.390165 0
+0.176777 -0.426777 0
+0.0883883 -0.463388 0
+5.55112e-17 -0.500000 0
+5.55112e-17 -0.500000 0
+-0.0883883 -0.463388 0
+-0.176777 -0.426777 0
+-0.265165 -0.390165 0
+-0.353553 -0.353553 0
+-0.353553 -0.353553 0
+-0.390165 -0.265165 0
+-0.426777 -0.176777 0
+-0.463388 -0.0883883 0
+-0.500000 -5.55112e-17 0
+-0.500000 -5.55112e-17 0
+-0.463388 0.0883883 0
+-0.426777 0.176777 0
+-0.390165 0.265165 0
+-0.353553 0.353553 0
+0.353553 0.353553 0
+0.390165 0.265165 0
+0.426777 0.176777 0
+0.463388 0.0883883 0
+0.500000 5.55112e-17 0
+0.500000 5.55112e-17 0
+0.463388 -0.0883883 0
+0.426777 -0.176777 0
+0.390165 -0.265165 0
+0.353553 -0.353553 0
+-0.353553 0.353553 0
+-0.265165 0.390165 0
+-0.176777 0.426777 0
+-0.0883883 0.463388 0
+-5.55112e-17 0.500000 0
+-5.55112e-17 0.500000 0
+0.0883883 0.463388 0
+0.176777 0.426777 0
+0.265165 0.390165 0
+0.353553 0.353553 0
+
+CELLS 32 96
+2      0       1
+2      1       2
+2      2       3
+2      3       4
+2      5       6
+2      6       7
+2      7       8
+2      8       9
+2      10      11
+2      11      12
+2      12      13
+2      13      14
+2      15      16
+2      16      17
+2      17      18
+2      18      19
+2      20      21
+2      21      22
+2      22      23
+2      23      24
+2      25      26
+2      26      27
+2      27      28
+2      28      29
+2      30      31
+2      31      32
+2      32      33
+2      33      34
+2      35      36
+2      36      37
+2      37      38
+2      38      39
+
+CELL_TYPES 32
+ 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
+POINT_DATA 40
+# vtk DataFile Version 3.0
+#This file was generated 
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 40 double
+0.353553 -0.353553 0
+0.277790 -0.415753 0
+0.191304 -0.461849 0
+0.0975546 -0.490409 0
+5.55112e-17 -0.500000 0
+5.55112e-17 -0.500000 0
+-0.0975546 -0.490409 0
+-0.191304 -0.461849 0
+-0.277790 -0.415753 0
+-0.353553 -0.353553 0
+-0.353553 -0.353553 0
+-0.415753 -0.277790 0
+-0.461849 -0.191304 0
+-0.490409 -0.0975546 0
+-0.500000 -5.55112e-17 0
+-0.500000 -5.55112e-17 0
+-0.490409 0.0975546 0
+-0.461849 0.191304 0
+-0.415753 0.277790 0
+-0.353553 0.353553 0
+0.353553 0.353553 0
+0.415753 0.277790 0
+0.461849 0.191304 0
+0.490409 0.0975546 0
+0.500000 5.55112e-17 0
+0.500000 5.55112e-17 0
+0.490409 -0.0975546 0
+0.461849 -0.191304 0
+0.415753 -0.277790 0
+0.353553 -0.353553 0
+-0.353553 0.353553 0
+-0.277790 0.415753 0
+-0.191304 0.461849 0
+-0.0975546 0.490409 0
+-5.55112e-17 0.500000 0
+-5.55112e-17 0.500000 0
+0.0975546 0.490409 0
+0.191304 0.461849 0
+0.277790 0.415753 0
+0.353553 0.353553 0
+
+CELLS 32 96
+2      0       1
+2      1       2
+2      2       3
+2      3       4
+2      5       6
+2      6       7
+2      7       8
+2      8       9
+2      10      11
+2      11      12
+2      12      13
+2      13      14
+2      15      16
+2      16      17
+2      17      18
+2      18      19
+2      20      21
+2      21      22
+2      22      23
+2      23      24
+2      25      26
+2      26      27
+2      27      28
+2      28      29
+2      30      31
+2      31      32
+2      32      33
+2      33      34
+2      35      36
+2      36      37
+2      37      38
+2      38      39
+
+CELL_TYPES 32
+ 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
+POINT_DATA 40

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.