From: Timo Heister Date: Fri, 16 Sep 2016 21:16:30 +0000 (-0400) Subject: add DofTools::write_gnuplot_dof_support_point_info X-Git-Tag: v8.5.0-rc1~648^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b66475714f45754946112583687f6f07b2edc198;p=dealii.git add DofTools::write_gnuplot_dof_support_point_info - add DofTools::write_gnuplot_dof_support_point_info - add tests --- diff --git a/doc/doxygen/images/support_point_dofs1.png b/doc/doxygen/images/support_point_dofs1.png new file mode 100644 index 0000000000..e96acf45f7 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 index 0000000000..75697ea53c Binary files /dev/null and b/doc/doxygen/images/support_point_dofs2.png differ diff --git a/doc/news/changes.h b/doc/news/changes.h index 4270543656..31e2a836ca 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -219,19 +219,27 @@ inconvenience this causes.
    +
  1. New: DoFTools::write_gnuplot_dof_support_point_info outputs + support point locations and dof indices to a format readable by + gnuplot. +
    + (Timo Heister, 2016/09/16) +
  2. +
  3. Fixed: The CMake macros DEAL_II_(ADD_TEST|SETUP_TARGET) now enforce a stricter CMAKE_BUILD_TYPE handling. This helps to avoid situations where targets with different build flavors might accidentally get linked against each other.
    (Matthias Maier, 2016/09/08) +
  4. Improved: the doxygen documentation now contains nicely formatted boxes containing the text message of each exception. Several messages haven been clarified and improved.
    (Timo Heister, 2016/09/06) -
  5. +
  6. Fixed: Reimplement copy_triangulation and load in dealii::parallel::shared::Triangulation, this avoids the loss of diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index ef07a09839..1816cca7a2 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -29,6 +29,8 @@ #include #include #include +#include + DEAL_II_NAMESPACE_OPEN @@ -51,7 +53,6 @@ namespace GridTools template 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 void @@ -2196,6 +2202,37 @@ namespace DoFTools const DoFHandlerType &dof_handler, std::map, 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): + *

    + * @image html support_point_dofs1.png + * @image html support_point_dofs2.png + *

    + */ + template + void + write_gnuplot_dof_support_point_info(std::ostream &out, + const std::map > &support_points); + + /** * Map a coupling table from the user friendly organization by components to * the organization by blocks. Specializations of this function for diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index aadfb77635..bc666ba5c9 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -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 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 + struct ComparisonHelper + { + /** + * Comparison operator. + * + * Return true if @lhs is considered less than @rhs. + */ + bool operator() (const Point &lhs, + const Point &rhs) const + { + for (unsigned int d=0; d + void + write_gnuplot_dof_support_point_info(std::ostream &out, + const std::map > &support_points) + { + AssertThrow (out, ExcIO()); + + typedef std::map< types::global_dof_index, Point > + dof_map_t; + + typedef std::map, std::vector, typename internal::ComparisonHelper > + 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 &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 &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 void diff --git a/source/dofs/dof_tools.inst.in b/source/dofs/dof_tools.inst.in index 8495a6ff4d..a11b5b3077 100644 --- a/source/dofs/dof_tools.inst.in +++ b/source/dofs/dof_tools.inst.in @@ -688,6 +688,12 @@ for (deal_II_dimension : DIMENSIONS) #endif + template + void + DoFTools::write_gnuplot_dof_support_point_info(std::ostream &, + const std::map > &); + + 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 index 0000000000..ebcc47f785 --- /dev/null +++ b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.cc @@ -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 +#include + +#include +#include +#include +#include + +#include +#include +#include + +template +void +test () +{ + std::ostream &out = deallog.get_file_stream(); + + Triangulation triangulation; + const int velocity_degree = 2; + + std::vector rep(dim, 1); + rep[0]=2; + GridGenerator::subdivided_hyper_rectangle(triangulation, + rep, + Point(), + Point(2.0,1.0)); + FESystem fe (FE_Q(velocity_degree), dim, + FE_Q(velocity_degree-1), 1); + + DoFHandler 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 > support_points; + DoFTools::map_dofs_to_support_points (MappingQ1(), + 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 index 0000000000..fadd5dff4e --- /dev/null +++ b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_01.output @@ -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 index 0000000000..623259d923 --- /dev/null +++ b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.cc @@ -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 +#include + +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +template +void +test () +{ + std::ostream &out = deallog.get_file_stream(); + + Point center; + + Triangulation triangulation; + GridGenerator::quarter_hyper_shell(triangulation, center, 0.5, 1.0); + + static const PolarManifold manifold(center); + triangulation.set_manifold (0, manifold); + triangulation.set_all_manifold_ids_on_boundary(0); + triangulation.set_all_manifold_ids(0); + + FE_Q fe (2); + + DoFHandler dof_handler (triangulation); + dof_handler.distribute_dofs (fe); + + MappingManifold mapping; + + std::map > 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 index 0000000000..f51627f212 --- /dev/null +++ b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_02.output @@ -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 index 0000000000..d80e908d9a --- /dev/null +++ b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.cc @@ -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 +#include + +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include + +template +void +test () +{ + std::ostream &out = deallog.get_file_stream(); + + Point center; + + Triangulation triangulation; + GridGenerator::quarter_hyper_shell(triangulation, center, 0.5, 1.0); + + static const PolarManifold manifold(center); + triangulation.set_manifold (0, manifold); + triangulation.set_all_manifold_ids_on_boundary(0); + triangulation.set_all_manifold_ids(0); + + FE_Q fe (2); + + DoFHandler dof_handler (triangulation); + dof_handler.distribute_dofs (fe); + + MappingQ mapping(3); + + std::map > 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 index 0000000000..e27f1e9d7e --- /dev/null +++ b/tests/dofs/dof_tools_write_gnuplot_dof_support_point_info_03.output @@ -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