]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add DofTools::write_gnuplot_dof_support_point_info 3126/head
authorTimo Heister <timo.heister@gmail.com>
Fri, 16 Sep 2016 21:16:30 +0000 (17:16 -0400)
committerTimo Heister <timo.heister@gmail.com>
Tue, 20 Sep 2016 19:27:54 +0000 (15:27 -0400)
- add DofTools::write_gnuplot_dof_support_point_info
- add tests

12 files changed:
doc/doxygen/images/support_point_dofs1.png [new file with mode: 0644]
doc/doxygen/images/support_point_dofs2.png [new file with mode: 0644]
doc/news/changes.h
include/deal.II/dofs/dof_tools.h
source/dofs/dof_tools.cc
source/dofs/dof_tools.inst.in
tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.cc [new file with mode: 0644]
tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.output [new file with mode: 0644]
tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.cc [new file with mode: 0644]
tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.output [new file with mode: 0644]
tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.cc [new file with mode: 0644]
tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.output [new file with mode: 0644]

diff --git a/doc/doxygen/images/support_point_dofs1.png b/doc/doxygen/images/support_point_dofs1.png
new file mode 100644 (file)
index 0000000..e96acf4
Binary files /dev/null and b/doc/doxygen/images/support_point_dofs1.png differ
diff --git a/doc/doxygen/images/support_point_dofs2.png b/doc/doxygen/images/support_point_dofs2.png
new file mode 100644 (file)
index 0000000..75697ea
Binary files /dev/null and b/doc/doxygen/images/support_point_dofs2.png differ
index 4270543656555a8faba620f184751719bd7059d0..31e2a836ca95d27aecf1b2d963b95a898a9d310c 100644 (file)
@@ -219,19 +219,27 @@ inconvenience this causes.
 
 <ol>
 
+ <li> New: DoFTools::write_gnuplot_dof_support_point_info outputs
+ support point locations and dof indices to a format readable by
+ gnuplot.
+ <br>
+ (Timo Heister, 2016/09/16)
+ </li>
+
  <li> Fixed: The CMake macros <code>DEAL_II_(ADD_TEST|SETUP_TARGET)</code>
  now enforce a stricter <code>CMAKE_BUILD_TYPE</code> handling. This helps
  to avoid situations where targets with different build flavors might
  accidentally get linked against each other.
  <br>
  (Matthias Maier, 2016/09/08)
+ </li>
 
  <li> Improved: the doxygen documentation now contains nicely formatted
  boxes containing the text message of each exception. Several messages
  haven been clarified and improved.
  <br>
  (Timo Heister, 2016/09/06)
-</li>
+ </li>
 
  <li> Fixed: Reimplement copy_triangulation and load in
  dealii::parallel::shared::Triangulation, this avoids the loss of
index ef07a0983958cdb63040a0020563f23c7fc5f167..1816cca7a2b562c169ed0d69acfd754d8f24bf48 100644 (file)
@@ -29,6 +29,8 @@
 #include <vector>
 #include <set>
 #include <map>
+#include <ostream>
+
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -51,7 +53,6 @@ namespace GridTools
   template <typename CellIterator> struct PeriodicFacePair;
 }
 
-//TODO: map_support_points_to_dofs should generate a multimap, rather than just a map, since several dofs may be located at the same support point
 
 /**
  * This is a collection of functions operating on, and manipulating the
@@ -2188,6 +2189,11 @@ namespace DoFTools
    * Just as with the function above, it is assumed that the finite element in
    * use here actually supports the notion of support points of all its
    * components.
+   *
+   * @todo This function should generate a multimap, rather than just a map,
+   * since several dofs may be located at the same support point. Currently,
+   * only the last value in the map returned by map_dofs_to_support_points() for
+   * each point will be returned.
    */
   template <typename DoFHandlerType, class Comp>
   void
@@ -2196,6 +2202,37 @@ namespace DoFTools
    const DoFHandlerType                                                            &dof_handler,
    std::map<Point<DoFHandlerType::space_dimension>, types::global_dof_index, Comp> &point_to_index_map);
 
