Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
data_out_base.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_data_out_base_h
17 #define dealii_data_out_base_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 #include <deal.II/base/mpi_stub.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/base/table.h>
26 
28 
30 
31 // To be able to serialize XDMFEntry
32 #include <boost/serialization/map.hpp>
33 
34 #include <limits>
35 #include <ostream>
36 #include <string>
37 #include <tuple>
38 #include <typeinfo>
39 #include <vector>
40 
42 
43 // Forward declarations
44 #ifndef DOXYGEN
45 class ParameterHandler;
46 class XDMFEntry;
47 #endif
48 
198 namespace DataOutBase
199 {
207  enum class CompressionLevel
208  {
216  best_speed,
229  plain_text
230  };
231 
232 
255  template <int dim, int spacedim = dim>
256  struct Patch
257  {
261  static const unsigned int space_dim = spacedim;
262 
280  std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell> vertices;
281 
292  std::array<unsigned int, GeometryInfo<dim>::faces_per_cell> neighbors;
293 
298  unsigned int patch_index;
299 
305  unsigned int n_subdivisions;
306 
329 
343 
348 
353  Patch();
354 
359  bool
360  operator==(const Patch &patch) const;
361 
368  std::size_t
369  memory_consumption() const;
370 
374  void
375  swap(Patch<dim, spacedim> &other_patch) noexcept;
376 
380  static const unsigned int no_neighbor = numbers::invalid_unsigned_int;
381 
392  int,
393  int,
394  << "It is not possible to have a structural dimension of " << arg1
395  << " to be larger than the space dimension of the surrounding"
396  << " space " << arg2);
398  };
399 
400 
401 
417  template <int spacedim>
418  struct Patch<0, spacedim>
419  {
423  static const unsigned int space_dim = spacedim;
424 
434 
439  static unsigned int neighbors[1];
440 
445  unsigned int patch_index;
446 
456  static const unsigned int n_subdivisions;
457 
479 
493 
500 
505  Patch();
506 
511  bool
512  operator==(const Patch &patch) const;
513 
520  std::size_t
521  memory_consumption() const;
522 
526  void
527  swap(Patch<0, spacedim> &other_patch) noexcept;
528 
532  static const unsigned int no_neighbor = numbers::invalid_unsigned_int;
533 
544  int,
545  int,
546  << "It is not possible to have a structural dimension of " << arg1
547  << " to be larger than the space dimension of the surrounding"
548  << " space " << arg2);
550  };
551 
552 
564  template <typename FlagsType>
566  {
574  static void
576 
584  void
586 
593  std::size_t
595  };
596 
597 
598  template <typename FlagsType>
599  void
601  {}
602 
603 
604  template <typename FlagsType>
605  void
607  {}
608 
609 
610  template <typename FlagsType>
611  std::size_t
613  {
614  return sizeof(FlagsType);
615  }
616 
617 
623  struct DXFlags : public OutputFlagsBase<DXFlags>
624  {
639 
644 
650 
654  DXFlags(const bool write_neighbors = false,
655  const bool int_binary = false,
656  const bool coordinates_binary = false,
657  const bool data_binary = false);
658 
663  static void
665 
672  void
674  };
675 
681  struct UcdFlags : public OutputFlagsBase<UcdFlags>
682  {
693 
697  UcdFlags(const bool write_preamble = false);
698 
703  static void
705 
712  void
714  };
715 
721  struct GnuplotFlags : public OutputFlagsBase<GnuplotFlags>
722  {
727  GnuplotFlags();
728 
732  GnuplotFlags(const std::vector<std::string> &space_dimension_labels);
733 
749  std::vector<std::string> space_dimension_labels;
750 
755  std::size_t
756  memory_consumption() const;
757 
763  "There should be at least one space dimension per spatial "
764  "dimension (extras are ignored).");
765  };
766 
773  struct PovrayFlags : public OutputFlagsBase<PovrayFlags>
774  {
780  bool smooth;
781 
788 
796 
800  PovrayFlags(const bool smooth = false,
801  const bool bicubic_patch = false,
802  const bool external_data = false);
803 
808  static void
810 
817  void
819  };
820 
821 
828  struct EpsFlags : public OutputFlagsBase<EpsFlags>
829  {
836  unsigned int height_vector;
837 
842  unsigned int color_vector;
843 
849  enum SizeType
850  {
854  height
855  };
856 
861 
869  unsigned int size;
870 
874  double line_width;
875 
881  double azimut_angle;
882 
906  double turn_angle;
907 
918  double z_scaling;
919 
926  bool draw_mesh;
927 
945 
955 
959  struct RgbValues
960  {
961  float red;
962  float green;
963  float blue;
964 
969  bool
970  is_grey() const;
971  };
972 
980  using ColorFunction = RgbValues (*)(const double value,
981  const double min_value,
982  const double max_value);
983 
990 
991 
1000  static RgbValues
1001  default_color_function(const double value,
1002  const double min_value,
1003  const double max_value);
1004 
1010  static RgbValues
1011  grey_scale_color_function(const double value,
1012  const double min_value,
1013  const double max_value);
1014 
1021  static RgbValues
1022  reverse_grey_scale_color_function(const double value,
1023  const double min_value,
1024  const double max_value);
1025 
1029  EpsFlags(const unsigned int height_vector = 0,
1030  const unsigned int color_vector = 0,
1031  const SizeType size_type = width,
1032  const unsigned int size = 300,
1033  const double line_width = 0.5,
1034  const double azimut_angle = 60,
1035  const double turn_angle = 30,
1036  const double z_scaling = 1.0,
1037  const bool draw_mesh = true,
1038  const bool draw_cells = true,
1039  const bool shade_cells = true,
1041 
1049  static void
1051 
1058  void
1059  parse_parameters(const ParameterHandler &prm);
1060  };
1061 
1068  struct GmvFlags : public OutputFlagsBase<GmvFlags>
1069  {};
1070 
1076  struct Hdf5Flags : public OutputFlagsBase<Hdf5Flags>
1077  {
1083 
1084  explicit Hdf5Flags(
1086  };
1087 
1093  struct TecplotFlags : public OutputFlagsBase<TecplotFlags>
1094  {
1099  const char *zone_name;
1100 
1107 
1111  TecplotFlags(const char *zone_name = nullptr,
1112  const double solution_time = -1.0);
1113 
1118  std::size_t
1119  memory_consumption() const;
1120  };
1121 
1127  struct VtkFlags : public OutputFlagsBase<VtkFlags>
1128  {
1139  double time;
1140 
1151  unsigned int cycle;
1152 
1160 
1169 
1178 
1184 
1203 
1225  std::map<std::string, std::string> physical_units;
1226 
1231  explicit VtkFlags(
1232  const double time = std::numeric_limits<double>::min(),
1233  const unsigned int cycle = std::numeric_limits<unsigned int>::min(),
1234  const bool print_date_and_time = true,
1236  const bool write_higher_order_cells = false,
1237  const std::map<std::string, std::string> &physical_units = {});
1238  };
1239 
1240 
1246  struct SvgFlags : public OutputFlagsBase<SvgFlags>
1247  {
1251  unsigned int height;
1252 
1257  unsigned int width;
1258 
1265  unsigned int height_vector;
1266 
1271 
1272  unsigned int line_thickness;
1273 
1277  bool margin;
1278 
1283 
1287  SvgFlags(const unsigned int height_vector = 0,
1288  const int azimuth_angle = 37,
1289  const int polar_angle = 45,
1290  const unsigned int line_thickness = 1,
1291  const bool margin = true,
1292  const bool draw_colorbar = true);
1293  };
1294 
1295 
1303  : public OutputFlagsBase<Deal_II_IntermediateFlags>
1304  {
1311  static const unsigned int format_version;
1312  };
1313 
1320  {
1331 
1337 
1341  DataOutFilterFlags(const bool filter_duplicate_vertices = false,
1342  const bool xdmf_hdf5_output = false);
1343 
1348  static void
1350 
1357  void
1358  parse_parameters(const ParameterHandler &prm);
1359 
1364  std::size_t
1366  };
1367 
1400  {
1401  public:
1405  DataOutFilter();
1406 
1412 
1418  template <int dim>
1419  void
1420  write_point(const unsigned int index, const Point<dim> &p);
1421 
1425  template <int dim>
1426  void
1427  write_cell(const unsigned int index,
1428  const unsigned int start,
1429  const std::array<unsigned int, dim> &offsets);
1430 
1435  void
1436  write_cell_single(const unsigned int index,
1437  const unsigned int start,
1438  const unsigned int n_points,
1439  const ReferenceCell &reference_cell);
1440 
1447  void
1448  write_data_set(const std::string &name,
1449  const unsigned int dimension,
1450  const unsigned int set_num,
1451  const Table<2, double> &data_vectors);
1452 
1457  void
1458  fill_node_data(std::vector<double> &node_data) const;
1459 
1464  void
1465  fill_cell_data(const unsigned int local_node_offset,
1466  std::vector<unsigned int> &cell_data) const;
1467 
1471  std::string
1472  get_data_set_name(const unsigned int set_num) const;
1473 
1477  unsigned int
1478  get_data_set_dim(const unsigned int set_num) const;
1479 
1484  const double *
1485  get_data_set(const unsigned int set_num) const;
1486 
1491  unsigned int
1492  n_nodes() const;
1493 
1498  unsigned int
1499  n_cells() const;
1500 
1505  unsigned int
1506  n_data_sets() const;
1507 
1511  void
1512  flush_points();
1513 
1517  void
1518  flush_cells();
1519 
1520 
1521  private:
1525  struct Point3Comp
1526  {
1527  bool
1528  operator()(const Point<3> &one, const Point<3> &two) const
1529  {
1530  /*
1531  * The return statement below is an optimized version of the following
1532  * code:
1533  *
1534  * for (unsigned int d=0; d<3; ++d)
1535  * {
1536  * if (one(d) < two(d))
1537  * return true;
1538  * else if (one(d) > two(d))
1539  * return false;
1540  * }
1541  * return false;
1542  */
1543 
1544  return (one(0) < two(0) ||
1545  (!(two(0) < one(0)) &&
1546  (one(1) < two(1) || (!(two(1) < one(1)) && one(2) < two(2)))));
1547  }
1548  };
1549 
1550  using Map3DPoint = std::multimap<Point<3>, unsigned int, Point3Comp>;
1551 
1556 
1563  unsigned int node_dim;
1564 
1569  unsigned int num_cells;
1570 
1575 
1579  std::map<unsigned int, unsigned int> filtered_points;
1580 
1584  std::map<unsigned int, unsigned int> filtered_cells;
1585 
1589  std::vector<std::string> data_set_names;
1590 
1594  std::vector<unsigned int> data_set_dims;
1595 
1599  std::vector<std::vector<double>> data_sets;
1600 
1604  void
1605  internal_add_cell(const unsigned int cell_index,
1606  const unsigned int pt_index);
1607  };
1608 
1609 
1614  {
1619 
1624 
1629 
1634 
1639 
1644 
1649 
1654 
1659 
1664 
1669 
1674 
1679 
1683  hdf5
1684  };
1685 
1686 
1690  template <int dim, int spacedim>
1691  void
1692  write_dx(
1693  const std::vector<Patch<dim, spacedim>> &patches,
1694  const std::vector<std::string> &data_names,
1695  const std::vector<
1696  std::tuple<unsigned int,
1697  unsigned int,
1698  std::string,
1700  &nonscalar_data_ranges,
1701  const DXFlags &flags,
1702  std::ostream &out);
1703 
1748  template <int spacedim>
1749  void
1750  write_eps(
1751  const std::vector<Patch<2, spacedim>> &patches,
1752  const std::vector<std::string> &data_names,
1753  const std::vector<
1754  std::tuple<unsigned int,
1755  unsigned int,
1756  std::string,
1758  &nonscalar_data_ranges,
1759  const EpsFlags &flags,
1760  std::ostream &out);
1761 
1767  template <int dim, int spacedim>
1768  void
1769  write_eps(
1770  const std::vector<Patch<dim, spacedim>> &patches,
1771  const std::vector<std::string> &data_names,
1772  const std::vector<
1773  std::tuple<unsigned int,
1774  unsigned int,
1775  std::string,
1777  &nonscalar_data_ranges,
1778  const EpsFlags &flags,
1779  std::ostream &out);
1780 
1781 
1791  template <int dim, int spacedim>
1792  void
1793  write_gmv(
1794  const std::vector<Patch<dim, spacedim>> &patches,
1795  const std::vector<std::string> &data_names,
1796  const std::vector<
1797  std::tuple<unsigned int,
1798  unsigned int,
1799  std::string,
1801  &nonscalar_data_ranges,
1802  const GmvFlags &flags,
1803  std::ostream &out);
1804 
1860  template <int dim, int spacedim>
1861  void
1862  write_gnuplot(
1863  const std::vector<Patch<dim, spacedim>> &patches,
1864  const std::vector<std::string> &data_names,
1865  const std::vector<
1866  std::tuple<unsigned int,
1867  unsigned int,
1868  std::string,
1870  &nonscalar_data_ranges,
1871  const GnuplotFlags &flags,
1872  std::ostream &out);
1873 
1919  template <int dim, int spacedim>
1920  void
1921  write_povray(
1922  const std::vector<Patch<dim, spacedim>> &patches,
1923  const std::vector<std::string> &data_names,
1924  const std::vector<
1925  std::tuple<unsigned int,
1926  unsigned int,
1927  std::string,
1929  &nonscalar_data_ranges,
1930  const PovrayFlags &flags,
1931  std::ostream &out);
1932 
1939  template <int dim, int spacedim>
1940  void
1941  write_tecplot(
1942  const std::vector<Patch<dim, spacedim>> &patches,
1943  const std::vector<std::string> &data_names,
1944  const std::vector<
1945  std::tuple<unsigned int,
1946  unsigned int,
1947  std::string,
1949  &nonscalar_data_ranges,
1950  const TecplotFlags &flags,
1951  std::ostream &out);
1952 
1967  template <int dim, int spacedim>
1968  void
1969  write_ucd(
1970  const std::vector<Patch<dim, spacedim>> &patches,
1971  const std::vector<std::string> &data_names,
1972  const std::vector<
1973  std::tuple<unsigned int,
1974  unsigned int,
1975  std::string,
1977  &nonscalar_data_ranges,
1978  const UcdFlags &flags,
1979  std::ostream &out);
1980 
2000  template <int dim, int spacedim>
2001  void
2002  write_vtk(
2003  const std::vector<Patch<dim, spacedim>> &patches,
2004  const std::vector<std::string> &data_names,
2005  const std::vector<
2006  std::tuple<unsigned int,
2007  unsigned int,
2008  std::string,
2010  &nonscalar_data_ranges,
2011  const VtkFlags &flags,
2012  std::ostream &out);
2013 
2014 
2038  template <int dim, int spacedim>
2039  void
2040  write_vtu(
2041  const std::vector<Patch<dim, spacedim>> &patches,
2042  const std::vector<std::string> &data_names,
2043  const std::vector<
2044  std::tuple<unsigned int,
2045  unsigned int,
2046  std::string,
2048  &nonscalar_data_ranges,
2049  const VtkFlags &flags,
2050  std::ostream &out);
2051 
2057  void
2058  write_vtu_header(std::ostream &out, const VtkFlags &flags);
2059 
2066  void
2067  write_vtu_footer(std::ostream &out);
2068 
2075  template <int dim, int spacedim>
2076  void
2078  const std::vector<Patch<dim, spacedim>> &patches,
2079  const std::vector<std::string> &data_names,
2080  const std::vector<
2081  std::tuple<unsigned int,
2082  unsigned int,
2083  std::string,
2085  &nonscalar_data_ranges,
2086  const VtkFlags &flags,
2087  std::ostream &out);
2088 
2130  void
2132  std::ostream &out,
2133  const std::vector<std::string> &piece_names,
2134  const std::vector<std::string> &data_names,
2135  const std::vector<
2136  std::tuple<unsigned int,
2137  unsigned int,
2138  std::string,
2140  &nonscalar_data_ranges,
2141  const VtkFlags &flags);
2142 
2191  void
2193  std::ostream &out,
2194  const std::vector<std::pair<double, std::string>> &times_and_names);
2195 
2207  void
2208  write_visit_record(std::ostream &out,
2209  const std::vector<std::string> &piece_names);
2210 
2238  void
2239  write_visit_record(std::ostream &out,
2240  const std::vector<std::vector<std::string>> &piece_names);
2241 
2274  void
2276  std::ostream &out,
2277  const std::vector<std::pair<double, std::vector<std::string>>>
2278  &times_and_piece_names);
2279 
2300  template <int spacedim>
2301  void
2302  write_svg(
2303  const std::vector<Patch<2, spacedim>> &patches,
2304  const std::vector<std::string> &data_names,
2305  const std::vector<
2306  std::tuple<unsigned int,
2307  unsigned int,
2308  std::string,
2310  &nonscalar_data_ranges,
2311  const SvgFlags &flags,
2312  std::ostream &out);
2313 
2351  template <int dim, int spacedim>
2352  void
2354  const std::vector<Patch<dim, spacedim>> &patches,
2355  const std::vector<std::string> &data_names,
2356  const std::vector<
2357  std::tuple<unsigned int,
2358  unsigned int,
2359  std::string,
2361  &nonscalar_data_ranges,
2362  const Deal_II_IntermediateFlags &flags,
2363  std::ostream &out);
2364 
2373  template <int dim, int spacedim>
2374  void
2376  const std::vector<Patch<dim, spacedim>> &patches,
2377  const std::vector<std::string> &data_names,
2378  const std::vector<
2379  std::tuple<unsigned int,
2380  unsigned int,
2381  std::string,
2383  &nonscalar_data_ranges,
2384  const Deal_II_IntermediateFlags &flags,
2385  const std::string &filename,
2386  const MPI_Comm comm,
2387  const CompressionLevel compression);
2388 
2393  template <int dim, int spacedim>
2394  void
2395  write_hdf5_parallel(const std::vector<Patch<dim, spacedim>> &patches,
2396  const DataOutFilter &data_filter,
2397  const DataOutBase::Hdf5Flags &flags,
2398  const std::string &filename,
2399  const MPI_Comm comm);
2400 
2408  template <int dim, int spacedim>
2409  void
2410  write_hdf5_parallel(const std::vector<Patch<dim, spacedim>> &patches,
2411  const DataOutFilter &data_filter,
2412  const DataOutBase::Hdf5Flags &flags,
2413  const bool write_mesh_file,
2414  const std::string &mesh_filename,
2415  const std::string &solution_filename,
2416  const MPI_Comm comm);
2417 
2424  template <int dim, int spacedim>
2425  void
2427  const std::vector<Patch<dim, spacedim>> &patches,
2428  const std::vector<std::string> &data_names,
2429  const std::vector<
2430  std::tuple<unsigned int,
2431  unsigned int,
2432  std::string,
2434  &nonscalar_data_ranges,
2435  DataOutFilter &filtered_data);
2436 
2448  std::pair<unsigned int, unsigned int>
2449  determine_intermediate_format_dimensions(std::istream &input);
2450 
2462  OutputFormat
2463  parse_output_format(const std::string &format_name);
2464 
2470  std::string
2472 
2492  std::string
2493  default_suffix(const OutputFormat output_format);
2494 
2504  int,
2505  int,
2506  << "The number of points in this data set is " << arg1
2507  << ", but we expected " << arg2
2508  << " in each space direction.");
2513  "You are trying to write graphical data into a file, but "
2514  "no data is available in the intermediate format that "
2515  "the DataOutBase functions require. Did you forget to "
2516  "call a function such as DataOut::build_patches()?");
2521  "The error code of one of the Tecplot functions was "
2522  "not zero as expected.");
2527  char *,
2528  << "There was an error opening Tecplot file " << arg1
2529  << " for output.");
2530 
2532 } // namespace DataOutBase
2533 
2534 
2535 
2642 template <int dim, int spacedim = dim>
2644 {
2645 public:
2649  DataOutInterface();
2650 
2655  virtual ~DataOutInterface() = default;
2656 
2661  void
2662  write_dx(std::ostream &out) const;
2663 
2668  void
2669  write_eps(std::ostream &out) const;
2670 
2675  void
2676  write_gmv(std::ostream &out) const;
2677 
2682  void
2683  write_gnuplot(std::ostream &out) const;
2684 
2689  void
2690  write_povray(std::ostream &out) const;
2691 
2696  void
2697  write_tecplot(std::ostream &out) const;
2698 
2703  void
2704  write_ucd(std::ostream &out) const;
2705 
2717  void
2718  write_vtk(std::ostream &out) const;
2719 
2734  void
2735  write_vtu(std::ostream &out) const;
2736 
2745  void
2746  write_vtu_in_parallel(const std::string &filename, const MPI_Comm comm) const;
2747 
2781  void
2782  write_pvtu_record(std::ostream &out,
2783  const std::vector<std::string> &piece_names) const;
2784 
2843  std::string
2845  const std::string &directory,
2846  const std::string &filename_without_extension,
2847  const unsigned int counter,
2848  const MPI_Comm mpi_communicator,
2849  const unsigned int n_digits_for_counter = numbers::invalid_unsigned_int,
2850  const unsigned int n_groups = 0) const;
2851 
2856  void
2857  write_svg(std::ostream &out) const;
2858 
2869  void
2870  write_deal_II_intermediate(std::ostream &out) const;
2871 
2878  void
2880  const std::string &filename,
2881  const MPI_Comm comm,
2882  const DataOutBase::CompressionLevel compression) const;
2883 
2889  XDMFEntry
2890  create_xdmf_entry(const DataOutBase::DataOutFilter &data_filter,
2891  const std::string &h5_filename,
2892  const double cur_time,
2893  const MPI_Comm comm) const;
2894 
2900  XDMFEntry
2901  create_xdmf_entry(const DataOutBase::DataOutFilter &data_filter,
2902  const std::string &h5_mesh_filename,
2903  const std::string &h5_solution_filename,
2904  const double cur_time,
2905  const MPI_Comm comm) const;
2906 
2931  void
2932  write_xdmf_file(const std::vector<XDMFEntry> &entries,
2933  const std::string &filename,
2934  const MPI_Comm comm) const;
2935 
2950  void
2952  const std::string &filename,
2953  const MPI_Comm comm) const;
2954 
2962  void
2964  const bool write_mesh_file,
2965  const std::string &mesh_filename,
2966  const std::string &solution_filename,
2967  const MPI_Comm comm) const;
2968 
2975  void
2976  write_filtered_data(DataOutBase::DataOutFilter &filtered_data) const;
2977 
2978 
2987  void
2988  write(std::ostream &out,
2989  const DataOutBase::OutputFormat output_format =
2991 
2996  void
2998 
2999 
3004  template <typename FlagType>
3005  void
3006  set_flags(const FlagType &flags);
3007 
3008 
3016  std::string
3017  default_suffix(const DataOutBase::OutputFormat output_format =
3019 
3034  static void
3036 
3044  void
3046 
3052  std::size_t
3053  memory_consumption() const;
3054 
3055 protected:
3063  virtual const std::vector<DataOutBase::Patch<dim, spacedim>> &
3064  get_patches() const = 0;
3065 
3070  virtual std::vector<std::string>
3071  get_dataset_names() const = 0;
3072 
3091  virtual std::vector<
3092  std::tuple<unsigned int,
3093  unsigned int,
3094  std::string,
3096  get_nonscalar_data_ranges() const;
3097 
3104  void
3105  validate_dataset_names() const;
3106 
3107 
3113  unsigned int default_subdivisions;
3114 
3115 private:
3122 
3128 
3134 
3140 
3146 
3152 
3158 
3164 
3170 
3176 
3182 
3188 };
3189 
3190 
3191 
3237 template <int dim, int spacedim = dim>
3238 class DataOutReader : public DataOutInterface<dim, spacedim>
3239 {
3240 public:
3246  void
3247  read(std::istream &in);
3248 
3254  void
3255  read_whole_parallel_file(std::istream &in);
3256 
3278  void
3279  merge(const DataOutReader<dim, spacedim> &other);
3280 
3285  "You are trying to merge two sets of patches for which the "
3286  "declared names of the variables do not match.");
3291  "You are trying to merge two sets of patches for which the "
3292  "number of subdivisions or the number of vector components "
3293  "do not match.");
3298  int,
3299  int,
3300  int,
3301  int,
3302  << "Either the dimensions <" << arg1 << "> and <" << arg2
3303  << "> or the space dimensions <" << arg3 << "> and <" << arg4
3304  << "> do not match!");
3305 
3306 protected:
3315  virtual const std::vector<::DataOutBase::Patch<dim, spacedim>> &
3316  get_patches() const override;
3317 
3324  virtual std::vector<std::string>
3325  get_dataset_names() const override;
3326 
3345  virtual std::vector<
3346  std::tuple<unsigned int,
3347  unsigned int,
3348  std::string,
3350  get_nonscalar_data_ranges() const override;
3351 
3352 private:
3357  std::vector<::DataOutBase::Patch<dim, spacedim>> patches;
3358  std::vector<std::string> dataset_names;
3359 
3364  std::vector<
3365  std::tuple<unsigned int,
3366  unsigned int,
3367  std::string,
3370 };
3371 
3372 
3373 
3382 {
3383 public:
3387  XDMFEntry();
3388 
3394  XDMFEntry(const std::string &filename,
3395  const double time,
3396  const std::uint64_t nodes,
3397  const std::uint64_t cells,
3398  const unsigned int dim,
3399  const ReferenceCell &cell_type);
3400 
3406  XDMFEntry(const std::string &filename,
3407  const double time,
3408  const std::uint64_t nodes,
3409  const std::uint64_t cells,
3410  const unsigned int dim);
3411 
3417  XDMFEntry(const std::string &mesh_filename,
3418  const std::string &solution_filename,
3419  const double time,
3420  const std::uint64_t nodes,
3421  const std::uint64_t cells,
3422  const unsigned int dim);
3423 
3428  XDMFEntry(const std::string &mesh_filename,
3429  const std::string &solution_filename,
3430  const double time,
3431  const std::uint64_t nodes,
3432  const std::uint64_t cells,
3433  const unsigned int dim,
3434  const ReferenceCell &cell_type);
3435 
3442  XDMFEntry(const std::string &mesh_filename,
3443  const std::string &solution_filename,
3444  const double time,
3445  const std::uint64_t nodes,
3446  const std::uint64_t cells,
3447  const unsigned int dim,
3448  const unsigned int spacedim);
3449 
3453  XDMFEntry(const std::string &mesh_filename,
3454  const std::string &solution_filename,
3455  const double time,
3456  const std::uint64_t nodes,
3457  const std::uint64_t cells,
3458  const unsigned int dim,
3459  const unsigned int spacedim,
3460  const ReferenceCell &cell_type);
3461 
3465  void
3466  add_attribute(const std::string &attr_name, const unsigned int dimension);
3467 
3473  template <class Archive>
3474  void
3475  serialize(Archive &ar, const unsigned int /*version*/)
3476  {
3479  }
3480 
3485  std::string
3486  get_xdmf_content(const unsigned int indent_level) const;
3487 
3495  std::string
3496  get_xdmf_content(const unsigned int indent_level,
3497  const ReferenceCell &reference_cell) const;
3498 
3499 private:
3503  bool valid;
3504 
3508  std::string h5_sol_filename;
3509 
3513  std::string h5_mesh_filename;
3514 
3518  double entry_time;
3519 
3523  std::uint64_t num_nodes;
3524 
3528  std::uint64_t num_cells;
3529 
3533  unsigned int dimension;
3534 
3539  unsigned int space_dimension;
3540 
3546 
3550  std::map<std::string, unsigned int> attribute_dims;
3551 };
3552 
3553 
3554 
3555 /* -------------------- inline functions ------------------- */
3556 
3557 namespace DataOutBase
3558 {
3559  inline bool
3561  {
3562  return (red == green) && (red == blue);
3563  }
3564 
3565 
3566  /* -------------------- template functions ------------------- */
3567 
3574  template <int dim, int spacedim>
3575  std::ostream &
3576  operator<<(std::ostream &out, const Patch<dim, spacedim> &patch);
3577 
3578 
3579 
3586  template <int dim, int spacedim>
3587  std::istream &
3588  operator>>(std::istream &in, Patch<dim, spacedim> &patch);
3589 } // namespace DataOutBase
3590 
3591 
3593 
3594 #endif
const double * get_data_set(const unsigned int set_num) const
std::map< unsigned int, unsigned int > filtered_points
std::string get_data_set_name(const unsigned int set_num) const
std::multimap< Point< 3 >, unsigned int, Point3Comp > Map3DPoint
void internal_add_cell(const unsigned int cell_index, const unsigned int pt_index)
void write_cell(const unsigned int index, const unsigned int start, const std::array< unsigned int, dim > &offsets)
std::vector< unsigned int > data_set_dims
unsigned int n_nodes() const
void fill_node_data(std::vector< double > &node_data) const
std::vector< std::vector< double > > data_sets
unsigned int get_data_set_dim(const unsigned int set_num) const
void write_data_set(const std::string &name, const unsigned int dimension, const unsigned int set_num, const Table< 2, double > &data_vectors)
std::map< unsigned int, unsigned int > filtered_cells
void write_cell_single(const unsigned int index, const unsigned int start, const unsigned int n_points, const ReferenceCell &reference_cell)
std::vector< std::string > data_set_names
void fill_cell_data(const unsigned int local_node_offset, std::vector< unsigned int > &cell_data) const
void write_point(const unsigned int index, const Point< dim > &p)
unsigned int n_cells() const
DataOutBase::DataOutFilterFlags flags
unsigned int n_data_sets() const
XDMFEntry create_xdmf_entry(const DataOutBase::DataOutFilter &data_filter, const std::string &h5_filename, const double cur_time, const MPI_Comm comm) const
virtual std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > get_nonscalar_data_ranges() const
unsigned int default_subdivisions
void parse_parameters(ParameterHandler &prm)
DataOutBase::Deal_II_IntermediateFlags deal_II_intermediate_flags
virtual std::vector< std::string > get_dataset_names() const =0
void write_filtered_data(DataOutBase::DataOutFilter &filtered_data) const
void write_pvtu_record(std::ostream &out, const std::vector< std::string > &piece_names) const
static void declare_parameters(ParameterHandler &prm)
DataOutBase::TecplotFlags tecplot_flags
void write_ucd(std::ostream &out) const
void write_povray(std::ostream &out) const
std::string default_suffix(const DataOutBase::OutputFormat output_format=DataOutBase::default_format) const
void write_xdmf_file(const std::vector< XDMFEntry > &entries, const std::string &filename, const MPI_Comm comm) const
DataOutBase::GmvFlags gmv_flags
std::size_t memory_consumption() const
void set_default_format(const DataOutBase::OutputFormat default_format)
DataOutBase::Hdf5Flags hdf5_flags
DataOutBase::PovrayFlags povray_flags
DataOutBase::UcdFlags ucd_flags
void write(std::ostream &out, const DataOutBase::OutputFormat output_format=DataOutBase::default_format) const
virtual const std::vector< DataOutBase::Patch< dim, spacedim > > & get_patches() const =0
DataOutBase::OutputFormat default_fmt
DataOutBase::GnuplotFlags gnuplot_flags
void write_gnuplot(std::ostream &out) const
virtual ~DataOutInterface()=default
void write_vtu(std::ostream &out) const
void write_hdf5_parallel(const DataOutBase::DataOutFilter &data_filter, const std::string &filename, const MPI_Comm comm) const
void write_tecplot(std::ostream &out) const
void write_deal_II_intermediate_in_parallel(const std::string &filename, const MPI_Comm comm, const DataOutBase::CompressionLevel compression) const
void write_svg(std::ostream &out) const
void write_vtu_in_parallel(const std::string &filename, const MPI_Comm comm) const
DataOutBase::SvgFlags svg_flags
void validate_dataset_names() const
void set_flags(const FlagType &flags)
void write_vtk(std::ostream &out) const
DataOutBase::VtkFlags vtk_flags
void write_gmv(std::ostream &out) const
std::string write_vtu_with_pvtu_record(const std::string &directory, const std::string &filename_without_extension, const unsigned int counter, const MPI_Comm mpi_communicator, const unsigned int n_digits_for_counter=numbers::invalid_unsigned_int, const unsigned int n_groups=0) const
DataOutBase::EpsFlags eps_flags
void write_eps(std::ostream &out) const
DataOutBase::DXFlags dx_flags
void write_dx(std::ostream &out) const
void write_deal_II_intermediate(std::ostream &out) const
void merge(const DataOutReader< dim, spacedim > &other)
void read_whole_parallel_file(std::istream &in)
std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > nonscalar_data_ranges
virtual const std::vector<::DataOutBase::Patch< dim, spacedim > > & get_patches() const override
void read(std::istream &in)
virtual std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > get_nonscalar_data_ranges() const override
virtual std::vector< std::string > get_dataset_names() const override
std::vector< std::string > dataset_names
std::vector<::DataOutBase::Patch< dim, spacedim > > patches
ReferenceCell cell_type
std::uint64_t num_nodes
unsigned int dimension
std::string h5_sol_filename
void serialize(Archive &ar, const unsigned int)
double entry_time
std::string h5_mesh_filename
std::string get_xdmf_content(const unsigned int indent_level) const
void add_attribute(const std::string &attr_name, const unsigned int dimension)
std::uint64_t num_cells
std::map< std::string, unsigned int > attribute_dims
unsigned int space_dimension
#define DEAL_II_DEPRECATED
Definition: config.h:178
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
unsigned int cell_index
static ::ExceptionBase & ExcIncompatiblePatchLists()
static ::ExceptionBase & ExcIncompatibleDimensions(int arg1, int arg2, int arg3, int arg4)
static ::ExceptionBase & ExcErrorOpeningTecplotFile(char *arg1)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition: exceptions.h:586
static ::ExceptionBase & ExcInvalidDatasetSize(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidCombinationOfDimensions(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:540
static ::ExceptionBase & ExcNotEnoughSpaceDimensionLabels()
static ::ExceptionBase & ExcIncompatibleDatasetNames()
static ::ExceptionBase & ExcTecplotAPIError()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:495
static ::ExceptionBase & ExcNoPatches()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:517
std::pair< unsigned int, unsigned int > determine_intermediate_format_dimensions(std::istream &input)
std::ostream & operator<<(std::ostream &out, const Patch< dim, spacedim > &patch)
void write_vtk(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const VtkFlags &flags, std::ostream &out)
void write_vtu(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const VtkFlags &flags, std::ostream &out)
void write_gnuplot(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const GnuplotFlags &flags, std::ostream &out)
void write_deal_II_intermediate_in_parallel(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const Deal_II_IntermediateFlags &flags, const std::string &filename, const MPI_Comm comm, const CompressionLevel compression)
void write_hdf5_parallel(const std::vector< Patch< dim, spacedim >> &patches, const DataOutFilter &data_filter, const DataOutBase::Hdf5Flags &flags, const std::string &filename, const MPI_Comm comm)
void write_pvtu_record(std::ostream &out, const std::vector< std::string > &piece_names, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const VtkFlags &flags)
void write_filtered_data(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, DataOutFilter &filtered_data)
void write_ucd(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const UcdFlags &flags, std::ostream &out)
void write_vtu_header(std::ostream &out, const VtkFlags &flags)
void write_pvd_record(std::ostream &out, const std::vector< std::pair< double, std::string >> &times_and_names)
void write_deal_II_intermediate(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const Deal_II_IntermediateFlags &flags, std::ostream &out)
void write_dx(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const DXFlags &flags, std::ostream &out)
void write_eps(const std::vector< Patch< 2, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const EpsFlags &flags, std::ostream &out)
void write_vtu_footer(std::ostream &out)
void write_tecplot(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const TecplotFlags &flags, std::ostream &out)
OutputFormat parse_output_format(const std::string &format_name)
void write_vtu_main(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const VtkFlags &flags, std::ostream &out)
void write_svg(const std::vector< Patch< 2, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const SvgFlags &flags, std::ostream &out)
std::istream & operator>>(std::istream &in, Patch< dim, spacedim > &patch)
std::string get_output_format_names()
void write_visit_record(std::ostream &out, const std::vector< std::string > &piece_names)
void write_povray(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const PovrayFlags &flags, std::ostream &out)
std::string default_suffix(const OutputFormat output_format)
void write_gmv(const std::vector< Patch< dim, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const GmvFlags &flags, std::ostream &out)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
static const types::blas_int one
static const unsigned int invalid_unsigned_int
Definition: types.h:221
void parse_parameters(const ParameterHandler &prm)
DXFlags(const bool write_neighbors=false, const bool int_binary=false, const bool coordinates_binary=false, const bool data_binary=false)
static void declare_parameters(ParameterHandler &prm)
static void declare_parameters(ParameterHandler &prm)
std::size_t memory_consumption() const
void parse_parameters(const ParameterHandler &prm)
DataOutFilterFlags(const bool filter_duplicate_vertices=false, const bool xdmf_hdf5_output=false)
bool operator()(const Point< 3 > &one, const Point< 3 > &two) const
static const unsigned int format_version
static void declare_parameters(ParameterHandler &prm)
static RgbValues default_color_function(const double value, const double min_value, const double max_value)
void parse_parameters(const ParameterHandler &prm)
ColorFunction color_function
RgbValues(*)(const double value, const double min_value, const double max_value) ColorFunction
static RgbValues grey_scale_color_function(const double value, const double min_value, const double max_value)
EpsFlags(const unsigned int height_vector=0, const unsigned int color_vector=0, const SizeType size_type=width, const unsigned int size=300, const double line_width=0.5, const double azimut_angle=60, const double turn_angle=30, const double z_scaling=1.0, const bool draw_mesh=true, const bool draw_cells=true, const bool shade_cells=true, const ColorFunction color_function=&default_color_function)
unsigned int color_vector
static RgbValues reverse_grey_scale_color_function(const double value, const double min_value, const double max_value)
@ width
Scale to given width.
@ height
Scale to given height.
unsigned int height_vector
std::vector< std::string > space_dimension_labels
std::size_t memory_consumption() const
DataOutBase::CompressionLevel compression_level
Hdf5Flags(const CompressionLevel compression_level=CompressionLevel::best_speed)
std::size_t memory_consumption() const
static void declare_parameters(ParameterHandler &prm)
void parse_parameters(const ParameterHandler &prm)
static const unsigned int n_subdivisions
static const ReferenceCell reference_cell
std::size_t memory_consumption() const
unsigned int patch_index
ReferenceCell reference_cell
Table< 2, float > data
static const unsigned int no_neighbor
bool operator==(const Patch &patch) const
void swap(Patch< dim, spacedim > &other_patch) noexcept
static const unsigned int space_dim
unsigned int n_subdivisions
std::array< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > vertices
std::array< unsigned int, GeometryInfo< dim >::faces_per_cell > neighbors
static void declare_parameters(ParameterHandler &prm)
PovrayFlags(const bool smooth=false, const bool bicubic_patch=false, const bool external_data=false)
void parse_parameters(const ParameterHandler &prm)
unsigned int height_vector
SvgFlags(const unsigned int height_vector=0, const int azimuth_angle=37, const int polar_angle=45, const unsigned int line_thickness=1, const bool margin=true, const bool draw_colorbar=true)
unsigned int line_thickness
std::size_t memory_consumption() const
TecplotFlags(const char *zone_name=nullptr, const double solution_time=-1.0)
void parse_parameters(const ParameterHandler &prm)
static void declare_parameters(ParameterHandler &prm)
UcdFlags(const bool write_preamble=false)
VtkFlags(const double time=std::numeric_limits< double >::min(), const unsigned int cycle=std::numeric_limits< unsigned int >::min(), const bool print_date_and_time=true, const CompressionLevel compression_level=CompressionLevel::best_speed, const bool write_higher_order_cells=false, const std::map< std::string, std::string > &physical_units={})
std::map< std::string, std::string > physical_units
static const DataOutBase::CompressionLevel default_compression
static const DataOutBase::CompressionLevel best_compression
static const DataOutBase::CompressionLevel best_speed
DataOutBase::CompressionLevel compression_level
static const DataOutBase::CompressionLevel no_compression
const MPI_Comm comm