Reference documentation for deal.II version GIT eef19498cf 2022-10-05 14:45:02+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 - 2022 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 <string>
36 #include <tuple>
37 #include <typeinfo>
38 #include <vector>
39 
40 // Only include the Tecplot API header if the appropriate files
41 // were detected by configure
42 #ifdef DEAL_II_HAVE_TECPLOT
43 # include <string.h>
44 
45 # include "TECIO.h"
46 #endif
47 
48 #include <ostream>
49 
51 
52 // Forward declarations
53 #ifndef DOXYGEN
54 class ParameterHandler;
55 class XDMFEntry;
56 #endif
57 
207 namespace DataOutBase
208 {
216  enum class CompressionLevel
217  {
225  best_speed,
239  plain_text
240  };
241 
242 
265  template <int dim, int spacedim = dim>
266  struct Patch
267  {
271  static const unsigned int space_dim = spacedim;
272 
286 
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);
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);
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 TecplotFlags : public OutputFlagsBase<TecplotFlags>
1077  {
1082  const char *zone_name;
1083 
1090 
1094  TecplotFlags(const char * zone_name = nullptr,
1095  const double solution_time = -1.0);
1096 
1101  std::size_t
1102  memory_consumption() const;
1103  };
1104 
1110  struct VtkFlags : public OutputFlagsBase<VtkFlags>
1111  {
1122  double time;
1123 
1134  unsigned int cycle;
1135 
1143 
1152 
1161 
1167 
1186 
1208  std::map<std::string, std::string> physical_units;
1209 
1214  explicit VtkFlags(
1215  const double time = std::numeric_limits<double>::min(),
1216  const unsigned int cycle = std::numeric_limits<unsigned int>::min(),
1217  const bool print_date_and_time = true,
1220  const bool write_higher_order_cells = false,
1221  const std::map<std::string, std::string> &physical_units = {});
1222  };
1223 
1224 
1230  struct SvgFlags : public OutputFlagsBase<SvgFlags>
1231  {
1235  unsigned int height;
1236 
1241  unsigned int width;
1242 
1249  unsigned int height_vector;
1250 
1255 
1256  unsigned int line_thickness;
1257 
1261  bool margin;
1262 
1267 
1271  SvgFlags(const unsigned int height_vector = 0,
1272  const int azimuth_angle = 37,
1273  const int polar_angle = 45,
1274  const unsigned int line_thickness = 1,
1275  const bool margin = true,
1276  const bool draw_colorbar = true);
1277  };
1278 
1279 
1287  : public OutputFlagsBase<Deal_II_IntermediateFlags>
1288  {
1295  static const unsigned int format_version;
1296  };
1297 
1304  {
1315 
1321 
1325  DataOutFilterFlags(const bool filter_duplicate_vertices = false,
1326  const bool xdmf_hdf5_output = false);
1327 
1332  static void
1334 
1341  void
1342  parse_parameters(const ParameterHandler &prm);
1343 
1348  std::size_t
1350  };
1351 
1384  {
1385  public:
1389  DataOutFilter();
1390 
1396 
1402  template <int dim>
1403  void
1404  write_point(const unsigned int index, const Point<dim> &p);
1405 
1409  template <int dim>
1410  void
1411  write_cell(const unsigned int index,
1412  const unsigned int start,
1413  const unsigned int d1,
1414  const unsigned int d2,
1415  const unsigned int d3);
1416 
1421  void
1422  write_cell_single(const unsigned int index,
1423  const unsigned int start,
1424  const unsigned int n_points,
1425  const ReferenceCell &reference_cell);
1426 
1433  void
1434  write_data_set(const std::string & name,
1435  const unsigned int dimension,
1436  const unsigned int set_num,
1437  const Table<2, double> &data_vectors);
1438 
1443  void
1444  fill_node_data(std::vector<double> &node_data) const;
1445 
1450  void
1451  fill_cell_data(const unsigned int local_node_offset,
1452  std::vector<unsigned int> &cell_data) const;
1453 
1457  std::string
1458  get_data_set_name(const unsigned int set_num) const;
1459 
1463  unsigned int
1464  get_data_set_dim(const unsigned int set_num) const;
1465 
1470  const double *
1471  get_data_set(const unsigned int set_num) const;
1472 
1477  unsigned int
1478  n_nodes() const;
1479 
1484  unsigned int
1485  n_cells() const;
1486 
1491  unsigned int
1492  n_data_sets() const;
1493 
1497  void
1498  flush_points();
1499 
1503  void
1504  flush_cells();
1505 
1506 
1507  private:
1511  struct Point3Comp
1512  {
1513  bool
1514  operator()(const Point<3> &one, const Point<3> &two) const
1515  {
1516  /*
1517  * The return statement below is an optimized version of the following
1518  * code:
1519  *
1520  * for (unsigned int d=0; d<3; ++d)
1521  * {
1522  * if (one(d) < two(d))
1523  * return true;
1524  * else if (one(d) > two(d))
1525  * return false;
1526  * }
1527  * return false;
1528  */
1529 
1530  return (one(0) < two(0) ||
1531  (!(two(0) < one(0)) &&
1532  (one(1) < two(1) || (!(two(1) < one(1)) && one(2) < two(2)))));
1533  }
1534  };
1535 
1536  using Map3DPoint = std::multimap<Point<3>, unsigned int, Point3Comp>;
1537 
1542 
1549  unsigned int node_dim;
1550 
1555  unsigned int num_cells;
1556 
1561 
1565  std::map<unsigned int, unsigned int> filtered_points;
1566 
1570  std::map<unsigned int, unsigned int> filtered_cells;
1571 
1575  std::vector<std::string> data_set_names;
1576 
1580  std::vector<unsigned int> data_set_dims;
1581 
1585  std::vector<std::vector<double>> data_sets;
1586 
1590  void
1591  internal_add_cell(const unsigned int cell_index,
1592  const unsigned int pt_index);
1593  };
1594 
1595 
1600  {
1605 
1610 
1615 
1620 
1625 
1630 
1635 
1640 
1645 
1653 
1658 
1663 
1668 
1673 
1677  hdf5
1678  };
1679 
1680 
1684  template <int dim, int spacedim>
1685  void
1686  write_dx(
1687  const std::vector<Patch<dim, spacedim>> &patches,
1688  const std::vector<std::string> & data_names,
1689  const std::vector<
1690  std::tuple<unsigned int,
1691  unsigned int,
1692  std::string,
1694  & nonscalar_data_ranges,
1695  const DXFlags &flags,
1696  std::ostream & out);
1697 
1742  template <int spacedim>
1743  void
1744  write_eps(
1745  const std::vector<Patch<2, spacedim>> &patches,
1746  const std::vector<std::string> & data_names,
1747  const std::vector<
1748  std::tuple<unsigned int,
1749  unsigned int,
1750  std::string,
1752  & nonscalar_data_ranges,
1753  const EpsFlags &flags,
1754  std::ostream & out);
1755 
1761  template <int dim, int spacedim>
1762  void
1763  write_eps(
1764  const std::vector<Patch<dim, spacedim>> &patches,
1765  const std::vector<std::string> & data_names,
1766  const std::vector<
1767  std::tuple<unsigned int,
1768  unsigned int,
1769  std::string,
1771  & nonscalar_data_ranges,
1772  const EpsFlags &flags,
1773  std::ostream & out);
1774 
1775 
1785  template <int dim, int spacedim>
1786  void
1787  write_gmv(
1788  const std::vector<Patch<dim, spacedim>> &patches,
1789  const std::vector<std::string> & data_names,
1790  const std::vector<
1791  std::tuple<unsigned int,
1792  unsigned int,
1793  std::string,
1795  & nonscalar_data_ranges,
1796  const GmvFlags &flags,
1797  std::ostream & out);
1798 
1854  template <int dim, int spacedim>
1855  void
1856  write_gnuplot(
1857  const std::vector<Patch<dim, spacedim>> &patches,
1858  const std::vector<std::string> & data_names,
1859  const std::vector<
1860  std::tuple<unsigned int,
1861  unsigned int,
1862  std::string,
1864  & nonscalar_data_ranges,
1865  const GnuplotFlags &flags,
1866  std::ostream & out);
1867 
1913  template <int dim, int spacedim>
1914  void
1915  write_povray(
1916  const std::vector<Patch<dim, spacedim>> &patches,
1917  const std::vector<std::string> & data_names,
1918  const std::vector<
1919  std::tuple<unsigned int,
1920  unsigned int,
1921  std::string,
1923  & nonscalar_data_ranges,
1924  const PovrayFlags &flags,
1925  std::ostream & out);
1926 
1933  template <int dim, int spacedim>
1934  void
1935  write_tecplot(
1936  const std::vector<Patch<dim, spacedim>> &patches,
1937  const std::vector<std::string> & data_names,
1938  const std::vector<
1939  std::tuple<unsigned int,
1940  unsigned int,
1941  std::string,
1943  & nonscalar_data_ranges,
1944  const TecplotFlags &flags,
1945  std::ostream & out);
1946 
1961  template <int dim, int spacedim>
1962  void
1963  write_ucd(
1964  const std::vector<Patch<dim, spacedim>> &patches,
1965  const std::vector<std::string> & data_names,
1966  const std::vector<
1967  std::tuple<unsigned int,
1968  unsigned int,
1969  std::string,
1971  & nonscalar_data_ranges,
1972  const UcdFlags &flags,
1973  std::ostream & out);
1974 
1994  template <int dim, int spacedim>
1995  void
1996  write_vtk(
1997  const std::vector<Patch<dim, spacedim>> &patches,
1998  const std::vector<std::string> & data_names,
1999  const std::vector<
2000  std::tuple<unsigned int,
2001  unsigned int,
2002  std::string,
2004  & nonscalar_data_ranges,
2005  const VtkFlags &flags,
2006  std::ostream & out);
2007 
2008 
2032  template <int dim, int spacedim>
2033  void
2034  write_vtu(
2035  const std::vector<Patch<dim, spacedim>> &patches,
2036  const std::vector<std::string> & data_names,
2037  const std::vector<
2038  std::tuple<unsigned int,
2039  unsigned int,
2040  std::string,
2042  & nonscalar_data_ranges,
2043  const VtkFlags &flags,
2044  std::ostream & out);
2045 
2051  void
2052  write_vtu_header(std::ostream &out, const VtkFlags &flags);
2053 
2060  void
2061  write_vtu_footer(std::ostream &out);
2062 
2069  template <int dim, int spacedim>
2070  void
2072  const std::vector<Patch<dim, spacedim>> &patches,
2073  const std::vector<std::string> & data_names,
2074  const std::vector<
2075  std::tuple<unsigned int,
2076  unsigned int,
2077  std::string,
2079  & nonscalar_data_ranges,
2080  const VtkFlags &flags,
2081  std::ostream & out);
2082 
2124  void
2126  std::ostream & out,
2127  const std::vector<std::string> &piece_names,
2128  const std::vector<std::string> &data_names,
2129  const std::vector<
2130  std::tuple<unsigned int,
2131  unsigned int,
2132  std::string,
2134  & nonscalar_data_ranges,
2135  const VtkFlags &flags);
2136 
2185  void
2187  std::ostream & out,
2188  const std::vector<std::pair<double, std::string>> &times_and_names);
2189 
2201  void
2202  write_visit_record(std::ostream & out,
2203  const std::vector<std::string> &piece_names);
2204 
2232  void
2233  write_visit_record(std::ostream & out,
2234  const std::vector<std::vector<std::string>> &piece_names);
2235 
2268  void
2270  std::ostream &out,
2271  const std::vector<std::pair<double, std::vector<std::string>>>
2272  &times_and_piece_names);
2273 
2294  template <int spacedim>
2295  void
2296  write_svg(
2297  const std::vector<Patch<2, spacedim>> &patches,
2298  const std::vector<std::string> & data_names,
2299  const std::vector<
2300  std::tuple<unsigned int,
2301  unsigned int,
2302  std::string,
2304  & nonscalar_data_ranges,
2305  const SvgFlags &flags,
2306  std::ostream & out);
2307 
2345  template <int dim, int spacedim>
2346  void
2348  const std::vector<Patch<dim, spacedim>> &patches,
2349  const std::vector<std::string> & data_names,
2350  const std::vector<
2351  std::tuple<unsigned int,
2352  unsigned int,
2353  std::string,
2355  & nonscalar_data_ranges,
2356  const Deal_II_IntermediateFlags &flags,
2357  std::ostream & out);
2358 
2367  template <int dim, int spacedim>
2368  void
2370  const std::vector<Patch<dim, spacedim>> &patches,
2371  const std::vector<std::string> & data_names,
2372  const std::vector<
2373  std::tuple<unsigned int,
2374  unsigned int,
2375  std::string,
2377  & nonscalar_data_ranges,
2378  const Deal_II_IntermediateFlags &flags,
2379  const std::string & filename,
2380  const MPI_Comm & comm,
2381  const CompressionLevel compression);
2382 
2387  template <int dim, int spacedim>
2388  void
2389  write_hdf5_parallel(const std::vector<Patch<dim, spacedim>> &patches,
2390  const DataOutFilter & data_filter,
2391  const std::string & filename,
2392  const MPI_Comm & comm);
2393 
2401  template <int dim, int spacedim>
2402  void
2403  write_hdf5_parallel(const std::vector<Patch<dim, spacedim>> &patches,
2404  const DataOutFilter & data_filter,
2405  const bool write_mesh_file,
2406  const std::string & mesh_filename,
2407  const std::string &solution_filename,
2408  const MPI_Comm & comm);
2409 
2416  template <int dim, int spacedim>
2417  void
2419  const std::vector<Patch<dim, spacedim>> &patches,
2420  const std::vector<std::string> & data_names,
2421  const std::vector<
2422  std::tuple<unsigned int,
2423  unsigned int,
2424  std::string,
2426  & nonscalar_data_ranges,
2427  DataOutFilter &filtered_data);
2428 
2440  std::pair<unsigned int, unsigned int>
2441  determine_intermediate_format_dimensions(std::istream &input);
2442 
2454  OutputFormat
2455  parse_output_format(const std::string &format_name);
2456 
2462  std::string
2464 
2485  std::string
2486  default_suffix(const OutputFormat output_format);
2487 
2497  int,
2498  int,
2499  << "The number of points in this data set is " << arg1
2500  << ", but we expected " << arg2
2501  << " in each space direction.");
2506  "You are trying to write graphical data into a file, but "
2507  "no data is available in the intermediate format that "
2508  "the DataOutBase functions require. Did you forget to "
2509  "call a function such as DataOut::build_patches()?");
2514  "The error code of one of the Tecplot functions was "
2515  "not zero as expected.");
2520  char *,
2521  << "There was an error opening Tecplot file " << arg1
2522  << " for output.");
2523 
2525 } // namespace DataOutBase
2526 
2527 
2528 
2635 template <int dim, int spacedim = dim>
2637 {
2638 public:
2642  DataOutInterface();
2643 
2648  virtual ~DataOutInterface() = default;
2649 
2654  void
2655  write_dx(std::ostream &out) const;
2656 
2661  void
2662  write_eps(std::ostream &out) const;
2663 
2668  void
2669  write_gmv(std::ostream &out) const;
2670 
2675  void
2676  write_gnuplot(std::ostream &out) const;
2677 
2682  void
2683  write_povray(std::ostream &out) const;
2684 
2689  void
2690  write_tecplot(std::ostream &out) const;
2691 
2696  void
2697  write_ucd(std::ostream &out) const;
2698 
2710  void
2711  write_vtk(std::ostream &out) const;
2712 
2727  void
2728  write_vtu(std::ostream &out) const;
2729 
2738  void
2739  write_vtu_in_parallel(const std::string &filename,
2740  const MPI_Comm & comm) const;
2741 
2775  void
2776  write_pvtu_record(std::ostream & out,
2777  const std::vector<std::string> &piece_names) const;
2778 
2837  std::string
2839  const std::string &directory,
2840  const std::string &filename_without_extension,
2841  const unsigned int counter,
2842  const MPI_Comm & mpi_communicator,
2843  const unsigned int n_digits_for_counter = numbers::invalid_unsigned_int,
2844  const unsigned int n_groups = 0) const;
2845 
2850  void
2851  write_svg(std::ostream &out) const;
2852 
2863  void
2864  write_deal_II_intermediate(std::ostream &out) const;
2865 
2872  void
2874  const std::string & filename,
2875  const MPI_Comm & comm,
2876  const DataOutBase::CompressionLevel compression) const;
2877 
2883  XDMFEntry
2884  create_xdmf_entry(const DataOutBase::DataOutFilter &data_filter,
2885  const std::string & h5_filename,
2886  const double cur_time,
2887  const MPI_Comm & comm) const;
2888 
2894  XDMFEntry
2895  create_xdmf_entry(const DataOutBase::DataOutFilter &data_filter,
2896  const std::string & h5_mesh_filename,
2897  const std::string & h5_solution_filename,
2898  const double cur_time,
2899  const MPI_Comm & comm) const;
2900 
2925  void
2926  write_xdmf_file(const std::vector<XDMFEntry> &entries,
2927  const std::string & filename,
2928  const MPI_Comm & comm) const;
2929 
2944  void
2946  const std::string & filename,
2947  const MPI_Comm & comm) const;
2948 
2956  void
2958  const bool write_mesh_file,
2959  const std::string & mesh_filename,
2960  const std::string & solution_filename,
2961  const MPI_Comm & comm) const;
2962 
2969  void
2970  write_filtered_data(DataOutBase::DataOutFilter &filtered_data) const;
2971 
2972 
2981  void
2982  write(std::ostream & out,
2983  const DataOutBase::OutputFormat output_format =
2985 
2990  void
2992 
2993 
2998  template <typename FlagType>
2999  void
3000  set_flags(const FlagType &flags);
3001 
3002 
3010  std::string
3011  default_suffix(const DataOutBase::OutputFormat output_format =
3013 
3028  static void
3030 
3038  void
3040 
3046  std::size_t
3047  memory_consumption() const;
3048 
3049 protected:
3057  virtual const std::vector<DataOutBase::Patch<dim, spacedim>> &
3058  get_patches() const = 0;
3059 
3064  virtual std::vector<std::string>
3065  get_dataset_names() const = 0;
3066 
3085  virtual std::vector<
3086  std::tuple<unsigned int,
3087  unsigned int,
3088  std::string,
3090  get_nonscalar_data_ranges() const;
3091 
3098  void
3099  validate_dataset_names() const;
3100 
3101 
3107  unsigned int default_subdivisions;
3108 
3109 private:
3116 
3122 
3128 
3134 
3140 
3146 
3152 
3158 
3164 
3170 
3176 };
3177 
3178 
3179 
3225 template <int dim, int spacedim = dim>
3226 class DataOutReader : public DataOutInterface<dim, spacedim>
3227 {
3228 public:
3234  void
3235  read(std::istream &in);
3236 
3242  void
3243  read_whole_parallel_file(std::istream &in);
3244 
3266  void
3267  merge(const DataOutReader<dim, spacedim> &other);
3268 
3273  "You are trying to merge two sets of patches for which the "
3274  "declared names of the variables do not match.");
3279  "You are trying to merge two sets of patches for which the "
3280  "number of subdivisions or the number of vector components "
3281  "do not match.");
3286  int,
3287  int,
3288  int,
3289  int,
3290  << "Either the dimensions <" << arg1 << "> and <" << arg2
3291  << "> or the space dimensions <" << arg3 << "> and <" << arg4
3292  << "> do not match!");
3293 
3294 protected:
3303  virtual const std::vector<::DataOutBase::Patch<dim, spacedim>> &
3304  get_patches() const override;
3305 
3312  virtual std::vector<std::string>
3313  get_dataset_names() const override;
3314 
3333  virtual std::vector<
3334  std::tuple<unsigned int,
3335  unsigned int,
3336  std::string,
3338  get_nonscalar_data_ranges() const override;
3339 
3340 private:
3345  std::vector<::DataOutBase::Patch<dim, spacedim>> patches;
3346  std::vector<std::string> dataset_names;
3347 
3352  std::vector<
3353  std::tuple<unsigned int,
3354  unsigned int,
3355  std::string,
3358 };
3359 
3360 
3361 
3370 {
3371 public:
3375  XDMFEntry();
3376 
3382  XDMFEntry(const std::string & filename,
3383  const double time,
3384  const std::uint64_t nodes,
3385  const std::uint64_t cells,
3386  const unsigned int dim,
3387  const ReferenceCell &cell_type);
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 
3405  XDMFEntry(const std::string & mesh_filename,
3406  const std::string & solution_filename,
3407  const double time,
3408  const std::uint64_t nodes,
3409  const std::uint64_t cells,
3410  const unsigned int dim);
3411 
3416  XDMFEntry(const std::string & mesh_filename,
3417  const std::string & solution_filename,
3418  const double time,
3419  const std::uint64_t nodes,
3420  const std::uint64_t cells,
3421  const unsigned int dim,
3422  const ReferenceCell &cell_type);
3423 
3430  XDMFEntry(const std::string & mesh_filename,
3431  const std::string & solution_filename,
3432  const double time,
3433  const std::uint64_t nodes,
3434  const std::uint64_t cells,
3435  const unsigned int dim,
3436  const unsigned int spacedim);
3437 
3441  XDMFEntry(const std::string & mesh_filename,
3442  const std::string & solution_filename,
3443  const double time,
3444  const std::uint64_t nodes,
3445  const std::uint64_t cells,
3446  const unsigned int dim,
3447  const unsigned int spacedim,
3448  const ReferenceCell &cell_type);
3449 
3453  void
3454  add_attribute(const std::string &attr_name, const unsigned int dimension);
3455 
3461  template <class Archive>
3462  void
3463  serialize(Archive &ar, const unsigned int /*version*/)
3464  {
3467  }
3468 
3473  std::string
3474  get_xdmf_content(const unsigned int indent_level) const;
3475 
3483  std::string
3484  get_xdmf_content(const unsigned int indent_level,
3485  const ReferenceCell &reference_cell) const;
3486 
3487 private:
3491  bool valid;
3492 
3496  std::string h5_sol_filename;
3497 
3501  std::string h5_mesh_filename;
3502 
3506  double entry_time;
3507 
3511  std::uint64_t num_nodes;
3512 
3516  std::uint64_t num_cells;
3517 
3521  unsigned int dimension;
3522 
3527  unsigned int space_dimension;
3528 
3534 
3538  std::map<std::string, unsigned int> attribute_dims;
3539 };
3540 
3541 
3542 
3543 /* -------------------- inline functions ------------------- */
3544 
3545 namespace DataOutBase
3546 {
3547  inline bool
3549  {
3550  return (red == green) && (red == blue);
3551  }
3552 
3553 
3554  /* -------------------- template functions ------------------- */
3555 
3562  template <int dim, int spacedim>
3563  std::ostream &
3564  operator<<(std::ostream &out, const Patch<dim, spacedim> &patch);
3565 
3566 
3567 
3574  template <int dim, int spacedim>
3575  std::istream &
3576  operator>>(std::istream &in, Patch<dim, spacedim> &patch);
3577 } // namespace DataOutBase
3578 
3579 
3581 
3582 #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 unsigned int d1, const unsigned int d2, const unsigned int d3)
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
virtual std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > get_nonscalar_data_ranges() const
unsigned int default_subdivisions
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
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)
void write_deal_II_intermediate_in_parallel(const std::string &filename, const MPI_Comm &comm, const DataOutBase::CompressionLevel compression) const
void write_hdf5_parallel(const DataOutBase::DataOutFilter &data_filter, const std::string &filename, const MPI_Comm &comm) const
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
XDMFEntry create_xdmf_entry(const DataOutBase::DataOutFilter &data_filter, const std::string &h5_filename, const double cur_time, const MPI_Comm &comm) const
void write_vtu_in_parallel(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::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_tecplot(std::ostream &out) const
void write_svg(std::ostream &out) const
DataOutBase::SvgFlags svg_flags
void write_xdmf_file(const std::vector< XDMFEntry > &entries, const std::string &filename, const MPI_Comm &comm) const
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
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:164
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:457
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:458
unsigned int cell_index
Definition: grid_tools.cc:1145
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:578
static ::ExceptionBase & ExcInvalidDatasetSize(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidCombinationOfDimensions(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
static ::ExceptionBase & ExcNotEnoughSpaceDimensionLabels()
static ::ExceptionBase & ExcIncompatibleDatasetNames()
static ::ExceptionBase & ExcTecplotAPIError()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcNoPatches()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
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)
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_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_hdf5_parallel(const std::vector< Patch< dim, spacedim >> &patches, const DataOutFilter &data_filter, const std::string &filename, const MPI_Comm &comm)
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:206
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
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
Point< spacedim > vertices[GeometryInfo< dim >::vertices_per_cell]
bool operator==(const Patch &patch) const
static const unsigned int space_dim
unsigned int n_subdivisions
void swap(Patch< dim, spacedim > &other_patch)
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)
std::map< std::string, std::string > physical_units
static const DataOutBase::CompressionLevel default_compression
static const DataOutBase::CompressionLevel best_compression
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_compression, const bool write_higher_order_cells=false, const std::map< std::string, std::string > &physical_units={})
static const DataOutBase::CompressionLevel best_speed
DataOutBase::CompressionLevel compression_level
static const DataOutBase::CompressionLevel no_compression
const MPI_Comm & comm