+  /**
+   * Generate text output readable by gnuplot with point data based on the
+   * given map @p support_points.  For each support point location, a string
+   * label containing a list of all DoFs from the map is generated.  The map
+   * can be generated with a call to map_dofs_to_support_points() and is useful
+   * to visualize location and global numbering of unknowns.
+   *
+   * An example for the format of each line in the output is:
+   * @code
+   * x [y] [z] "dof1, dof2"
+   * @endcode
+   * where x, y, and z (present only in corresponding dimension) are the
+   * coordinates of the support point, followed by a list of DoF numbers.
+   *
+   * The points with labels can be plotted as follows in gnuplot:
+   * @code
+   * plot "./points.gpl" using 1:2:3 with labels point offset 1,1
+   * @endcode
+   *
+   * Examples (this also includes the grid written separately using GridOut):
+   * <p ALIGN="center">
+   * @image html support_point_dofs1.png
+   * @image html support_point_dofs2.png
+   * </p>
+   */
+  template <int spacedim>
+  void
+  write_gnuplot_dof_support_point_info(std::ostream &out,
+                                       const std::map<types::global_dof_index, Point<spacedim> > &support_points);
+
+
   /**
    * Map a coupling table from the user friendly organization by components to
    * the organization by blocks. Specializations of this function for
index aadfb77635499638d8ff502ac0fa6d834c7894dd..bc666ba5c969389df14cb03a20214440fd34c2dc 100644 (file)
@@ -54,6 +54,40 @@ namespace DoFTools
 {
   namespace internal
   {
+    /**
+     * Comparison functor struct to compare two Points and return if one is
+     * "less" than the other one. This can be used to use Point<dim> as a key in
+     * std::map.
+     *
+     * Comparison is done by comparing values in each dimension in ascending
+     * order (first x, then y, etc.). Note that comparisons are done without an
+     * epsilon, so points need to have identical floating point components to be
+     * considered equal.
+     */
+    template <int dim,typename Number=double>
+    struct ComparisonHelper
+    {
+      /**
+        * Comparison operator.
+        *
+        * Return true if @lhs is considered less than @rhs.
+        */
+      bool operator() (const Point<dim, Number> &lhs,
+                       const Point<dim, Number> &rhs) const
+      {
+        for (unsigned int d=0; d<dim; ++d)
+          {
+            if (lhs[d] == rhs[d])
+              continue;
+            return lhs[d]<rhs[d];
+          }
+        return false;
+      }
+
+    };
+
+
+
     // return an array that for each dof on the reference cell
     // lists the corresponding vector component.
     //
@@ -1998,6 +2032,50 @@ namespace DoFTools
                                           support_points);
   }
 
+  template <int spacedim>
+  void
+  write_gnuplot_dof_support_point_info(std::ostream &out,
+                                       const std::map<types::global_dof_index, Point<spacedim> > &support_points)
+  {
+    AssertThrow (out, ExcIO());
+
+    typedef std::map< types::global_dof_index, Point<spacedim> >
+    dof_map_t;
+
+    typedef std::map<Point<spacedim>, std::vector<types::global_dof_index>, typename internal::ComparisonHelper<spacedim> >
+    point_map_t;
+
+    point_map_t point_map;
+
+    // convert to map point -> list of DoFs
+    for (typename dof_map_t::const_iterator it = support_points.begin();
+         it!=support_points.end();
+         ++it)
+      {
+        std::vector<types::global_dof_index> &v = point_map[it->second];
+        v.push_back(it->first);
+      }
+
+    // print the newly created map:
+    for (typename point_map_t::iterator it = point_map.begin();
+         it!=point_map.end();
+         ++it)
+      {
+        out << it->first << " \"";
+        const std::vector<types::global_dof_index> &v = it->second;
+        for (unsigned int i=0; i < v.size(); ++i)
+          {
+            if (i>0)
+              out << ", ";
+            out << v[i];
+          }
+
+        out << "\"\n";
+      }
+
+    out << std::flush;
+  }
+
 
   template<int dim, int spacedim>
   void
index 8495a6ff4d192bad27a752b970e39e51b8f8aa3e..a11b5b3077b87e0b6b8a20563ffcf9a56176d5bd 100644 (file)
@@ -688,6 +688,12 @@ for (deal_II_dimension : DIMENSIONS)
 
 #endif
 
+    template
+    void
+    DoFTools::write_gnuplot_dof_support_point_info(std::ostream &,
+            const std::map<types::global_dof_index, Point<deal_II_dimension> > &);
+
+
     template
     void
     DoFTools::convert_couplings_to_blocks (
diff --git a/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.cc b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.cc
new file mode 100644 (file)
index 0000000..ebcc47f
--- /dev/null
@@ -0,0 +1,74 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 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 DoFTools::write_gnuplot_dof_support_point_info with two
+// cells and a Stokes Taylor-Hood element.
+
+#include "../tests.h"
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.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_q1.h>
+
+template <int dim>
+void
+test ()
+{
+  std::ostream &out = deallog.get_file_stream();
+
+  Triangulation<dim> triangulation;
+  const int velocity_degree = 2;
+
+  std::vector<unsigned int> rep(dim, 1);
+  rep[0]=2;
+  GridGenerator::subdivided_hyper_rectangle(triangulation,
+                                            rep,
+                                            Point<dim>(),
+                                            Point<dim>(2.0,1.0));
+  FESystem<dim> fe (FE_Q<dim>(velocity_degree), dim,
+                    FE_Q<dim>(velocity_degree-1), 1);
+
+  DoFHandler<dim> dof_handler (triangulation);
+  dof_handler.distribute_dofs (fe);
+
+  out << "plot '-' using 1:2 with lines, '-' with labels point pt 2 offset 1,1" << std::endl;
+  GridOut().write_gnuplot (triangulation, deallog.get_file_stream());
+  out << "e" << std::endl;
+
+  std::map<types::global_dof_index, Point<dim> > support_points;
+  DoFTools::map_dofs_to_support_points (MappingQ1<dim>(),
+                                        dof_handler,
+                                        support_points);
+  DoFTools::write_gnuplot_dof_support_point_info(deallog.get_file_stream(),
+                                                 support_points);
+  out << "e" << std::endl;
+}
+
+
+
+int main()
+{
+  initlog();
+  test<2>();
+}
+
diff --git a/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.output b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.output
new file mode 100644 (file)
index 0000000..fadd5df
--- /dev/null
@@ -0,0 +1,33 @@
+
+plot '-' using 1:2 with lines, '-' with labels point pt 2 offset 1,1
+0.00000 0.00000 0 0
+1.00000 0.00000 0 0
+1.00000 1.00000 0 0
+0.00000 1.00000 0 0
+0.00000 0.00000 0 0
+
+
+1.00000 0.00000 0 0
+2.00000 0.00000 0 0
+2.00000 1.00000 0 0
+1.00000 1.00000 0 0
+1.00000 0.00000 0 0
+
+
+e
+0.00000 0.00000 "0, 1, 2"
+0.00000 0.500000 "12, 13"
+0.00000 1.00000 "6, 7, 8"
+0.500000 0.00000 "16, 17"
+0.500000 0.500000 "20, 21"
+0.500000 1.00000 "18, 19"
+1.00000 0.00000 "3, 4, 5"
+1.00000 0.500000 "14, 15"
+1.00000 1.00000 "9, 10, 11"
+1.50000 0.00000 "30, 31"
+1.50000 0.500000 "34, 35"
+1.50000 1.00000 "32, 33"
+2.00000 0.00000 "22, 23, 24"
+2.00000 0.500000 "28, 29"
+2.00000 1.00000 "25, 26, 27"
+e
diff --git a/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.cc b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.cc
new file mode 100644 (file)
index 0000000..623259d
--- /dev/null
@@ -0,0 +1,79 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 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 DoFTools::write_gnuplot_dof_support_point_info with a mapping
+
+#include "../tests.h"
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/manifold_lib.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_q1.h>
+#include <deal.II/fe/mapping_manifold.h>
+
+template <int dim>
+void
+test ()
+{
+  std::ostream &out = deallog.get_file_stream();
+
+  Point<dim> center;
+
+  Triangulation<dim> triangulation;
+  GridGenerator::quarter_hyper_shell(triangulation, center, 0.5, 1.0);
+
+  static const PolarManifold<dim,dim> manifold(center);
+  triangulation.set_manifold (0, manifold);
+  triangulation.set_all_manifold_ids_on_boundary(0);
+  triangulation.set_all_manifold_ids(0);
+
+  FE_Q<dim> fe (2);
+
+  DoFHandler<dim> dof_handler (triangulation);
+  dof_handler.distribute_dofs (fe);
+
+  MappingManifold<dim> mapping;
+
+  std::map<types::global_dof_index, Point<dim> > support_points;
+  DoFTools::map_dofs_to_support_points (mapping,
+                                        dof_handler,
+                                        support_points);
+
+  out << "set view equal xy" << std::endl
+      << "plot '-' using 1:2 with lines, '-' with labels point pt 2 offset 1,1" << std::endl;
+  GridOut().write_gnuplot (triangulation, deallog.get_file_stream());
+  out << "e" << std::endl;
+
+  DoFTools::write_gnuplot_dof_support_point_info(deallog.get_file_stream(),
+                                                 support_points);
+  out << "e" << std::endl;
+}
+
+
+
+int main()
+{
+  initlog();
+  test<2>();
+}
+
diff --git a/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.output b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.output
new file mode 100644 (file)
index 0000000..f51627f
--- /dev/null
@@ -0,0 +1,47 @@
+
+set view equal xy
+plot '-' using 1:2 with lines, '-' with labels point pt 2 offset 1,1
+1.00000 0.00000 0 0
+0.866025 0.500000 0 0
+0.433013 0.250000 0 0
+0.500000 0.00000 0 0
+1.00000 0.00000 0 0
+
+
+0.866025 0.500000 0 0
+0.500000 0.866025 0 0
+0.250000 0.433013 0 0
+0.433013 0.250000 0 0
+0.866025 0.500000 0 0
+
+
+0.500000 0.866025 0 0
+0.00000 1.00000 0 0
+0.00000 0.500000 0 0
+0.250000 0.433013 0 0
+0.500000 0.866025 0 0
+
+
+e
+3.06162e-17 0.500000 "16"
+4.59243e-17 0.750000 "17"
+6.12323e-17 1.00000 "15"
+0.129410 0.482963 "19"
+0.194114 0.724444 "20"
+0.250000 0.433013 "10"
+0.258819 0.965926 "18"
+0.353553 0.353553 "13"
+0.375000 0.649519 "11"
+0.433013 0.250000 "3"
+0.482963 0.129410 "7"
+0.500000 0.00000 "2"
+0.500000 0.866025 "9"
+0.530330 0.530330 "14"
+0.649519 0.375000 "5"
+0.707107 0.707107 "12"
+0.724444 0.194114 "8"
+0.750000 0.00000 "4"
+0.866025 0.500000 "1"
+0.965926 0.258819 "6"
+1.00000 0.00000 "0"
+e
diff --git a/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.cc b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.cc
new file mode 100644 (file)
index 0000000..d80e908
--- /dev/null
@@ -0,0 +1,80 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 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 DoFTools::write_gnuplot_dof_support_point_info with a mapping
+
+#include "../tests.h"
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/manifold_lib.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_q1.h>
+#include <deal.II/fe/mapping_q.h>
+#include <deal.II/fe/mapping_manifold.h>
+
+template <int dim>
+void
+test ()
+{
+  std::ostream &out = deallog.get_file_stream();
+
+  Point<dim> center;
+
+  Triangulation<dim> triangulation;
+  GridGenerator::quarter_hyper_shell(triangulation, center, 0.5, 1.0);
+
+  static const PolarManifold<dim,dim> manifold(center);
+  triangulation.set_manifold (0, manifold);
+  triangulation.set_all_manifold_ids_on_boundary(0);
+  triangulation.set_all_manifold_ids(0);
+
+  FE_Q<dim> fe (2);
+
+  DoFHandler<dim> dof_handler (triangulation);
+  dof_handler.distribute_dofs (fe);
+
+  MappingQ<dim> mapping(3);
+
+  std::map<types::global_dof_index, Point<dim> > support_points;
+  DoFTools::map_dofs_to_support_points (mapping,
+                                        dof_handler,
+                                        support_points);
+
+  out << "set view equal xy" << std::endl
+      << "plot '-' using 1:2 with lines, '-' with labels point pt 2 offset 1,1" << std::endl;
+  GridOut().write_gnuplot (triangulation, deallog.get_file_stream());
+  out << "e" << std::endl;
+
+  DoFTools::write_gnuplot_dof_support_point_info(deallog.get_file_stream(),
+                                                 support_points);
+  out << "e" << std::endl;
+}
+
+
+
+int main()
+{
+  initlog();
+  test<2>();
+}
+
diff --git a/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.output b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.output
new file mode 100644 (file)
index 0000000..e27f1e9
--- /dev/null
@@ -0,0 +1,47 @@
+
+set view equal xy
+plot '-' using 1:2 with lines, '-' with labels point pt 2 offset 1,1
+1.00000 0.00000 0 0
+0.866025 0.500000 0 0
+0.433013 0.250000 0 0
+0.500000 0.00000 0 0
+1.00000 0.00000 0 0
+
+
+0.866025 0.500000 0 0
+0.500000 0.866025 0 0
+0.250000 0.433013 0 0
+0.433013 0.250000 0 0
+0.866025 0.500000 0 0
+
+
+0.500000 0.866025 0 0
+0.00000 1.00000 0 0
+0.00000 0.500000 0 0
+0.250000 0.433013 0 0
+0.500000 0.866025 0 0
+
+
+e
+0.00000 0.500000 "16"
+0.00000 1.00000 "15"
+5.74034e-17 0.750000 "17"
+0.129404 0.482944 "19"
+0.194107 0.724416 "20"
+0.250000 0.433013 "10"
+0.258809 0.965888 "18"
+0.353540 0.353540 "13"
+0.375000 0.649519 "11"
+0.433013 0.250000 "3"
+0.482944 0.129404 "7"
+0.500000 0.00000 "2"
+0.500000 0.866025 "9"
+0.530309 0.530309 "14"
+0.649519 0.375000 "5"
+0.707079 0.707079 "12"
+0.724416 0.194107 "8"
+0.750000 0.00000 "4"
+0.866025 0.500000 "1"
+0.965888 0.258809 "6"
+1.00000 0.00000 "0"
+e

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.