40#ifdef DEAL_II_WITH_ZLIB
44#ifdef DEAL_II_WITH_HDF5
48#include <boost/iostreams/copy.hpp>
49#include <boost/iostreams/device/back_inserter.hpp>
50#include <boost/iostreams/filtering_stream.hpp>
51#ifdef DEAL_II_WITH_ZLIB
52# include <boost/iostreams/filter/zlib.hpp>
66 <<
"Unexpected input: expected line\n <" << arg1
67 <<
">\nbut got\n <" << arg2 <<
">");
73#ifdef DEAL_II_WITH_ZLIB
74 constexpr bool deal_ii_with_zlib =
true;
76 constexpr bool deal_ii_with_zlib =
false;
80#ifdef DEAL_II_WITH_ZLIB
91 return Z_NO_COMPRESSION;
95 return Z_BEST_COMPRESSION;
97 return Z_DEFAULT_COMPRESSION;
100 return Z_NO_COMPRESSION;
104# ifdef DEAL_II_WITH_MPI
115 return boost::iostreams::zlib::no_compression;
117 return boost::iostreams::zlib::best_speed;
119 return boost::iostreams::zlib::best_compression;
121 return boost::iostreams::zlib::default_compression;
124 return boost::iostreams::zlib::no_compression;
134 template <
typename T>
136 compress_array(
const std::vector<T> & data,
139#ifdef DEAL_II_WITH_ZLIB
140 if (data.size() != 0)
142 const std::size_t uncompressed_size = (data.size() *
sizeof(T));
150 std::numeric_limits<std::uint32_t>::max(),
154 auto compressed_data_length = compressBound(uncompressed_size);
156 std::numeric_limits<std::uint32_t>::max(),
159 std::vector<unsigned char> compressed_data(compressed_data_length);
161 int err = compress2(&compressed_data[0],
162 &compressed_data_length,
163 reinterpret_cast<const Bytef *
>(data.data()),
165 get_zlib_compression_level(compression_level));
170 compressed_data.resize(compressed_data_length);
173 const std::uint32_t compression_header[4] = {
175 static_cast<std::uint32_t
>(uncompressed_size),
176 static_cast<std::uint32_t
>(
178 static_cast<std::uint32_t
>(
179 compressed_data_length)};
181 const auto header_start =
182 reinterpret_cast<const unsigned char *
>(&compression_header[0]);
185 {header_start, header_start + 4 *
sizeof(std::uint32_t)}) +
192 (void)compression_level;
194 ExcMessage(
"This function can only be called if cmake found "
195 "a working libz installation."));
210 template <
typename T>
212 vtu_stringize_array(
const std::vector<T> & data,
216 if (deal_ii_with_zlib &&
220 return compress_array(data, compression_level);
224 std::ostringstream stream;
225 stream.precision(precision);
226 for (
const T &el : data)
241 struct ParallelIntermediateHeader
244 std::uint64_t version;
245 std::uint64_t compression;
246 std::uint64_t dimension;
247 std::uint64_t space_dimension;
248 std::uint64_t n_ranks;
249 std::uint64_t n_patches;
297 SvgCell::operator<(
const SvgCell &e)
const
338 EpsCell2d::operator<(
const EpsCell2d &e)
const
357 template <
int dim,
int spacedim,
typename Number =
double>
358 std::unique_ptr<Table<2, Number>>
359 create_global_data_table(
const std::vector<Patch<dim, spacedim>> &patches)
362 if (patches.size() == 0)
363 return std::make_unique<Table<2, Number>>();
372 const unsigned int n_data_sets = patches[0].points_are_available ?
373 (patches[0].data.n_rows() - spacedim) :
374 patches[0].data.n_rows();
375 const unsigned int n_data_points =
376 std::accumulate(patches.begin(),
379 [](
const unsigned int count,
381 return count + patch.data.n_cols();
384 std::unique_ptr<Table<2, Number>> global_data_table =
385 std::make_unique<Table<2, Number>>(n_data_sets, n_data_points);
388 unsigned int next_value = 0;
389 for (
const auto &patch : patches)
391 const unsigned int n_subdivisions = patch.n_subdivisions;
392 (void)n_subdivisions;
394 Assert((patch.data.n_rows() == n_data_sets &&
395 !patch.points_are_available) ||
396 (patch.data.n_rows() == n_data_sets + spacedim &&
397 patch.points_are_available),
399 (n_data_sets + spacedim) :
401 patch.data.n_rows()));
402 Assert(patch.reference_cell != ReferenceCells::get_hypercube<dim>() ||
403 (n_data_sets == 0) ||
404 (patch.data.n_cols() ==
405 Utilities::fixed_power<dim>(n_subdivisions + 1)),
407 n_subdivisions + 1));
409 for (
unsigned int i = 0; i < patch.data.n_cols(); ++i, ++next_value)
410 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
411 (*global_data_table)[data_set][next_value] =
412 patch.data(data_set, i);
416 return global_data_table;
424 , node_dim(
numbers::invalid_unsigned_int)
432 , node_dim(
numbers::invalid_unsigned_int)
445 for (
unsigned int d = 0; d < dim; ++d)
449 unsigned int internal_ind;
460 internal_ind = it->second;
470 const unsigned int pt_index)
488 for (
unsigned int d = 0; d <
node_dim; ++d)
489 node_data[
node_dim * existing_point.second + d] =
490 existing_point.first(d);
498 std::vector<unsigned int> &cell_data)
const
504 cell_data[filtered_cell.first] =
505 filtered_cell.second + local_node_offset;
574 const unsigned int start,
575 const std::array<unsigned int, dim> &offsets)
579 const unsigned int base_entry =
592 const unsigned int d1 = offsets[0];
601 const unsigned int d1 = offsets[0];
602 const unsigned int d2 = offsets[1];
613 const unsigned int d1 = offsets[0];
614 const unsigned int d2 = offsets[1];
615 const unsigned int d3 = offsets[2];
637 const unsigned int start,
638 const unsigned int n_points,
643 const unsigned int base_entry = index * n_points;
645 static const std::array<unsigned int, 5> table = {{0, 1, 3, 2, 4}};
647 for (
unsigned int i = 0; i < n_points; ++i)
658 const unsigned int dimension,
659 const unsigned int set_num,
662 unsigned int new_dim;
681 for (
unsigned int d = 0; d < new_dim; ++d)
684 data_sets.back()[r * new_dim + d] = data_vectors(set_num + d, i);
700 const char *gmv_cell_type[4] = {
"",
"line 2",
"quad 4",
"hex 8"};
702 const char *ucd_cell_type[4] = {
"pt",
"line",
"quad",
"hex"};
704 const char *tecplot_cell_type[4] = {
"",
"lineseg",
"quadrilateral",
"brick"};
719 template <
int dim,
int spacedim>
720 std::array<unsigned int, 3>
722 const bool write_higher_order_cells)
724 std::array<unsigned int, 3> vtk_cell_id = {
730 if (write_higher_order_cells)
733 vtk_cell_id[2] = patch.
data.n_cols();
739 patch.
data.n_cols() == 6)
743 vtk_cell_id[2] = patch.
data.n_cols();
746 patch.
data.n_cols() == 10)
750 vtk_cell_id[2] = patch.
data.n_cols();
775 template <
int dim,
int spacedim>
777 get_equispaced_location(
779 const std::initializer_list<unsigned int> &lattice_location,
780 const unsigned int n_subdivisions)
787 const unsigned int xstep = (dim > 0 ? *(lattice_location.begin() + 0) : 0);
788 const unsigned int ystep = (dim > 1 ? *(lattice_location.begin() + 1) : 0);
789 const unsigned int zstep = (dim > 2 ? *(lattice_location.begin() + 2) : 0);
797 unsigned int point_no = 0;
802 point_no += (n_subdivisions + 1) * (n_subdivisions + 1) * zstep;
806 point_no += (n_subdivisions + 1) * ystep;
820 for (
unsigned int d = 0;
d < spacedim; ++
d)
821 node[d] = patch.
data(patch.
data.size(0) - spacedim + d, point_no);
833 const double stepsize = 1. / n_subdivisions;
834 const double xfrac = xstep * stepsize;
840 const double yfrac = ystep * stepsize;
842 node += ((patch.
vertices[3] * xfrac) +
843 (patch.
vertices[2] * (1 - xfrac))) *
847 const double zfrac = zstep * stepsize;
849 node += (((patch.
vertices[5] * xfrac) +
850 (patch.
vertices[4] * (1 - xfrac))) *
853 (patch.
vertices[6] * (1 - xfrac))) *
866 template <
int dim,
int spacedim>
869 const unsigned int node_index)
874 unsigned int point_no_actual = node_index;
879 static const std::array<unsigned int, 5> table = {{0, 1, 3, 2, 4}};
880 point_no_actual = table[node_index];
888 for (
unsigned int d = 0;
d < spacedim; ++
d)
890 patch.
data(patch.
data.size(0) - spacedim + d, point_no_actual);
903 return patch.
vertices[point_no_actual];
914 template <
int dim,
int spacedim>
915 std::tuple<unsigned int, unsigned int>
916 count_nodes_and_cells(
919 unsigned int n_nodes = 0;
921 for (
const auto &patch : patches)
925 "The reference cell for this patch is set to 'Invalid', "
926 "but that is clearly not a valid choice. Did you forget "
927 "to set the reference cell for the patch?"));
942 return std::make_tuple(n_nodes, n_cells);
952 template <
int dim,
int spacedim>
953 std::tuple<unsigned int, unsigned int, unsigned int>
954 count_nodes_and_cells_and_points(
956 const bool write_higher_order_cells)
958 unsigned int n_nodes = 0;
960 unsigned int n_points_and_n_cells = 0;
962 for (
const auto &patch : patches)
968 if (write_higher_order_cells)
974 n_points_and_n_cells +=
983 const unsigned int n_subcells =
986 n_points_and_n_cells +=
992 n_nodes += patch.
data.n_cols();
994 n_points_and_n_cells += patch.
data.n_cols() + 1;
998 return std::make_tuple(n_nodes, n_cells, n_points_and_n_cells);
1006 template <
typename FlagsType>
1013 StreamBase(std::ostream &stream,
const FlagsType &flags)
1014 : selected_component(
numbers::invalid_unsigned_int)
1025 write_point(
const unsigned int,
const Point<dim> &)
1028 ExcMessage(
"The derived class you are using needs to "
1029 "reimplement this function if you want to call "
1049 write_cell(
const unsigned int ,
1050 const unsigned int ,
1051 std::array<unsigned int, dim> & )
1054 ExcMessage(
"The derived class you are using needs to "
1055 "reimplement this function if you want to call "
1066 write_cell_single(
const unsigned int index,
1067 const unsigned int start,
1068 const unsigned int n_points,
1074 (void)reference_cell;
1077 ExcMessage(
"The derived class you are using needs to "
1078 "reimplement this function if you want to call "
1096 template <
typename T>
1110 unsigned int selected_component;
1117 std::ostream &stream;
1122 const FlagsType flags;
1128 class DXStream :
public StreamBase<DataOutBase::DXFlags>
1135 write_point(
const unsigned int index,
const Point<dim> &);
1147 write_cell(
const unsigned int index,
1148 const unsigned int start,
1149 const std::array<unsigned int, dim> &offsets);
1157 template <
typename data>
1159 write_dataset(
const unsigned int index,
const std::vector<data> &values);
1165 class GmvStream :
public StreamBase<DataOutBase::GmvFlags>
1172 write_point(
const unsigned int index,
const Point<dim> &);
1184 write_cell(
const unsigned int index,
1185 const unsigned int start,
1186 const std::array<unsigned int, dim> &offsets);
1192 class TecplotStream :
public StreamBase<DataOutBase::TecplotFlags>
1199 write_point(
const unsigned int index,
const Point<dim> &);
1211 write_cell(
const unsigned int index,
1212 const unsigned int start,
1213 const std::array<unsigned int, dim> &offsets);
1219 class UcdStream :
public StreamBase<DataOutBase::UcdFlags>
1226 write_point(
const unsigned int index,
const Point<dim> &);
1240 write_cell(
const unsigned int index,
1241 const unsigned int start,
1242 const std::array<unsigned int, dim> &offsets);
1250 template <
typename data>
1252 write_dataset(
const unsigned int index,
const std::vector<data> &values);
1258 class VtkStream :
public StreamBase<DataOutBase::VtkFlags>
1265 write_point(
const unsigned int index,
const Point<dim> &);
1277 write_cell(
const unsigned int index,
1278 const unsigned int start,
1279 const std::array<unsigned int, dim> &offsets);
1285 write_cell_single(
const unsigned int index,
1286 const unsigned int start,
1287 const unsigned int n_points,
1299 write_high_order_cell(
const unsigned int start,
1300 const std::vector<unsigned> &connectivity);
1313 DXStream::write_point(
const unsigned int,
const Point<dim> &p)
1315 if (flags.coordinates_binary)
1318 for (
unsigned int d = 0;
d < dim; ++
d)
1320 stream.write(
reinterpret_cast<const char *
>(data), dim *
sizeof(*data));
1324 for (
unsigned int d = 0;
d < dim; ++
d)
1325 stream << p(d) <<
'\t';
1339 std::array<unsigned int, GeometryInfo<0>::vertices_per_cell>
1340 set_node_numbers(
const unsigned int ,
1341 const std::array<unsigned int, 0> & )
1349 std::array<unsigned int, GeometryInfo<1>::vertices_per_cell>
1350 set_node_numbers(
const unsigned int start,
1351 const std::array<unsigned int, 1> &offsets)
1353 std::array<unsigned int, GeometryInfo<1>::vertices_per_cell> nodes;
1355 nodes[1] = start + offsets[0];
1361 std::array<unsigned int, GeometryInfo<2>::vertices_per_cell>
1362 set_node_numbers(
const unsigned int start,
1363 const std::array<unsigned int, 2> &offsets)
1366 const unsigned int d1 = offsets[0];
1367 const unsigned int d2 = offsets[1];
1369 std::array<unsigned int, GeometryInfo<2>::vertices_per_cell> nodes;
1371 nodes[1] = start + d1;
1372 nodes[2] = start + d2;
1373 nodes[3] = start + d2 + d1;
1379 std::array<unsigned int, GeometryInfo<3>::vertices_per_cell>
1380 set_node_numbers(
const unsigned int start,
1381 const std::array<unsigned int, 3> &offsets)
1383 const unsigned int d1 = offsets[0];
1384 const unsigned int d2 = offsets[1];
1385 const unsigned int d3 = offsets[2];
1387 std::array<unsigned int, GeometryInfo<3>::vertices_per_cell> nodes;
1389 nodes[1] = start + d1;
1390 nodes[2] = start + d2;
1391 nodes[3] = start + d2 + d1;
1392 nodes[4] = start + d3;
1393 nodes[5] = start + d3 + d1;
1394 nodes[6] = start + d3 + d2;
1395 nodes[7] = start + d3 + d2 + d1;
1404 DXStream::write_cell(
const unsigned int,
1405 const unsigned int start,
1406 const std::array<unsigned int, dim> &offsets)
1409 DataOutBaseImplementation::set_node_numbers(start, offsets);
1411 if (flags.int_binary)
1413 std::array<unsigned int, GeometryInfo<dim>::vertices_per_cell> temp;
1414 for (
unsigned int i = 0; i < nodes.size(); ++i)
1416 stream.write(
reinterpret_cast<const char *
>(temp.data()),
1417 temp.size() *
sizeof(temp[0]));
1421 for (
unsigned int i = 0; i < nodes.size() - 1; ++i)
1423 stream << nodes[GeometryInfo<dim>::dx_to_deal[nodes.size() - 1]]
1430 template <
typename data>
1432 DXStream::write_dataset(
const unsigned int,
const std::vector<data> &values)
1434 if (flags.data_binary)
1436 stream.write(
reinterpret_cast<const char *
>(
values.data()),
1437 values.size() *
sizeof(data));
1441 for (
unsigned int i = 0; i <
values.size(); ++i)
1442 stream <<
'\t' << values[i];
1458 GmvStream::write_point(
const unsigned int,
const Point<dim> &p)
1462 stream << p(selected_component) <<
' ';
1469 GmvStream::write_cell(
const unsigned int,
1470 const unsigned int s,
1471 const std::array<unsigned int, dim> &offsets)
1474 const unsigned int start = s + 1;
1475 stream << gmv_cell_type[dim] <<
'\n';
1487 const unsigned int d1 = offsets[0];
1489 stream <<
'\t' << start + d1;
1495 const unsigned int d1 = offsets[0];
1496 const unsigned int d2 = offsets[1];
1498 stream <<
'\t' << start + d1;
1499 stream <<
'\t' << start + d2 + d1 <<
'\t' << start + d2;
1505 const unsigned int d1 = offsets[0];
1506 const unsigned int d2 = offsets[1];
1507 const unsigned int d3 = offsets[2];
1509 stream <<
'\t' << start + d1;
1510 stream <<
'\t' << start + d2 + d1 <<
'\t' << start + d2;
1511 stream <<
'\t' << start + d3 <<
'\t' << start + d3 + d1 <<
'\t'
1512 << start + d3 + d2 + d1 <<
'\t' << start + d3 + d2;
1524 TecplotStream::TecplotStream(std::ostream & out,
1532 TecplotStream::write_point(
const unsigned int,
const Point<dim> &p)
1536 stream << p(selected_component) <<
'\n';
1543 TecplotStream::write_cell(
const unsigned int,
1544 const unsigned int s,
1545 const std::array<unsigned int, dim> &offsets)
1547 const unsigned int start = s + 1;
1559 const unsigned int d1 = offsets[0];
1561 stream <<
'\t' << start + d1;
1567 const unsigned int d1 = offsets[0];
1568 const unsigned int d2 = offsets[1];
1570 stream <<
'\t' << start + d1;
1571 stream <<
'\t' << start + d2 + d1 <<
'\t' << start + d2;
1577 const unsigned int d1 = offsets[0];
1578 const unsigned int d2 = offsets[1];
1579 const unsigned int d3 = offsets[2];
1581 stream <<
'\t' << start + d1;
1582 stream <<
'\t' << start + d2 + d1 <<
'\t' << start + d2;
1583 stream <<
'\t' << start + d3 <<
'\t' << start + d3 + d1 <<
'\t'
1584 << start + d3 + d2 + d1 <<
'\t' << start + d3 + d2;
1603 UcdStream::write_point(
const unsigned int index,
const Point<dim> &p)
1605 stream <<
index + 1 <<
" ";
1607 for (
unsigned int i = 0; i < dim; ++i)
1608 stream << p(i) <<
' ';
1610 for (
unsigned int i = dim; i < 3; ++i)
1619 UcdStream::write_cell(
const unsigned int index,
1620 const unsigned int start,
1621 const std::array<unsigned int, dim> &offsets)
1624 DataOutBaseImplementation::set_node_numbers(start, offsets);
1627 stream <<
index + 1 <<
"\t0 " << ucd_cell_type[dim];
1628 for (
unsigned int i = 0; i < nodes.size(); ++i)
1635 template <
typename data>
1637 UcdStream::write_dataset(
const unsigned int index,
1638 const std::vector<data> &values)
1640 stream <<
index + 1;
1641 for (
unsigned int i = 0; i <
values.size(); ++i)
1642 stream <<
'\t' << values[i];
1657 VtkStream::write_point(
const unsigned int,
const Point<dim> &p)
1662 for (
unsigned int i = dim; i < 3; ++i)
1671 VtkStream::write_cell(
const unsigned int,
1672 const unsigned int start,
1673 const std::array<unsigned int, dim> &offsets)
1675 stream << GeometryInfo<dim>::vertices_per_cell <<
'\t';
1687 const unsigned int d1 = offsets[0];
1689 stream <<
'\t' << start + d1;
1695 const unsigned int d1 = offsets[0];
1696 const unsigned int d2 = offsets[1];
1698 stream <<
'\t' << start + d1;
1699 stream <<
'\t' << start + d2 + d1 <<
'\t' << start + d2;
1705 const unsigned int d1 = offsets[0];
1706 const unsigned int d2 = offsets[1];
1707 const unsigned int d3 = offsets[2];
1709 stream <<
'\t' << start + d1;
1710 stream <<
'\t' << start + d2 + d1 <<
'\t' << start + d2;
1711 stream <<
'\t' << start + d3 <<
'\t' << start + d3 + d1 <<
'\t'
1712 << start + d3 + d2 + d1 <<
'\t' << start + d3 + d2;
1725 VtkStream::write_cell_single(
const unsigned int index,
1726 const unsigned int start,
1727 const unsigned int n_points,
1732 static const std::array<unsigned int, 5> table = {{0, 1, 3, 2, 4}};
1734 stream <<
'\t' << n_points;
1735 for (
unsigned int i = 0; i < n_points; ++i)
1744 VtkStream::write_high_order_cell(
const unsigned int start,
1745 const std::vector<unsigned> &connectivity)
1747 stream << connectivity.size();
1748 for (
const auto &c : connectivity)
1749 stream <<
'\t' << start + c;
1761 template <
int dim,
int spacedim>
1765 template <
int dim,
int spacedim>
1769 template <
int dim,
int spacedim>
1771 : patch_index(no_neighbor)
1773 , points_are_available(false)
1787 template <
int dim,
int spacedim>
1795 const double epsilon = 3e-16;
1813 if (data.n_rows() != patch.
data.n_rows())
1816 if (data.n_cols() != patch.
data.n_cols())
1819 for (
unsigned int i = 0; i < data.n_rows(); ++i)
1820 for (
unsigned int j = 0; j < data.n_cols(); ++j)
1821 if (data[i][j] != patch.
data[i][j])
1829 template <
int dim,
int spacedim>
1835 sizeof(neighbors) /
sizeof(neighbors[0]) *
1841 sizeof(reference_cell));
1846 template <
int dim,
int spacedim>
1851 std::swap(neighbors, other_patch.
neighbors);
1854 data.swap(other_patch.
data);
1861 template <
int spacedim>
1865 template <
int spacedim>
1869 template <
int spacedim>
1873 template <
int spacedim>
1876 template <
int spacedim>
1880 template <
int spacedim>
1882 : patch_index(no_neighbor)
1883 , points_are_available(false)
1890 template <
int spacedim>
1894 const unsigned int dim = 0;
1897 const double epsilon = 3e-16;
1908 if (data.n_rows() != patch.
data.n_rows())
1911 if (data.n_cols() != patch.
data.n_cols())
1914 for (
unsigned int i = 0; i < data.n_rows(); ++i)
1915 for (
unsigned int j = 0; j < data.n_cols(); ++j)
1916 if (data[i][j] != patch.
data[i][j])
1924 template <
int spacedim>
1936 template <
int spacedim>
1942 data.swap(other_patch.
data);
1949 : write_preamble(write_preamble)
1964 : space_dimension_labels(labels)
1978 const bool bicubic_patch,
1979 const bool external_data)
1981 , bicubic_patch(bicubic_patch)
1982 , external_data(external_data)
1987 const bool xdmf_hdf5_output)
1988 : filter_duplicate_vertices(filter_duplicate_vertices)
1989 , xdmf_hdf5_output(xdmf_hdf5_output)
1997 "Filter duplicate vertices",
2000 "Whether to remove duplicate vertex values. deal.II duplicates "
2001 "vertices once for each adjacent cell so that it can output "
2002 "discontinuous quantities for which there may be more than one "
2003 "value for each vertex position. Setting this flag to "
2004 "'true' will merge all of these values by selecting a "
2005 "random one and outputting this as 'the' value for the vertex. "
2006 "As long as the data to be output corresponds to continuous "
2007 "fields, merging vertices has no effect. On the other hand, "
2008 "if the data to be output corresponds to discontinuous fields "
2009 "(either because you are using a discontinuous finite element, "
2010 "or because you are using a DataPostprocessor that yields "
2011 "discontinuous data, or because the data to be output has been "
2012 "produced by entirely different means), then the data in the "
2013 "output file no longer faithfully represents the underlying data "
2014 "because the discontinuous field has been replaced by a "
2015 "continuous one. Note also that the filtering can not occur "
2016 "on processor boundaries. Thus, a filtered discontinuous field "
2017 "looks like a continuous field inside of a subdomain, "
2018 "but like a discontinuous field at the subdomain boundary."
2020 "In any case, filtering results in drastically smaller output "
2021 "files (smaller by about a factor of 2^dim).");
2026 "Whether the data will be used in an XDMF/HDF5 combination.");
2041 const bool int_binary,
2042 const bool coordinates_binary,
2043 const bool data_binary)
2044 : write_neighbors(write_neighbors)
2045 , int_binary(int_binary)
2046 , coordinates_binary(coordinates_binary)
2047 , data_binary(data_binary)
2048 , data_double(false)
2058 "A boolean field indicating whether neighborship "
2059 "information between cells is to be written to the "
2060 "OpenDX output file");
2064 "Output format of integer numbers, which is "
2065 "either a text representation (ascii) or binary integer "
2066 "values of 32 or 64 bits length");
2070 "Output format of vertex coordinates, which is "
2071 "either a text representation (ascii) or binary "
2072 "floating point values of 32 or 64 bits length");
2076 "Output format of data values, which is "
2077 "either a text representation (ascii) or binary "
2078 "floating point values of 32 or 64 bits length");
2098 "A flag indicating whether a comment should be "
2099 "written to the beginning of the output file "
2100 "indicating date and time of creation as well "
2101 "as the creating program");
2115 const int azimuth_angle,
2116 const int polar_angle,
2117 const unsigned int line_thickness,
2119 const bool draw_colorbar)
2122 , height_vector(height_vector)
2123 , azimuth_angle(azimuth_angle)
2124 , polar_angle(polar_angle)
2125 , line_thickness(line_thickness)
2127 , draw_colorbar(draw_colorbar)
2138 "A flag indicating whether POVRAY should use smoothed "
2139 "triangles instead of the usual ones");
2143 "Whether POVRAY should use bicubic patches");
2147 "Whether camera and lighting information should "
2148 "be put into an external file \"data.inc\" or into "
2149 "the POVRAY input file");
2165 const unsigned int color_vector,
2167 const unsigned int size,
2168 const double line_width,
2169 const double azimut_angle,
2170 const double turn_angle,
2171 const double z_scaling,
2172 const bool draw_mesh,
2173 const bool draw_cells,
2174 const bool shade_cells,
2176 : height_vector(height_vector)
2177 , color_vector(color_vector)
2178 , size_type(size_type)
2180 , line_width(line_width)
2181 , azimut_angle(azimut_angle)
2182 , turn_angle(turn_angle)
2183 , z_scaling(z_scaling)
2184 , draw_mesh(draw_mesh)
2185 , draw_cells(draw_cells)
2186 , shade_cells(shade_cells)
2187 , color_function(color_function)
2226 double sum = xmax + xmin;
2227 double sum13 = xmin + 3 * xmax;
2228 double sum22 = 2 * xmin + 2 * xmax;
2229 double sum31 = 3 * xmin + xmax;
2230 double dif = xmax - xmin;
2231 double rezdif = 1.0 / dif;
2235 if (x < (sum31) / 4)
2237 else if (x < (sum22) / 4)
2239 else if (x < (sum13) / 4)
2250 rgb_values.
green = 0;
2251 rgb_values.
blue = (x - xmin) * 4. * rezdif;
2255 rgb_values.
green = (4 * x - 3 * xmin - xmax) * rezdif;
2256 rgb_values.
blue = (sum22 - 4. * x) * rezdif;
2259 rgb_values.
red = (4 * x - 2 * sum) * rezdif;
2260 rgb_values.
green = (xmin + 3 * xmax - 4 * x) * rezdif;
2261 rgb_values.
blue = 0;
2265 rgb_values.
green = (4 * x - xmin - 3 * xmax) * rezdif;
2266 rgb_values.
blue = (4. * x - sum13) * rezdif;
2273 rgb_values.
red = rgb_values.
green = rgb_values.
blue = 1;
2287 (x - xmin) / (xmax - xmin);
2300 1 - (x - xmin) / (xmax - xmin);
2312 "Number of the input vector that is to be used to "
2313 "generate height information");
2317 "Number of the input vector that is to be used to "
2318 "generate color information");
2322 "Whether width or height should be scaled to match "
2327 "The size (width or height) to which the eps output "
2328 "file is to be scaled");
2332 "The width in which the postscript renderer is to "
2337 "Angle of the viewing position against the vertical "
2342 "Angle of the viewing direction against the y-axis");
2346 "Scaling for the z-direction relative to the scaling "
2347 "used in x- and y-directions");
2351 "Whether the mesh lines, or only the surface should be "
2356 "Whether only the mesh lines, or also the interior of "
2357 "cells should be plotted. If this flag is false, then "
2358 "one can see through the mesh");
2362 "Whether the interior of cells shall be shaded");
2366 "default|grey scale|reverse grey scale"),
2367 "Name of a color function used to colorize mesh lines "
2368 "and/or cell interiors");
2378 if (prm.
get(
"Scale to width or height") ==
"width")
2390 if (prm.
get(
"Color function") ==
"default")
2392 else if (prm.
get(
"Color function") ==
"grey scale")
2394 else if (prm.
get(
"Color function") ==
"reverse grey scale")
2404 : compression_level(compression_level)
2409 : zone_name(zone_name)
2410 , solution_time(solution_time)
2424 const unsigned int cycle,
2425 const bool print_date_and_time,
2427 const bool write_higher_order_cells,
2428 const std::map<std::string, std::string> &physical_units)
2431 , print_date_and_time(print_date_and_time)
2432 , compression_level(compression_level)
2433 , write_higher_order_cells(write_higher_order_cells)
2434 , physical_units(physical_units)
2442 if (format_name ==
"none")
2445 if (format_name ==
"dx")
2448 if (format_name ==
"ucd")
2451 if (format_name ==
"gnuplot")
2454 if (format_name ==
"povray")
2457 if (format_name ==
"eps")
2460 if (format_name ==
"gmv")
2463 if (format_name ==
"tecplot")
2466 if (format_name ==
"vtk")
2469 if (format_name ==
"vtu")
2472 if (format_name ==
"deal.II intermediate")
2475 if (format_name ==
"hdf5")
2479 ExcMessage(
"The given file format name is not recognized: <" +
2480 format_name +
">"));
2491 return "none|dx|ucd|gnuplot|povray|eps|gmv|tecplot|vtk|vtu|hdf5|svg|deal.II intermediate";
2499 switch (output_format)
2541 template <
int dim,
int spacedim>
2542 std::vector<Point<spacedim>>
2546 static const std::array<unsigned int, 5> table = {{0, 1, 3, 2, 4}};
2548 std::vector<Point<spacedim>> node_positions;
2549 for (
const auto &patch : patches)
2552 if (patch.
reference_cell != ReferenceCells::get_hypercube<dim>())
2554 for (
unsigned int point_no = 0; point_no < patch.
data.n_cols();
2556 node_positions.emplace_back(get_node_location(
2565 const unsigned int n = n_subdivisions + 1;
2570 node_positions.emplace_back(
2571 get_equispaced_location(patch, {}, n_subdivisions));
2574 for (
unsigned int i1 = 0; i1 < n; ++i1)
2575 node_positions.emplace_back(
2576 get_equispaced_location(patch, {i1}, n_subdivisions));
2579 for (
unsigned int i2 = 0; i2 < n; ++i2)
2580 for (
unsigned int i1 = 0; i1 < n; ++i1)
2581 node_positions.emplace_back(get_equispaced_location(
2582 patch, {i1, i2}, n_subdivisions));
2585 for (
unsigned int i3 = 0; i3 < n; ++i3)
2586 for (
unsigned int i2 = 0; i2 < n; ++i2)
2587 for (
unsigned int i1 = 0; i1 < n; ++i1)
2588 node_positions.emplace_back(get_equispaced_location(
2589 patch, {i1, i2, i3}, n_subdivisions));
2598 return node_positions;
2602 template <
int dim,
int spacedim,
typename StreamType>
2608 const std::vector<Point<spacedim>> node_positions =
2612 for (
const auto &node : node_positions)
2613 out.write_point(count++, node);
2619 template <
int dim,
int spacedim,
typename StreamType>
2624 unsigned int count = 0;
2625 unsigned int first_vertex_of_patch = 0;
2626 for (
const auto &patch : patches)
2629 if (patch.
reference_cell != ReferenceCells::get_hypercube<dim>())
2631 out.write_cell_single(count++,
2632 first_vertex_of_patch,
2633 patch.
data.n_cols(),
2635 first_vertex_of_patch += patch.
data.n_cols();
2640 const unsigned int n = n_subdivisions + 1;
2646 const unsigned int offset = first_vertex_of_patch;
2647 out.template write_cell<0>(count++, offset, {});
2653 constexpr unsigned int d1 = 1;
2655 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
2657 const unsigned int offset =
2658 first_vertex_of_patch + i1 * d1;
2659 out.template write_cell<1>(count++, offset, {{d1}});
2667 constexpr unsigned int d1 = 1;
2668 const unsigned int d2 = n;
2670 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
2671 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
2673 const unsigned int offset =
2674 first_vertex_of_patch + i2 * d2 + i1 * d1;
2675 out.template write_cell<2>(count++,
2685 constexpr unsigned int d1 = 1;
2686 const unsigned int d2 = n;
2687 const unsigned int d3 = n * n;
2689 for (
unsigned int i3 = 0; i3 < n_subdivisions; ++i3)
2690 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
2691 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
2693 const unsigned int offset = first_vertex_of_patch +
2696 out.template write_cell<3>(count++,
2708 first_vertex_of_patch +=
2709 Utilities::fixed_power<dim>(n_subdivisions + 1);
2718 template <
int dim,
int spacedim,
typename StreamType>
2722 const bool legacy_format)
2725 unsigned int first_vertex_of_patch = 0;
2727 std::vector<unsigned> connectivity;
2729 for (
const auto &patch : patches)
2731 if (patch.
reference_cell != ReferenceCells::get_hypercube<dim>())
2733 connectivity.resize(patch.
data.n_cols());
2735 for (
unsigned int i = 0; i < patch.
data.n_cols(); ++i)
2736 connectivity[i] = i;
2738 out.template write_high_order_cell<dim>(first_vertex_of_patch,
2741 first_vertex_of_patch += patch.
data.n_cols();
2746 const unsigned int n = n_subdivisions + 1;
2748 connectivity.resize(Utilities::fixed_power<dim>(n));
2755 ExcMessage(
"Point-like cells should not be possible "
2756 "when writing higher-order cells."));
2761 for (
unsigned int i1 = 0; i1 < n_subdivisions + 1; ++i1)
2763 const unsigned int local_index = i1;
2764 const unsigned int connectivity_index =
2766 .template vtk_lexicographic_to_node_index<1>(
2767 {{i1}}, {{n_subdivisions}}, legacy_format);
2768 connectivity[connectivity_index] = local_index;
2775 for (
unsigned int i2 = 0; i2 < n_subdivisions + 1; ++i2)
2776 for (
unsigned int i1 = 0; i1 < n_subdivisions + 1; ++i1)
2778 const unsigned int local_index = i2 * n + i1;
2779 const unsigned int connectivity_index =
2781 .template vtk_lexicographic_to_node_index<2>(
2783 {{n_subdivisions, n_subdivisions}},
2785 connectivity[connectivity_index] = local_index;
2792 for (
unsigned int i3 = 0; i3 < n_subdivisions + 1; ++i3)
2793 for (
unsigned int i2 = 0; i2 < n_subdivisions + 1; ++i2)
2794 for (
unsigned int i1 = 0; i1 < n_subdivisions + 1; ++i1)
2796 const unsigned int local_index =
2797 i3 * n * n + i2 * n + i1;
2798 const unsigned int connectivity_index =
2800 .template vtk_lexicographic_to_node_index<3>(
2806 connectivity[connectivity_index] = local_index;
2817 out.template write_high_order_cell<dim>(first_vertex_of_patch,
2821 first_vertex_of_patch += Utilities::fixed_power<dim>(n);
2829 template <
int dim,
int spacedim,
class StreamType>
2832 unsigned int n_data_sets,
2833 const bool double_precision,
2837 unsigned int count = 0;
2839 for (
const auto &patch : patches)
2842 const unsigned int n = n_subdivisions + 1;
2846 (patch.
data.n_rows() == n_data_sets + spacedim &&
2849 (n_data_sets + spacedim) :
2851 patch.data.n_rows()));
2852 Assert(patch.
data.n_cols() == Utilities::fixed_power<dim>(n),
2855 std::vector<float> floats(n_data_sets);
2856 std::vector<double> doubles(n_data_sets);
2859 for (
unsigned int i = 0; i < Utilities::fixed_power<dim>(n);
2861 if (double_precision)
2863 for (
unsigned int data_set = 0; data_set < n_data_sets;
2865 doubles[data_set] = patch.
data(data_set, i);
2866 out.write_dataset(count, doubles);
2870 for (
unsigned int data_set = 0; data_set < n_data_sets;
2872 floats[data_set] = patch.
data(data_set, i);
2873 out.write_dataset(count, floats);
2898 camera_vertical[0] = camera_horizontal[1] * camera_direction[2] -
2899 camera_horizontal[2] * camera_direction[1];
2900 camera_vertical[1] = camera_horizontal[2] * camera_direction[0] -
2901 camera_horizontal[0] * camera_direction[2];
2902 camera_vertical[2] = camera_horizontal[0] * camera_direction[1] -
2903 camera_horizontal[1] * camera_direction[0];
2907 phi /= (point[0] - camera_position[0]) * camera_direction[0] +
2908 (point[1] - camera_position[1]) * camera_direction[1] +
2909 (point[2] - camera_position[2]) * camera_direction[2];
2913 camera_position[0] + phi * (point[0] - camera_position[0]);
2915 camera_position[1] + phi * (point[1] - camera_position[1]);
2917 camera_position[2] + phi * (point[2] - camera_position[2]);
2920 projection_decomposition[0] = (projection[0] - camera_position[0] -
2921 camera_focus * camera_direction[0]) *
2922 camera_horizontal[0];
2923 projection_decomposition[0] += (projection[1] - camera_position[1] -
2924 camera_focus * camera_direction[1]) *
2925 camera_horizontal[1];
2926 projection_decomposition[0] += (projection[2] - camera_position[2] -
2927 camera_focus * camera_direction[2]) *
2928 camera_horizontal[2];
2930 projection_decomposition[1] = (projection[0] - camera_position[0] -
2931 camera_focus * camera_direction[0]) *
2933 projection_decomposition[1] += (projection[1] - camera_position[1] -
2934 camera_focus * camera_direction[1]) *
2936 projection_decomposition[1] += (projection[2] - camera_position[2] -
2937 camera_focus * camera_direction[2]) *
2940 return projection_decomposition;
2949 svg_get_gradient_parameters(
Point<3> points[])
2955 for (
int i = 0; i < 2; ++i)
2957 for (
int j = 0; j < 2 - i; ++j)
2959 if (points[j][2] > points[j + 1][2])
2962 points[j] = points[j + 1];
2963 points[j + 1] = temp;
2970 v_inter = points[1];
2977 A[0][0] = v_max[0] - v_min[0];
2978 A[0][1] = v_inter[0] - v_min[0];
2979 A[1][0] = v_max[1] - v_min[1];
2980 A[1][1] = v_inter[1] - v_min[1];
2986 bool col_change =
false;
2995 double temp = A[1][0];
3000 for (
unsigned int k = 0; k < 1; ++k)
3002 for (
unsigned int i = k + 1; i < 2; ++i)
3004 x = A[i][k] / A[k][k];
3006 for (
unsigned int j = k + 1; j < 2; ++j)
3007 A[i][j] = A[i][j] - A[k][j] * x;
3009 b[i] = b[i] - b[k] * x;
3013 b[1] =
b[1] / A[1][1];
3015 for (
int i = 0; i >= 0; i--)
3019 for (
unsigned int j = i + 1; j < 2; ++j)
3020 sum = sum - A[i][j] * b[j];
3022 b[i] = sum / A[i][i];
3032 double c =
b[0] * (v_max[2] - v_min[2]) + b[1] * (v_inter[2] - v_min[2]) +
3036 A[0][0] = v_max[0] - v_min[0];
3037 A[0][1] = v_inter[0] - v_min[0];
3038 A[1][0] = v_max[1] - v_min[1];
3039 A[1][1] = v_inter[1] - v_min[1];
3041 b[0] = 1.0 - v_min[0];
3053 double temp = A[1][0];
3058 for (
unsigned int k = 0; k < 1; ++k)
3060 for (
unsigned int i = k + 1; i < 2; ++i)
3062 x = A[i][k] / A[k][k];
3064 for (
unsigned int j = k + 1; j < 2; ++j)
3065 A[i][j] = A[i][j] - A[k][j] * x;
3067 b[i] = b[i] - b[k] * x;
3071 b[1] =
b[1] / A[1][1];
3073 for (
int i = 0; i >= 0; i--)
3077 for (
unsigned int j = i + 1; j < 2; ++j)
3078 sum = sum - A[i][j] * b[j];
3080 b[i] =
sum / A[i][i];
3090 gradient[0] =
b[0] * (v_max[2] - v_min[2]) +
3091 b[1] * (v_inter[2] - v_min[2]) - c + v_min[2];
3094 A[0][0] = v_max[0] - v_min[0];
3095 A[0][1] = v_inter[0] - v_min[0];
3096 A[1][0] = v_max[1] - v_min[1];
3097 A[1][1] = v_inter[1] - v_min[1];
3100 b[1] = 1.0 - v_min[1];
3111 double temp = A[1][0];
3116 for (
unsigned int k = 0; k < 1; ++k)
3118 for (
unsigned int i = k + 1; i < 2; ++i)
3120 x = A[i][k] / A[k][k];
3122 for (
unsigned int j = k + 1; j < 2; ++j)
3123 A[i][j] = A[i][j] - A[k][j] * x;
3125 b[i] = b[i] - b[k] * x;
3129 b[1] =
b[1] / A[1][1];
3131 for (
int i = 0; i >= 0; i--)
3135 for (
unsigned int j = i + 1; j < 2; ++j)
3136 sum = sum - A[i][j] * b[j];
3138 b[i] =
sum / A[i][i];
3148 gradient[1] =
b[0] * (v_max[2] - v_min[2]) +
3149 b[1] * (v_inter[2] - v_min[2]) - c + v_min[2];
3155 gradient[1] * (v_min[1] - v_max[1]);
3159 gradient_parameters[0] = v_min[0];
3160 gradient_parameters[1] = v_min[1];
3165 gradient_parameters[4] = v_min[2];
3166 gradient_parameters[5] = v_max[2];
3168 return gradient_parameters;
3174 template <
int dim,
int spacedim>
3178 const std::vector<std::string> & data_names,
3180 std::tuple<
unsigned int,
3193#ifndef DEAL_II_WITH_MPI
3202 if (patches.size() == 0)
3206 const unsigned int n_data_sets = data_names.size();
3208 UcdStream ucd_out(out, flags);
3211 unsigned int n_nodes;
3212 unsigned int n_cells;
3213 std::tie(n_nodes, n_cells) = count_nodes_and_cells(patches);
3219 <<
"# This file was generated by the deal.II library." <<
'\n'
3223 <<
"# For a description of the UCD format see the AVS Developer's guide."
3229 out << n_nodes <<
' ' << n_cells <<
' ' << n_data_sets <<
' ' << 0
3242 if (n_data_sets != 0)
3244 out << n_data_sets <<
" ";
3245 for (
unsigned int i = 0; i < n_data_sets; ++i)
3250 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
3251 out << data_names[data_set]
3255 write_data(patches, n_data_sets,
true, ucd_out);
3265 template <
int dim,
int spacedim>
3269 const std::vector<std::string> & data_names,
3271 std::tuple<
unsigned int,
3283#ifndef DEAL_II_WITH_MPI
3292 if (patches.size() == 0)
3296 DXStream dx_out(out, flags);
3299 unsigned int offset = 0;
3301 const unsigned int n_data_sets = data_names.size();
3304 unsigned int n_nodes;
3305 unsigned int n_cells;
3306 std::tie(n_nodes, n_cells) = count_nodes_and_cells(patches);
3309 out <<
"object \"vertices\" class array type float rank 1 shape "
3310 << spacedim <<
" items " << n_nodes;
3314 out <<
" lsb ieee data 0" <<
'\n';
3315 offset += n_nodes * spacedim *
sizeof(float);
3319 out <<
" data follows" <<
'\n';
3328 out <<
"object \"cells\" class array type int rank 1 shape "
3333 out <<
" lsb binary data " << offset <<
'\n';
3334 offset += n_cells *
sizeof(
int);
3338 out <<
" data follows" <<
'\n';
3344 out <<
"attribute \"element type\" string \"";
3351 out <<
"\"" <<
'\n' <<
"attribute \"ref\" string \"positions\"" <<
'\n';
3358 out <<
"object \"neighbors\" class array type int rank 1 shape "
3362 for (
const auto &patch : patches)
3365 const unsigned int n1 = (dim > 0) ? n : 1;
3366 const unsigned int n2 = (dim > 1) ? n : 1;
3367 const unsigned int n3 = (dim > 2) ? n : 1;
3368 const unsigned int x_minus = (dim > 0) ? 0 : 0;
3369 const unsigned int x_plus = (dim > 0) ? 1 : 0;
3370 const unsigned int y_minus = (dim > 1) ? 2 : 0;
3371 const unsigned int y_plus = (dim > 1) ? 3 : 0;
3372 const unsigned int z_minus = (dim > 2) ? 4 : 0;
3373 const unsigned int z_plus = (dim > 2) ? 5 : 0;
3374 unsigned int cells_per_patch = Utilities::fixed_power<dim>(n);
3375 unsigned int dx = 1;
3376 unsigned int dy = n;
3377 unsigned int dz = n * n;
3379 const unsigned int patch_start =
3382 for (
unsigned int i3 = 0; i3 < n3; ++i3)
3383 for (
unsigned int i2 = 0; i2 < n2; ++i2)
3384 for (
unsigned int i1 = 0; i1 < n1; ++i1)
3386 const unsigned int nx = i1 *
dx;
3387 const unsigned int ny = i2 * dy;
3388 const unsigned int nz = i3 * dz;
3400 const unsigned int nn = patch.
neighbors[x_minus];
3404 << (nn * cells_per_patch + ny + nz +
dx * (n - 1));
3410 out <<
'\t' << patch_start + nx -
dx + ny + nz;
3415 const unsigned int nn = patch.
neighbors[x_plus];
3418 out << (nn * cells_per_patch + ny + nz);
3424 out <<
'\t' << patch_start + nx +
dx + ny + nz;
3431 const unsigned int nn = patch.
neighbors[y_minus];
3435 << (nn * cells_per_patch + nx + nz + dy * (n - 1));
3441 out <<
'\t' << patch_start + nx + ny - dy + nz;
3446 const unsigned int nn = patch.
neighbors[y_plus];
3449 out << (nn * cells_per_patch + nx + nz);
3455 out <<
'\t' << patch_start + nx + ny + dy + nz;
3463 const unsigned int nn = patch.
neighbors[z_minus];
3467 << (nn * cells_per_patch + nx + ny + dz * (n - 1));
3473 out <<
'\t' << patch_start + nx + ny + nz - dz;
3478 const unsigned int nn = patch.
neighbors[z_plus];
3481 out << (nn * cells_per_patch + nx + ny);
3487 out <<
'\t' << patch_start + nx + ny + nz + dz;
3495 if (n_data_sets != 0)
3497 out <<
"object \"data\" class array type float rank 1 shape "
3498 << n_data_sets <<
" items " << n_nodes;
3502 out <<
" lsb ieee data " << offset <<
'\n';
3503 offset += n_data_sets * n_nodes *
3504 ((flags.
data_double) ?
sizeof(
double) :
sizeof(float));
3508 out <<
" data follows" <<
'\n';
3513 out <<
"attribute \"dep\" string \"positions\"" <<
'\n';
3517 out <<
"object \"data\" class constantarray type float rank 0 items "
3518 << n_nodes <<
" data follows" <<
'\n'
3524 out <<
"object \"deal data\" class field" <<
'\n'
3525 <<
"component \"positions\" value \"vertices\"" <<
'\n'
3526 <<
"component \"connections\" value \"cells\"" <<
'\n'
3527 <<
"component \"data\" value \"data\"" <<
'\n';
3530 out <<
"component \"neighbors\" value \"neighbors\"" <<
'\n';
3537 out <<
"end" <<
'\n';
3555 template <
int dim,
int spacedim>
3559 const std::vector<std::string> & data_names,
3561 std::tuple<
unsigned int,
3570#ifndef DEAL_II_WITH_MPI
3580 if (patches.size() == 0)
3584 const unsigned int n_data_sets = data_names.size();
3588 out <<
"# This file was generated by the deal.II library." <<
'\n'
3592 <<
"# For a description of the GNUPLOT format see the GNUPLOT manual."
3599 for (
unsigned int spacedim_n = 0; spacedim_n < spacedim; ++spacedim_n)
3604 for (
const auto &data_name : data_names)
3605 out <<
'<' << data_name <<
"> ";
3611 for (
const auto &patch : patches)
3614 const unsigned int n_points_per_direction = n_subdivisions + 1;
3616 Assert((patch.
data.n_rows() == n_data_sets &&
3618 (patch.
data.n_rows() == n_data_sets + spacedim &&
3621 (n_data_sets + spacedim) :
3623 patch.
data.n_rows()));
3625 auto output_point_data =
3626 [&out, &patch, n_data_sets](
const unsigned int point_index)
mutable {
3627 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
3628 out << patch.data(data_set, point_index) <<
' ';
3637 Assert(patch.data.n_cols() == 1,
3639 n_subdivisions + 1));
3643 out << get_equispaced_location(patch, {}, n_subdivisions)
3645 output_point_data(0);
3655 Assert(patch.data.n_cols() ==
3656 Utilities::fixed_power<dim>(n_points_per_direction),
3658 n_subdivisions + 1));
3660 for (
unsigned int i1 = 0; i1 < n_points_per_direction; ++i1)
3663 out << get_equispaced_location(patch, {i1}, n_subdivisions)
3666 output_point_data(i1);
3679 Assert(patch.data.n_cols() == Utilities::fixed_power<dim>(
3680 n_points_per_direction),
3682 n_subdivisions + 1));
3684 for (
unsigned int i2 = 0; i2 < n_points_per_direction; ++i2)
3686 for (
unsigned int i1 = 0; i1 < n_points_per_direction;
3690 out << get_equispaced_location(patch,
3695 output_point_data(i1 + i2 * n_points_per_direction);
3720 out << get_node_location(patch, 0) <<
' ';
3721 output_point_data(0);
3724 out << get_node_location(patch, 1) <<
' ';
3725 output_point_data(1);
3729 out << get_node_location(patch, 2) <<
' ';
3730 output_point_data(2);
3733 out << get_node_location(patch, 2) <<
' ';
3734 output_point_data(2);
3754 Assert(patch.data.n_cols() == Utilities::fixed_power<dim>(
3755 n_points_per_direction),
3757 n_subdivisions + 1));
3762 for (
unsigned int i3 = 0; i3 < n_points_per_direction; ++i3)
3763 for (
unsigned int i2 = 0; i2 < n_points_per_direction;
3765 for (
unsigned int i1 = 0; i1 < n_points_per_direction;
3770 get_equispaced_location(patch,
3774 if (i1 < n_subdivisions)
3777 out << this_point <<
' ';
3778 output_point_data(i1 +
3779 i2 * n_points_per_direction +
3780 i3 * n_points_per_direction *
3781 n_points_per_direction);
3785 out << get_equispaced_location(patch,
3790 output_point_data((i1 + 1) +
3791 i2 * n_points_per_direction +
3792 i3 * n_points_per_direction *
3793 n_points_per_direction);
3797 out <<
'\n' <<
'\n';
3801 if (i2 < n_subdivisions)
3804 out << this_point <<
' ';
3805 output_point_data(i1 +
3806 i2 * n_points_per_direction +
3807 i3 * n_points_per_direction *
3808 n_points_per_direction);
3812 out << get_equispaced_location(patch,
3818 i1 + (i2 + 1) * n_points_per_direction +
3819 i3 * n_points_per_direction *
3820 n_points_per_direction);
3824 out <<
'\n' <<
'\n';
3828 if (i3 < n_subdivisions)
3831 out << this_point <<
' ';
3832 output_point_data(i1 +
3833 i2 * n_points_per_direction +
3834 i3 * n_points_per_direction *
3835 n_points_per_direction);
3839 out << get_equispaced_location(patch,
3845 i1 + i2 * n_points_per_direction +
3846 (i3 + 1) * n_points_per_direction *
3847 n_points_per_direction);
3850 out <<
'\n' <<
'\n';
3859 for (
const unsigned int v : {0, 1, 2, 0, 3, 2})
3861 out << get_node_location(patch, v) <<
' ';
3862 output_point_data(v);
3867 for (
const unsigned int v : {3, 1})
3869 out << get_node_location(patch, v) <<
' ';
3870 output_point_data(v);
3880 for (
const unsigned int v : {0, 1, 3, 2, 0, 4, 1})
3882 out << get_node_location(patch, v) <<
' ';
3883 output_point_data(v);
3888 for (
const unsigned int v : {2, 4, 3})
3890 out << get_node_location(patch, v) <<
' ';
3891 output_point_data(v);
3905 for (
const unsigned int v : {0, 1, 2, 0, 3, 4, 5, 3})
3907 out << get_node_location(patch, v) <<
' ';
3908 output_point_data(v);
3913 for (
const unsigned int v : {1, 4})
3915 out << get_node_location(patch, v) <<
' ';
3916 output_point_data(v);
3921 for (
const unsigned int v : {2, 5})
3923 out << get_node_location(patch, v) <<
' ';
3924 output_point_data(v);
3949 template <
int dim,
int spacedim>
3951 do_write_povray(
const std::vector<Patch<dim, spacedim>> &,
3952 const std::vector<std::string> &,
3953 const PovrayFlags &,
3957 ExcMessage(
"Writing files in POVRAY format is only supported "
3958 "for two-dimensional meshes."));
3964 do_write_povray(
const std::vector<Patch<2, 2>> &patches,
3965 const std::vector<std::string> &data_names,
3966 const PovrayFlags & flags,
3971#ifndef DEAL_II_WITH_MPI
3980 if (patches.size() == 0)
3983 constexpr int dim = 2;
3985 constexpr int spacedim = 2;
3987 const unsigned int n_data_sets = data_names.size();
3993 <<
"/* This file was generated by the deal.II library." <<
'\n'
3997 <<
" For a description of the POVRAY format see the POVRAY manual."
4002 out <<
"#include \"colors.inc\" " <<
'\n'
4003 <<
"#include \"textures.inc\" " <<
'\n';
4007 if (flags.external_data)
4008 out <<
"#include \"data.inc\" " <<
'\n';
4014 <<
"camera {" <<
'\n'
4015 <<
" location <1,4,-7>" <<
'\n'
4016 <<
" look_at <0,0,0>" <<
'\n'
4017 <<
" angle 30" <<
'\n'
4022 <<
"light_source {" <<
'\n'
4023 <<
" <1,4,-7>" <<
'\n'
4024 <<
" color Grey" <<
'\n'
4027 <<
"light_source {" <<
'\n'
4028 <<
" <0,20,0>" <<
'\n'
4029 <<
" color White" <<
'\n'
4036 double hmin = patches[0].data(0, 0);
4037 double hmax = patches[0].data(0, 0);
4039 for (
const auto &patch : patches)
4043 Assert((patch.
data.n_rows() == n_data_sets &&
4045 (patch.
data.n_rows() == n_data_sets + spacedim &&
4048 (n_data_sets + spacedim) :
4050 patch.data.n_rows()));
4052 Utilities::fixed_power<dim>(n_subdivisions + 1),
4054 n_subdivisions + 1));
4056 for (
unsigned int i = 0; i < n_subdivisions + 1; ++i)
4057 for (
unsigned int j = 0; j < n_subdivisions + 1; ++j)
4059 const int dl = i * (n_subdivisions + 1) + j;
4060 if (patch.
data(0, dl) < hmin)
4061 hmin = patch.
data(0, dl);
4062 if (patch.
data(0, dl) > hmax)
4063 hmax = patch.
data(0, dl);
4067 out <<
"#declare HMIN=" << hmin <<
";" <<
'\n'
4068 <<
"#declare HMAX=" << hmax <<
";" <<
'\n'
4071 if (!flags.external_data)
4074 out <<
"#declare Tex=texture{" <<
'\n'
4075 <<
" pigment {" <<
'\n'
4076 <<
" gradient y" <<
'\n'
4077 <<
" scale y*(HMAX-HMIN)*" << 0.1 <<
'\n'
4078 <<
" color_map {" <<
'\n'
4079 <<
" [0.00 color Light_Purple] " <<
'\n'
4080 <<
" [0.95 color Light_Purple] " <<
'\n'
4081 <<
" [1.00 color White] " <<
'\n'
4086 if (!flags.bicubic_patch)
4089 out <<
'\n' <<
"mesh {" <<
'\n';
4093 for (
const auto &patch : patches)
4096 const unsigned int n = n_subdivisions + 1;
4097 const unsigned int d1 = 1;
4098 const unsigned int d2 = n;
4100 Assert((patch.
data.n_rows() == n_data_sets &&
4102 (patch.
data.n_rows() == n_data_sets + spacedim &&
4105 (n_data_sets + spacedim) :
4107 patch.data.n_rows()));
4108 Assert(patch.
data.n_cols() == Utilities::fixed_power<dim>(n),
4110 n_subdivisions + 1));
4113 std::vector<Point<spacedim>> ver(n * n);
4115 for (
unsigned int i2 = 0; i2 < n; ++i2)
4116 for (
unsigned int i1 = 0; i1 < n; ++i1)
4119 ver[i1 * d1 + i2 * d2] =
4120 get_equispaced_location(patch, {i1, i2}, n_subdivisions);
4124 if (!flags.bicubic_patch)
4127 std::vector<Point<3>> nrml;
4137 for (
unsigned int i = 0; i < n; ++i)
4138 for (
unsigned int j = 0; j < n; ++j)
4140 const unsigned int il = (i == 0) ? i : (i - 1);
4141 const unsigned int ir =
4142 (i == n_subdivisions) ? i : (i + 1);
4143 const unsigned int jl = (j == 0) ? j : (j - 1);
4144 const unsigned int jr =
4145 (j == n_subdivisions) ? j : (j + 1);
4148 ver[ir * d1 + j * d2](0) - ver[il * d1 + j * d2](0);
4149 h1(1) = patch.
data(0, ir * d1 + j * d2) -
4150 patch.
data(0, il * d1 + j * d2);
4152 ver[ir * d1 + j * d2](1) - ver[il * d1 + j * d2](1);
4155 ver[i * d1 + jr * d2](0) - ver[i * d1 + jl * d2](0);
4156 h2(1) = patch.
data(0, i * d1 + jr * d2) -
4157 patch.
data(0, i * d1 + jl * d2);
4159 ver[i * d1 + jr * d2](1) - ver[i * d1 + jl * d2](1);
4161 nrml[i * d1 + j * d2](0) =
4162 h1(1) * h2(2) - h1(2) * h2(1);
4163 nrml[i * d1 + j * d2](1) =
4164 h1(2) * h2(0) - h1(0) * h2(2);
4165 nrml[i * d1 + j * d2](2) =
4166 h1(0) * h2(1) - h1(1) * h2(0);
4171 std::pow(nrml[i * d1 + j * d2](1), 2.) +
4172 std::pow(nrml[i * d1 + j * d2](2), 2.));
4174 if (nrml[i * d1 + j * d2](1) < 0)
4177 for (
unsigned int k = 0; k < 3; ++k)
4178 nrml[i * d1 + j * d2](k) /=
norm;
4183 for (
unsigned int i = 0; i < n_subdivisions; ++i)
4184 for (
unsigned int j = 0; j < n_subdivisions; ++j)
4187 const int dl = i * d1 + j * d2;
4193 out <<
"smooth_triangle {" <<
'\n'
4194 <<
"\t<" << ver[dl](0) <<
"," << patch.
data(0, dl)
4195 <<
"," << ver[dl](1) <<
">, <" << nrml[dl](0)
4196 <<
", " << nrml[dl](1) <<
", " << nrml[dl](2)
4198 out <<
" \t<" << ver[dl + d1](0) <<
","
4199 << patch.
data(0, dl + d1) <<
"," << ver[dl + d1](1)
4200 <<
">, <" << nrml[dl + d1](0) <<
", "
4201 << nrml[dl + d1](1) <<
", " << nrml[dl + d1](2)
4203 out <<
"\t<" << ver[dl + d1 + d2](0) <<
","
4204 << patch.
data(0, dl + d1 + d2) <<
","
4205 << ver[dl + d1 + d2](1) <<
">, <"
4206 << nrml[dl + d1 + d2](0) <<
", "
4207 << nrml[dl + d1 + d2](1) <<
", "
4208 << nrml[dl + d1 + d2](2) <<
">}" <<
'\n';
4211 out <<
"smooth_triangle {" <<
'\n'
4212 <<
"\t<" << ver[dl](0) <<
"," << patch.
data(0, dl)
4213 <<
"," << ver[dl](1) <<
">, <" << nrml[dl](0)
4214 <<
", " << nrml[dl](1) <<
", " << nrml[dl](2)
4216 out <<
"\t<" << ver[dl + d1 + d2](0) <<
","
4217 << patch.
data(0, dl + d1 + d2) <<
","
4218 << ver[dl + d1 + d2](1) <<
">, <"
4219 << nrml[dl + d1 + d2](0) <<
", "
4220 << nrml[dl + d1 + d2](1) <<
", "
4221 << nrml[dl + d1 + d2](2) <<
">," <<
'\n';
4222 out <<
"\t<" << ver[dl + d2](0) <<
","
4223 << patch.
data(0, dl + d2) <<
"," << ver[dl + d2](1)
4224 <<
">, <" << nrml[dl + d2](0) <<
", "
4225 << nrml[dl + d2](1) <<
", " << nrml[dl + d2](2)
4231 out <<
"triangle {" <<
'\n'
4232 <<
"\t<" << ver[dl](0) <<
"," << patch.
data(0, dl)
4233 <<
"," << ver[dl](1) <<
">," <<
'\n';
4234 out <<
"\t<" << ver[dl + d1](0) <<
","
4235 << patch.
data(0, dl + d1) <<
"," << ver[dl + d1](1)
4237 out <<
"\t<" << ver[dl + d1 + d2](0) <<
","
4238 << patch.
data(0, dl + d1 + d2) <<
","
4239 << ver[dl + d1 + d2](1) <<
">}" <<
'\n';
4242 out <<
"triangle {" <<
'\n'
4243 <<
"\t<" << ver[dl](0) <<
"," << patch.
data(0, dl)
4244 <<
"," << ver[dl](1) <<
">," <<
'\n';
4245 out <<
"\t<" << ver[dl + d1 + d2](0) <<
","
4246 << patch.
data(0, dl + d1 + d2) <<
","
4247 << ver[dl + d1 + d2](1) <<
">," <<
'\n';
4248 out <<
"\t<" << ver[dl + d2](0) <<
","
4249 << patch.
data(0, dl + d2) <<
"," << ver[dl + d2](1)
4257 Assert(n_subdivisions == 3,
4260 <<
"bicubic_patch {" <<
'\n'
4261 <<
" type 0" <<
'\n'
4262 <<
" flatness 0" <<
'\n'
4263 <<
" u_steps 0" <<
'\n'
4264 <<
" v_steps 0" <<
'\n';
4265 for (
int i = 0; i < 16; ++i)
4267 out <<
"\t<" << ver[i](0) <<
"," << patch.
data(0, i) <<
","
4268 << ver[i](1) <<
">";
4273 out <<
" texture {Tex}" <<
'\n' <<
"}" <<
'\n';
4277 if (!flags.bicubic_patch)
4280 out <<
" texture {Tex}" <<
'\n' <<
"}" <<
'\n' <<
'\n';
4292 template <
int dim,
int spacedim>
4296 const std::vector<std::string> & data_names,
4298 std::tuple<
unsigned int,
4305 do_write_povray(patches, data_names, flags, out);
4310 template <
int dim,
int spacedim>
4314 const std::vector<std::string> & ,
4316 std::tuple<
unsigned int,
4328 template <
int spacedim>
4332 const std::vector<std::string> & ,
4334 std::tuple<
unsigned int,
4343#ifndef DEAL_II_WITH_MPI
4352 if (patches.size() == 0)
4362 std::multiset<EpsCell2d> cells;
4366 float min_color_value = std::numeric_limits<float>::max();
4367 float max_color_value = std::numeric_limits<float>::min();
4371 double heights[4] = {0, 0, 0, 0};
4375 for (
const auto &patch : patches)
4378 const unsigned int n = n_subdivisions + 1;
4379 const unsigned int d1 = 1;
4380 const unsigned int d2 = n;
4382 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
4383 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
4387 get_equispaced_location(patch, {i1, i2}, n_subdivisions);
4389 get_equispaced_location(patch, {i1 + 1, i2}, n_subdivisions);
4391 get_equispaced_location(patch, {i1, i2 + 1}, n_subdivisions);
4392 points[3] = get_equispaced_location(patch,
4400 patch.
data.n_rows() == 0,
4403 patch.
data.n_rows()));
4405 patch.
data.n_rows() != 0 ?
4409 heights[1] = patch.
data.n_rows() != 0 ?
4411 (i1 + 1) * d1 + i2 * d2) *
4414 heights[2] = patch.
data.n_rows() != 0 ?
4416 i1 * d1 + (i2 + 1) * d2) *
4419 heights[3] = patch.
data.n_rows() != 0 ?
4421 (i1 + 1) * d1 + (i2 + 1) * d2) *
4428 for (
unsigned int i = 0; i < 4; ++i)
4429 heights[i] = points[i](2);
4454 for (
unsigned int vertex = 0; vertex < 4; ++vertex)
4456 const double x = points[vertex](0), y = points[vertex](1),
4457 z = -heights[vertex];
4459 eps_cell.vertices[vertex](0) = -cz * x + sz * y;
4460 eps_cell.vertices[vertex](1) =
4461 -cx * sz * x - cx * cz * y - sx * z;
4485 (points[0] + points[1] + points[2] + points[3]) / 4;
4486 const double center_height =
4487 -(heights[0] + heights[1] + heights[2] + heights[3]) / 4;
4490 eps_cell.depth = -sx * sz * center_point(0) -
4491 sx * cz * center_point(1) + cx * center_height;
4496 patch.
data.n_rows() == 0,
4499 patch.
data.n_rows()));
4500 const double color_values[4] = {
4501 patch.
data.n_rows() != 0 ?
4505 patch.
data.n_rows() != 0 ?
4509 patch.
data.n_rows() != 0 ?
4513 patch.
data.n_rows() != 0 ?
4515 (i1 + 1) * d1 + (i2 + 1) * d2) :
4519 eps_cell.color_value = (color_values[0] + color_values[1] +
4520 color_values[3] + color_values[2]) /
4525 std::min(min_color_value, eps_cell.color_value);
4527 std::max(max_color_value, eps_cell.color_value);
4531 cells.insert(eps_cell);
4537 double x_min = cells.begin()->vertices[0](0);
4538 double x_max = x_min;
4539 double y_min = cells.begin()->vertices[0](1);
4540 double y_max = y_min;
4542 for (
const auto &cell : cells)
4543 for (
const auto &vertex : cell.vertices)
4545 x_min =
std::min(x_min, vertex(0));
4546 x_max =
std::max(x_max, vertex(0));
4547 y_min =
std::min(y_min, vertex(1));
4548 y_max =
std::max(y_max, vertex(1));
4553 const double scale =
4557 const Point<2> offset(x_min, y_min);
4562 out <<
"%!PS-Adobe-2.0 EPSF-1.2" <<
'\n'
4563 <<
"%%Title: deal.II Output" <<
'\n'
4564 <<
"%%Creator: the deal.II library" <<
'\n'
4567 <<
"%%BoundingBox: "
4571 <<
static_cast<unsigned int>((x_max - x_min) * scale + 0.5) <<
' '
4572 <<
static_cast<unsigned int>((y_max - y_min) * scale + 0.5) <<
'\n';
4581 out <<
"/m {moveto} bind def" <<
'\n'
4582 <<
"/l {lineto} bind def" <<
'\n'
4583 <<
"/s {setrgbcolor} bind def" <<
'\n'
4584 <<
"/sg {setgray} bind def" <<
'\n'
4585 <<
"/lx {lineto closepath stroke} bind def" <<
'\n'
4586 <<
"/lf {lineto closepath fill} bind def" <<
'\n';
4588 out <<
"%%EndProlog" <<
'\n' <<
'\n';
4590 out << flags.
line_width <<
" setlinewidth" <<
'\n';
4598 if (max_color_value == min_color_value)
4599 max_color_value = min_color_value + 1;
4603 for (
const auto &cell : cells)
4616 out << rgb_values.
red <<
" sg ";
4618 out << rgb_values.
red <<
' ' << rgb_values.
green <<
' '
4619 << rgb_values.
blue <<
" s ";
4624 out << (cell.vertices[0] - offset) * scale <<
" m "
4625 << (cell.vertices[1] - offset) * scale <<
" l "
4626 << (cell.vertices[3] - offset) * scale <<
" l "
4627 << (cell.vertices[2] - offset) * scale <<
" lf" <<
'\n';
4632 << (cell.vertices[0] - offset) * scale <<
" m "
4633 << (cell.vertices[1] - offset) * scale <<
" l "
4634 << (cell.vertices[3] - offset) * scale <<
" l "
4635 << (cell.vertices[2] - offset) * scale <<
" lx" <<
'\n';
4637 out <<
"showpage" <<
'\n';
4646 template <
int dim,
int spacedim>
4650 const std::vector<std::string> & data_names,
4652 std::tuple<
unsigned int,
4668#ifndef DEAL_II_WITH_MPI
4677 if (patches.size() == 0)
4681 GmvStream gmv_out(out, flags);
4682 const unsigned int n_data_sets = data_names.size();
4685 Assert((patches[0].data.n_rows() == n_data_sets &&
4686 !patches[0].points_are_available) ||
4687 (patches[0].data.n_rows() == n_data_sets + spacedim &&
4688 patches[0].points_are_available),
4690 (n_data_sets + spacedim) :
4692 patches[0].data.n_rows()));
4696 out <<
"gmvinput ascii" <<
'\n' <<
'\n';
4699 unsigned int n_nodes;
4700 unsigned int n_cells;
4701 std::tie(n_nodes, n_cells) = count_nodes_and_cells(patches);
4716 [&patches]() {
return create_global_data_table(patches); });
4722 out <<
"nodes " << n_nodes <<
'\n';
4723 for (
unsigned int d = 0; d < spacedim; ++d)
4725 gmv_out.selected_component = d;
4731 for (
unsigned int d = spacedim; d < 3; ++d)
4733 for (
unsigned int i = 0; i < n_nodes; ++i)
4740 out <<
"cells " << n_cells <<
'\n';
4745 out <<
"variable" <<
'\n';
4749 std::move(*create_global_data_table_task.
return_value());
4753 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
4755 out << data_names[data_set] <<
" 1" <<
'\n';
4756 std::copy(data_vectors[data_set].begin(),
4757 data_vectors[data_set].end(),
4758 std::ostream_iterator<double>(out,
" "));
4759 out <<
'\n' <<
'\n';
4765 out <<
"endvars" <<
'\n';
4768 out <<
"endgmv" <<
'\n';
4779 template <
int dim,
int spacedim>
4783 const std::vector<std::string> & data_names,
4785 std::tuple<
unsigned int,
4799#ifndef DEAL_II_WITH_MPI
4808 if (patches.size() == 0)
4812 TecplotStream tecplot_out(out, flags);
4814 const unsigned int n_data_sets = data_names.size();
4817 Assert((patches[0].data.n_rows() == n_data_sets &&
4818 !patches[0].points_are_available) ||
4819 (patches[0].data.n_rows() == n_data_sets + spacedim &&
4820 patches[0].points_are_available),
4822 (n_data_sets + spacedim) :
4824 patches[0].data.n_rows()));
4827 unsigned int n_nodes;
4828 unsigned int n_cells;
4829 std::tie(n_nodes, n_cells) = count_nodes_and_cells(patches);
4835 <<
"# This file was generated by the deal.II library." <<
'\n'
4839 <<
"# For a description of the Tecplot format see the Tecplot documentation."
4844 out <<
"Variables=";
4852 out <<
"\"x\", \"y\"";
4855 out <<
"\"x\", \"y\", \"z\"";
4861 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
4862 out <<
", \"" << data_names[data_set] <<
"\"";
4868 out <<
"t=\"" << flags.
zone_name <<
"\" ";
4871 out <<
"strandid=1, solutiontime=" << flags.
solution_time <<
", ";
4873 out <<
"f=feblock, n=" << n_nodes <<
", e=" << n_cells
4874 <<
", et=" << tecplot_cell_type[dim] <<
'\n';
4891 [&patches]() {
return create_global_data_table(patches); });
4897 for (
unsigned int d = 0; d < spacedim; ++d)
4899 tecplot_out.selected_component = d;
4910 std::move(*create_global_data_table_task.
return_value());
4913 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
4915 std::copy(data_vectors[data_set].begin(),
4916 data_vectors[data_set].end(),
4917 std::ostream_iterator<double>(out,
"\n"));
4932 template <
int dim,
int spacedim>
4936 const std::vector<std::string> & data_names,
4938 std::tuple<
unsigned int,
4942 & nonscalar_data_ranges,
4948#ifndef DEAL_II_WITH_MPI
4957 if (patches.size() == 0)
4961 VtkStream vtk_out(out, flags);
4963 const unsigned int n_data_sets = data_names.size();
4965 if (patches[0].points_are_available)
4977 out <<
"# vtk DataFile Version 3.0" <<
'\n'
4978 <<
"#This file was generated by the deal.II library";
4986 out <<
'\n' <<
"ASCII" <<
'\n';
4988 out <<
"DATASET UNSTRUCTURED_GRID\n" <<
'\n';
4995 const unsigned int n_metadata =
4996 ((flags.
cycle != std::numeric_limits<unsigned int>::min() ? 1 : 0) +
4997 (flags.
time != std::numeric_limits<double>::min() ? 1 : 0));
5000 out <<
"FIELD FieldData " << n_metadata <<
'\n';
5002 if (flags.
cycle != std::numeric_limits<unsigned int>::min())
5004 out <<
"CYCLE 1 1 int\n" << flags.
cycle <<
'\n';
5006 if (flags.
time != std::numeric_limits<double>::min())
5008 out <<
"TIME 1 1 double\n" << flags.
time <<
'\n';
5014 unsigned int n_nodes;
5015 unsigned int n_cells;
5016 unsigned int n_points_and_n_cells;
5017 std::tie(n_nodes, n_cells, n_points_and_n_cells) =
5033 [&patches]() {
return create_global_data_table(patches); });
5039 out <<
"POINTS " << n_nodes <<
" double" <<
'\n';
5044 out <<
"CELLS " << n_cells <<
' ' << n_points_and_n_cells <<
'\n';
5052 out <<
"CELL_TYPES " << n_cells <<
'\n';
5056 for (
const auto &patch : patches)
5058 const auto vtk_cell_id =
5061 for (
unsigned int i = 0; i < vtk_cell_id[1]; ++i)
5062 out <<
' ' << vtk_cell_id[0];
5071 std::move(*create_global_data_table_task.
return_value());
5076 out <<
"POINT_DATA " << n_nodes <<
'\n';
5080 std::vector<bool> data_set_written(n_data_sets,
false);
5081 for (
const auto &nonscalar_data_range : nonscalar_data_ranges)
5086 "The VTK writer does not currently support outputting "
5087 "tensor data. Use the VTU writer instead."));
5090 std::get<0>(nonscalar_data_range),
5092 std::get<0>(nonscalar_data_range)));
5093 AssertThrow(std::get<1>(nonscalar_data_range) < n_data_sets,
5097 AssertThrow(std::get<1>(nonscalar_data_range) + 1 -
5098 std::get<0>(nonscalar_data_range) <=
5101 "Can't declare a vector with more than 3 components "
5105 for (
unsigned int i = std::get<0>(nonscalar_data_range);
5106 i <= std::get<1>(nonscalar_data_range);
5108 data_set_written[i] =
true;
5114 if (!std::get<2>(nonscalar_data_range).empty())
5115 out << std::get<2>(nonscalar_data_range);
5118 for (
unsigned int i = std::get<0>(nonscalar_data_range);
5119 i < std::get<1>(nonscalar_data_range);
5121 out << data_names[i] <<
"__";
5122 out << data_names[std::get<1>(nonscalar_data_range)];
5125 out <<
" double" <<
'\n';
5128 for (
unsigned int n = 0; n < n_nodes; ++n)
5130 switch (std::get<1>(nonscalar_data_range) -
5131 std::get<0>(nonscalar_data_range))
5134 out << data_vectors(std::get<0>(nonscalar_data_range), n)
5139 out << data_vectors(std::get<0>(nonscalar_data_range), n)
5141 << data_vectors(std::get<0>(nonscalar_data_range) + 1, n)
5145 out << data_vectors(std::get<0>(nonscalar_data_range), n)
5147 << data_vectors(std::get<0>(nonscalar_data_range) + 1, n)
5149 << data_vectors(std::get<0>(nonscalar_data_range) + 2, n)
5162 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
5163 if (data_set_written[data_set] ==
false)
5165 out <<
"SCALARS " << data_names[data_set] <<
" double 1" <<
'\n'
5166 <<
"LOOKUP_TABLE default" <<
'\n';
5167 std::copy(data_vectors[data_set].begin(),
5168 data_vectors[data_set].end(),
5169 std::ostream_iterator<double>(out,
" "));
5185 out <<
"<?xml version=\"1.0\" ?> \n";
5187 out <<
"# vtk DataFile Version 3.0" <<
'\n'
5188 <<
"#This file was generated by the deal.II library";
5199 out <<
"<VTKFile type=\"UnstructuredGrid\" version=\"2.2\"";
5201 out <<
"<VTKFile type=\"UnstructuredGrid\" version=\"0.1\"";
5202 if (deal_ii_with_zlib &&
5204 out <<
" compressor=\"vtkZLibDataCompressor\"";
5205#ifdef DEAL_II_WORDS_BIGENDIAN
5206 out <<
" byte_order=\"BigEndian\"";
5208 out <<
" byte_order=\"LittleEndian\"";
5212 out <<
"<UnstructuredGrid>";
5222 out <<
" </UnstructuredGrid>\n";
5223 out <<
"</VTKFile>\n";
5228 template <
int dim,
int spacedim>
5232 const std::vector<std::string> & data_names,
5234 std::tuple<
unsigned int,
5238 & nonscalar_data_ranges,
5243 write_vtu_main(patches, data_names, nonscalar_data_ranges, flags, out);
5250 template <
int dim,
int spacedim>
5254 const std::vector<std::string> & data_names,
5256 std::tuple<
unsigned int,
5260 & nonscalar_data_ranges,
5274 unit.second.find(
'\"') == std::string::npos,
5276 "A physical unit you provided, <" + unit.second +
5277 ">, contained a quotation mark character. This is not allowed."));
5280#ifndef DEAL_II_WITH_MPI
5289 if (patches.size() == 0)
5294 out <<
"<Piece NumberOfPoints=\"0\" NumberOfCells=\"0\" >\n"
5296 <<
"<DataArray type=\"UInt8\" Name=\"types\"></DataArray>\n"
5298 <<
" <PointData Scalars=\"scalars\">\n";
5299 std::vector<bool> data_set_written(data_names.size(),
false);
5300 for (
const auto &nonscalar_data_range : nonscalar_data_ranges)
5303 for (
unsigned int i = std::get<0>(nonscalar_data_range);
5304 i <= std::get<1>(nonscalar_data_range);
5306 data_set_written[i] =
true;
5310 out <<
" <DataArray type=\"Float32\" Name=\"";
5312 if (!std::get<2>(nonscalar_data_range).empty())
5313 out << std::get<2>(nonscalar_data_range);
5316 for (
unsigned int i = std::get<0>(nonscalar_data_range);
5317 i < std::get<1>(nonscalar_data_range);
5319 out << data_names[i] <<
"__";
5320 out << data_names[std::get<1>(nonscalar_data_range)];
5323 out <<
"\" NumberOfComponents=\"3\"></DataArray>\n";
5326 for (
unsigned int data_set = 0; data_set < data_names.size();
5328 if (data_set_written[data_set] ==
false)
5330 out <<
" <DataArray type=\"Float32\" Name=\""
5331 << data_names[data_set] <<
"\"></DataArray>\n";
5334 out <<
" </PointData>\n";
5335 out <<
"</Piece>\n";
5349 const unsigned int n_metadata =
5350 ((flags.
cycle != std::numeric_limits<unsigned int>::min() ? 1 : 0) +
5351 (flags.
time != std::numeric_limits<double>::min() ? 1 : 0));
5353 out <<
"<FieldData>\n";
5355 if (flags.
cycle != std::numeric_limits<unsigned int>::min())
5358 <<
"<DataArray type=\"Float32\" Name=\"CYCLE\" NumberOfTuples=\"1\" format=\"ascii\">"
5359 << flags.
cycle <<
"</DataArray>\n";
5361 if (flags.
time != std::numeric_limits<double>::min())
5364 <<
"<DataArray type=\"Float32\" Name=\"TIME\" NumberOfTuples=\"1\" format=\"ascii\">"
5365 << flags.
time <<
"</DataArray>\n";
5369 out <<
"</FieldData>\n";
5373 const unsigned int n_data_sets = data_names.size();
5376 if (patches[0].points_are_available)
5385 const char *ascii_or_binary =
5386 (deal_ii_with_zlib &&
5393 unsigned int n_nodes;
5394 unsigned int n_cells;
5395 std::tie(n_nodes, n_cells, std::ignore) =
5403 const auto stringize_vertex_information = [&patches,
5407 ascii_or_binary]() {
5408 std::ostringstream o;
5410 o <<
" <DataArray type=\"Float32\" NumberOfComponents=\"3\" format=\""
5411 << ascii_or_binary <<
"\">\n";
5412 const std::vector<Point<spacedim>> node_positions =
5417 std::vector<float> node_coordinates_3d;
5418 node_coordinates_3d.reserve(node_positions.size() * 3);
5419 for (
const auto &node_position : node_positions)
5421 for (
unsigned int d = 0; d < 3; ++d)
5423 node_coordinates_3d.emplace_back(node_position[d]);
5425 node_coordinates_3d.emplace_back(0.0f);
5427 o << vtu_stringize_array(node_coordinates_3d,
5431 o <<
" </DataArray>\n";
5432 o <<
" </Points>\n\n";
5441 const auto stringize_cell_to_vertex_information = [&patches,
5445 out.precision()]() {
5446 std::ostringstream o;
5449 o <<
" <DataArray type=\"Int32\" Name=\"connectivity\" format=\""
5450 << ascii_or_binary <<
"\">\n";
5452 std::vector<int32_t> cells;
5455 unsigned int first_vertex_of_patch = 0;
5457 for (
const auto &patch : patches)
5469 const unsigned int n_points = patch.
data.n_cols();
5470 Assert((dim == 2 && n_points == 6) ||
5471 (dim == 3 && n_points == 10),
5474 if (deal_ii_with_zlib &&
5478 for (
unsigned int i = 0; i < n_points; ++i)
5479 cells.push_back(first_vertex_of_patch + i);
5483 for (
unsigned int i = 0; i < n_points; ++i)
5484 o <<
'\t' << first_vertex_of_patch + i;
5488 first_vertex_of_patch += n_points;
5493 else if (patch.
reference_cell != ReferenceCells::get_hypercube<dim>())
5497 const unsigned int n_points = patch.
data.n_cols();
5499 if (deal_ii_with_zlib &&
5503 for (
unsigned int i = 0; i < n_points; ++i)
5505 first_vertex_of_patch +
5510 for (
unsigned int i = 0; i < n_points; ++i)
5512 << (first_vertex_of_patch +
5517 first_vertex_of_patch += n_points;
5522 const unsigned int n_points_per_direction = n_subdivisions + 1;
5524 std::vector<unsigned> local_vertex_order;
5528 const auto flush_current_cell = [&flags,
5531 first_vertex_of_patch,
5532 &local_vertex_order]() {
5533 if (deal_ii_with_zlib &&
5537 for (
const auto &c : local_vertex_order)
5538 cells.push_back(first_vertex_of_patch + c);
5542 for (
const auto &c : local_vertex_order)
5543 o <<
'\t' << first_vertex_of_patch + c;
5547 local_vertex_order.clear();
5552 local_vertex_order.reserve(Utilities::fixed_power<dim>(2));
5558 local_vertex_order.emplace_back(0);
5559 flush_current_cell();
5565 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
5567 const unsigned int starting_offset = i1;
5568 local_vertex_order.emplace_back(starting_offset);
5569 local_vertex_order.emplace_back(starting_offset +
5571 flush_current_cell();
5578 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
5579 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
5581 const unsigned int starting_offset =
5582 i2 * n_points_per_direction + i1;
5583 local_vertex_order.emplace_back(
5585 local_vertex_order.emplace_back(
5586 starting_offset + 1);
5587 local_vertex_order.emplace_back(
5588 starting_offset + n_points_per_direction + 1);
5589 local_vertex_order.emplace_back(
5590 starting_offset + n_points_per_direction);
5591 flush_current_cell();
5598 for (
unsigned int i3 = 0; i3 < n_subdivisions; ++i3)
5599 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
5600 for (
unsigned int i1 = 0; i1 < n_subdivisions;
5603 const unsigned int starting_offset =
5604 i3 * n_points_per_direction *
5605 n_points_per_direction +
5606 i2 * n_points_per_direction + i1;
5607 local_vertex_order.emplace_back(
5609 local_vertex_order.emplace_back(
5610 starting_offset + 1);
5611 local_vertex_order.emplace_back(
5612 starting_offset + n_points_per_direction +
5614 local_vertex_order.emplace_back(
5615 starting_offset + n_points_per_direction);
5616 local_vertex_order.emplace_back(
5617 starting_offset + n_points_per_direction *
5618 n_points_per_direction);
5619 local_vertex_order.emplace_back(
5621 n_points_per_direction *
5622 n_points_per_direction +
5624 local_vertex_order.emplace_back(
5626 n_points_per_direction *
5627 n_points_per_direction +
5628 n_points_per_direction + 1);
5629 local_vertex_order.emplace_back(
5631 n_points_per_direction *
5632 n_points_per_direction +
5633 n_points_per_direction);
5634 flush_current_cell();
5645 local_vertex_order.resize(
5646 Utilities::fixed_power<dim>(n_points_per_direction));
5654 "Point-like cells should not be possible "
5655 "when writing higher-order cells."));
5660 for (
unsigned int i1 = 0; i1 < n_subdivisions + 1;
5663 const unsigned int local_index = i1;
5664 const unsigned int connectivity_index =
5666 .template vtk_lexicographic_to_node_index<1>(
5670 local_vertex_order[connectivity_index] =
5672 flush_current_cell();
5679 for (
unsigned int i2 = 0; i2 < n_subdivisions + 1;
5681 for (
unsigned int i1 = 0; i1 < n_subdivisions + 1;
5684 const unsigned int local_index =
5685 i2 * n_points_per_direction + i1;
5686 const unsigned int connectivity_index =
5688 .template vtk_lexicographic_to_node_index<
5690 {{n_subdivisions, n_subdivisions}},
5692 local_vertex_order[connectivity_index] =
5695 flush_current_cell();
5701 for (
unsigned int i3 = 0; i3 < n_subdivisions + 1;
5703 for (
unsigned int i2 = 0; i2 < n_subdivisions + 1;
5705 for (
unsigned int i1 = 0; i1 < n_subdivisions + 1;
5708 const unsigned int local_index =
5709 i3 * n_points_per_direction *
5710 n_points_per_direction +
5711 i2 * n_points_per_direction + i1;
5712 const unsigned int connectivity_index =
5714 .template vtk_lexicographic_to_node_index<
5720 local_vertex_order[connectivity_index] =
5724 flush_current_cell();
5734 first_vertex_of_patch +=
5743 o << vtu_stringize_array(cells,
5748 o <<
" </DataArray>\n";
5768 const auto stringize_cell_offset_and_type_information =
5773 output_precision = out.precision()]() {
5774 std::ostringstream o;
5776 o <<
" <DataArray type=\"Int32\" Name=\"offsets\" format=\""
5777 << ascii_or_binary <<
"\">\n";
5779 std::vector<int32_t> offsets;
5780 offsets.reserve(n_cells);
5784 std::vector<unsigned int> cell_types;
5785 cell_types.reserve(n_cells);
5787 unsigned int first_vertex_of_patch = 0;
5789 for (
const auto &patch : patches)
5791 const auto vtk_cell_id =
5794 for (
unsigned int i = 0; i < vtk_cell_id[1]; ++i)
5796 cell_types.push_back(vtk_cell_id[0]);
5797 first_vertex_of_patch += vtk_cell_id[2];
5798 offsets.push_back(first_vertex_of_patch);
5802 o << vtu_stringize_array(offsets,
5806 o <<
" </DataArray>\n";
5808 o <<
" <DataArray type=\"UInt8\" Name=\"types\" format=\""
5809 << ascii_or_binary <<
"\">\n";
5811 if (deal_ii_with_zlib &&
5814 std::vector<uint8_t> cell_types_uint8_t(cell_types.size());
5815 for (
unsigned int i = 0; i < cell_types.size(); ++i)
5816 cell_types_uint8_t[i] =
static_cast<std::uint8_t
>(cell_types[i]);
5818 o << vtu_stringize_array(cell_types_uint8_t,
5824 o << vtu_stringize_array(cell_types,
5830 o <<
" </DataArray>\n";
5840 const auto stringize_nonscalar_data_range =
5846 output_precision = out.precision()](
const Table<2, float> &data_vectors,
5847 const auto & range) {
5848 std::ostringstream o;
5850 const auto first_component = std::get<0>(range);
5851 const auto last_component = std::get<1>(range);
5852 const auto &name = std::get<2>(range);
5853 const bool is_tensor =
5854 (std::get<3>(range) ==
5856 const unsigned int n_components = (is_tensor ? 9 : 3);
5863 AssertThrow((last_component + 1 - first_component <= 9),
5865 "Can't declare a tensor with more than 9 components "
5866 "in VTK/VTU format."));
5870 AssertThrow((last_component + 1 - first_component <= 3),
5872 "Can't declare a vector with more than 3 components "
5873 "in VTK/VTU format."));
5878 o <<
" <DataArray type=\"Float32\" Name=\"";
5884 for (
unsigned int i = first_component; i < last_component; ++i)
5885 o << data_names[i] <<
"__";
5886 o << data_names[last_component];
5889 o <<
"\" NumberOfComponents=\"" << n_components <<
"\" format=\""
5890 << ascii_or_binary <<
"\"";
5909 std::vector<float> data;
5910 data.reserve(n_nodes * n_components);
5912 for (
unsigned int n = 0; n < n_nodes; ++n)
5916 switch (last_component - first_component)
5919 data.push_back(data_vectors(first_component, n));
5925 data.push_back(data_vectors(first_component, n));
5926 data.push_back(data_vectors(first_component + 1, n));
5931 data.push_back(data_vectors(first_component, n));
5932 data.push_back(data_vectors(first_component + 1, n));
5933 data.push_back(data_vectors(first_component + 2, n));
5946 const unsigned int size = last_component - first_component + 1;
5950 vtk_data[0][0] = data_vectors(first_component, n);
5955 for (
unsigned int c = 0; c < size; ++c)
5959 vtk_data[ind[0]][ind[1]] =
5960 data_vectors(first_component + c, n);
5966 for (
unsigned int c = 0; c < size; ++c)
5970 vtk_data[ind[0]][ind[1]] =
5971 data_vectors(first_component + c, n);
5982 for (
unsigned int i = 0; i < 3; ++i)
5983 for (
unsigned int j = 0; j < 3; ++j)
5984 data.push_back(vtk_data[i][j]);
5988 o << vtu_stringize_array(data,
5992 o <<
" </DataArray>\n";
5997 const auto stringize_scalar_data_set =
6001 output_precision = out.precision()](
const Table<2, float> &data_vectors,
6002 const unsigned int data_set) {
6003 std::ostringstream o;
6005 o <<
" <DataArray type=\"Float32\" Name=\"" << data_names[data_set]
6006 <<
"\" format=\"" << ascii_or_binary <<
"\"";
6010 o <<
" units=\"" << flags.
physical_units.at(data_names[data_set])
6015 const std::vector<float> data(data_vectors[data_set].begin(),
6016 data_vectors[data_set].end());
6017 o << vtu_stringize_array(data,
6021 o <<
" </DataArray>\n";
6040 return create_global_data_table<dim, spacedim, float>(patches);
6054 std::move(*create_global_data_table_task.
return_value());
6060 std::vector<bool> data_set_handled(n_data_sets,
false);
6061 for (
const auto &range : nonscalar_data_ranges)
6064 const auto first_component = std::get<0>(range);
6065 const auto last_component = std::get<1>(range);
6066 for (
unsigned int i = first_component; i <= last_component; ++i)
6067 data_set_handled[i] =
true;
6070 return stringize_nonscalar_data_range(data_vectors, range);
6075 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
6076 if (data_set_handled[data_set] ==
false)
6079 return stringize_scalar_data_set(data_vectors, data_set);
6085 out <<
"<Piece NumberOfPoints=\"" << n_nodes <<
"\" NumberOfCells=\""
6086 << n_cells <<
"\" >\n";
6087 for (
const auto &s : mesh_tasks.return_values())
6089 out <<
" <PointData Scalars=\"scalars\">\n";
6092 out <<
" </PointData>\n";
6093 out <<
" </Piece>\n";
6107 const std::vector<std::string> &piece_names,
6108 const std::vector<std::string> &data_names,
6110 std::tuple<
unsigned int,
6114 & nonscalar_data_ranges,
6127 unit.second.find(
'\"') == std::string::npos,
6129 "A physical unit you provided, <" + unit.second +
6130 ">, contained a quotation mark character. This is not allowed."));
6133 const unsigned int n_data_sets = data_names.size();
6135 out <<
"<?xml version=\"1.0\"?>\n";
6138 out <<
"#This file was generated by the deal.II library"
6143 <<
"<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n";
6144 out <<
" <PUnstructuredGrid GhostLevel=\"0\">\n";
6145 out <<
" <PPointData Scalars=\"scalars\">\n";
6148 std::vector<bool> data_set_written(n_data_sets,
false);
6149 for (
const auto &nonscalar_data_range : nonscalar_data_ranges)
6151 const auto first_component = std::get<0>(nonscalar_data_range);
6152 const auto last_component = std::get<1>(nonscalar_data_range);
6153 const bool is_tensor =
6154 (std::get<3>(nonscalar_data_range) ==
6156 const unsigned int n_components = (is_tensor ? 9 : 3);
6163 AssertThrow((last_component + 1 - first_component <= 9),
6165 "Can't declare a tensor with more than 9 components "
6170 Assert((last_component + 1 - first_component <= 3),
6172 "Can't declare a vector with more than 3 components "
6177 for (
unsigned int i = std::get<0>(nonscalar_data_range);
6178 i <= std::get<1>(nonscalar_data_range);
6180 data_set_written[i] =
true;
6184 out <<
" <PDataArray type=\"Float32\" Name=\"";
6186 const std::string &name = std::get<2>(nonscalar_data_range);
6191 for (
unsigned int i = std::get<0>(nonscalar_data_range);
6192 i < std::get<1>(nonscalar_data_range);
6194 out << data_names[i] <<
"__";
6195 out << data_names[std::get<1>(nonscalar_data_range)];
6198 out <<
"\" NumberOfComponents=\"" << n_components
6199 <<
"\" format=\"ascii\"";
6211 data_names[std::get<1>(nonscalar_data_range)]) !=
6215 data_names[std::get<1>(nonscalar_data_range)])
6223 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
6224 if (data_set_written[data_set] ==
false)
6226 out <<
" <PDataArray type=\"Float32\" Name=\""
6227 << data_names[data_set] <<
"\" format=\"ascii\"";
6231 out <<
" units=\"" << flags.
physical_units.at(data_names[data_set])
6237 out <<
" </PPointData>\n";
6239 out <<
" <PPoints>\n";
6240 out <<
" <PDataArray type=\"Float32\" NumberOfComponents=\"3\"/>\n";
6241 out <<
" </PPoints>\n";
6243 for (
const auto &piece_name : piece_names)
6244 out <<
" <Piece Source=\"" << piece_name <<
"\"/>\n";
6246 out <<
" </PUnstructuredGrid>\n";
6247 out <<
"</VTKFile>\n";
6260 const std::vector<std::pair<double, std::string>> ×_and_names)
6264 out <<
"<?xml version=\"1.0\"?>\n";
6267 out <<
"#This file was generated by the deal.II library"
6272 <<
"<VTKFile type=\"Collection\" version=\"0.1\" ByteOrder=\"LittleEndian\">\n";
6273 out <<
" <Collection>\n";
6275 std::streamsize ss = out.precision();
6278 for (
const auto &time_and_name : times_and_names)
6279 out <<
" <DataSet timestep=\"" << time_and_name.first
6280 <<
"\" group=\"\" part=\"0\" file=\"" << time_and_name.second
6283 out <<
" </Collection>\n";
6284 out <<
"</VTKFile>\n";
6296 const std::vector<std::string> &piece_names)
6298 out <<
"!NBLOCKS " << piece_names.size() <<
'\n';
6299 for (
const auto &piece_name : piece_names)
6300 out << piece_name <<
'\n';
6309 const std::vector<std::vector<std::string>> &piece_names)
6313 if (piece_names.size() == 0)
6316 const double nblocks = piece_names[0].size();
6318 ExcMessage(
"piece_names should be a vector of nonempty vectors."));
6320 out <<
"!NBLOCKS " << nblocks <<
'\n';
6321 for (
const auto &domain : piece_names)
6323 Assert(domain.size() == nblocks,
6325 "piece_names should be a vector of equal sized vectors."));
6326 for (
const auto &subdomain : domain)
6327 out << subdomain <<
'\n';
6338 const std::vector<std::pair<
double, std::vector<std::string>>>
6339 ×_and_piece_names)
6343 if (times_and_piece_names.size() == 0)
6346 const double nblocks = times_and_piece_names[0].second.size();
6350 "time_and_piece_names should contain nonempty vectors of filenames for every timestep."));
6352 for (
const auto &domain : times_and_piece_names)
6353 out <<
"!TIME " << domain.first <<
'\n';
6355 out <<
"!NBLOCKS " << nblocks <<
'\n';
6356 for (
const auto &domain : times_and_piece_names)
6358 Assert(domain.second.size() == nblocks,
6360 "piece_names should be a vector of equal sized vectors."));
6361 for (
const auto &subdomain : domain.second)
6362 out << subdomain <<
'\n';
6370 template <
int dim,
int spacedim>
6374 const std::vector<std::string> &,
6376 std::tuple<
unsigned int,
6386 template <
int spacedim>
6390 const std::vector<std::string> & ,
6392 std::tuple<
unsigned int,
6400 const unsigned int height = flags.
height;
6401 unsigned int width = flags.
width;
6404 unsigned int margin_in_percent = 0;
6406 margin_in_percent = 5;
6410 double x_dimension, y_dimension, z_dimension;
6412 const auto &first_patch = patches[0];
6414 unsigned int n_subdivisions = first_patch.n_subdivisions;
6415 unsigned int n = n_subdivisions + 1;
6416 const unsigned int d1 = 1;
6417 const unsigned int d2 = n;
6420 std::array<Point<spacedim>, 4> projected_points;
6423 std::array<Point<2>, 4> projection_decompositions;
6426 get_equispaced_location(first_patch, {0, 0}, n_subdivisions);
6428 if (first_patch.data.n_rows() != 0)
6433 double x_min = projected_point[0];
6434 double x_max = x_min;
6435 double y_min = projected_point[1];
6436 double y_max = y_min;
6437 double z_min = first_patch.data.n_rows() != 0 ?
6440 double z_max = z_min;
6443 for (
const auto &patch : patches)
6446 n = n_subdivisions + 1;
6448 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
6450 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
6452 projected_points[0] =
6453 get_equispaced_location(patch, {i1, i2}, n_subdivisions);
6454 projected_points[1] =
6455 get_equispaced_location(patch, {i1 + 1, i2}, n_subdivisions);
6456 projected_points[2] =
6457 get_equispaced_location(patch, {i1, i2 + 1}, n_subdivisions);
6458 projected_points[3] = get_equispaced_location(patch,
6462 x_min =
std::min(x_min, projected_points[0][0]);
6463 x_min =
std::min(x_min, projected_points[1][0]);
6464 x_min =
std::min(x_min, projected_points[2][0]);
6465 x_min =
std::min(x_min, projected_points[3][0]);
6467 x_max =
std::max(x_max, projected_points[0][0]);
6468 x_max =
std::max(x_max, projected_points[1][0]);
6469 x_max =
std::max(x_max, projected_points[2][0]);
6470 x_max =
std::max(x_max, projected_points[3][0]);
6472 y_min =
std::min(y_min, projected_points[0][1]);
6473 y_min =
std::min(y_min, projected_points[1][1]);
6474 y_min =
std::min(y_min, projected_points[2][1]);
6475 y_min =
std::min(y_min, projected_points[3][1]);
6477 y_max =
std::max(y_max, projected_points[0][1]);
6478 y_max =
std::max(y_max, projected_points[1][1]);
6479 y_max =
std::max(y_max, projected_points[2][1]);
6480 y_max =
std::max(y_max, projected_points[3][1]);
6483 patch.
data.n_rows() == 0,
6486 patch.
data.n_rows()));
6488 z_min = std::min<double>(z_min,
6490 i1 * d1 + i2 * d2));
6491 z_min = std::min<double>(z_min,
6493 (i1 + 1) * d1 + i2 * d2));
6494 z_min = std::min<double>(z_min,
6496 i1 * d1 + (i2 + 1) * d2));
6498 std::min<double>(z_min,
6500 (i1 + 1) * d1 + (i2 + 1) * d2));
6502 z_max = std::max<double>(z_max,
6504 i1 * d1 + i2 * d2));
6505 z_max = std::max<double>(z_max,
6507 (i1 + 1) * d1 + i2 * d2));
6508 z_max = std::max<double>(z_max,
6510 i1 * d1 + (i2 + 1) * d2));
6512 std::max<double>(z_max,
6514 (i1 + 1) * d1 + (i2 + 1) * d2));
6519 x_dimension = x_max - x_min;
6520 y_dimension = y_max - y_min;
6521 z_dimension = z_max - z_min;
6528 float camera_focus = 0;
6531 camera_position[0] = 0.;
6532 camera_position[1] = 0.;
6533 camera_position[2] = z_min + 2. * z_dimension;
6535 camera_direction[0] = 0.;
6536 camera_direction[1] = 0.;
6537 camera_direction[2] = -1.;
6539 camera_horizontal[0] = 1.;
6540 camera_horizontal[1] = 0.;
6541 camera_horizontal[2] = 0.;
6543 camera_focus = .5 * z_dimension;
6549 const float angle_factor = 3.14159265f / 180.f;
6552 camera_position_temp[1] =
6555 camera_position_temp[2] =
6559 camera_direction_temp[1] =
6562 camera_direction_temp[2] =
6566 camera_horizontal_temp[1] =
6569 camera_horizontal_temp[2] =
6573 camera_position[1] = camera_position_temp[1];
6574 camera_position[2] = camera_position_temp[2];
6576 camera_direction[1] = camera_direction_temp[1];
6577 camera_direction[2] = camera_direction_temp[2];
6579 camera_horizontal[1] = camera_horizontal_temp[1];
6580 camera_horizontal[2] = camera_horizontal_temp[2];
6583 camera_position_temp[0] =
6586 camera_position_temp[1] =
6590 camera_direction_temp[0] =
6593 camera_direction_temp[1] =
6597 camera_horizontal_temp[0] =
6600 camera_horizontal_temp[1] =
6604 camera_position[0] = camera_position_temp[0];
6605 camera_position[1] = camera_position_temp[1];
6607 camera_direction[0] = camera_direction_temp[0];
6608 camera_direction[1] = camera_direction_temp[1];
6610 camera_horizontal[0] = camera_horizontal_temp[0];
6611 camera_horizontal[1] = camera_horizontal_temp[1];
6614 camera_position[0] = x_min + .5 * x_dimension;
6615 camera_position[1] = y_min + .5 * y_dimension;
6617 camera_position[0] += (z_min + 2. * z_dimension) *
6620 camera_position[1] -= (z_min + 2. * z_dimension) *
6626 double x_min_perspective, y_min_perspective;
6627 double x_max_perspective, y_max_perspective;
6628 double x_dimension_perspective, y_dimension_perspective;
6630 n_subdivisions = first_patch.n_subdivisions;
6631 n = n_subdivisions + 1;
6636 get_equispaced_location(first_patch, {0, 0}, n_subdivisions);
6638 if (first_patch.data.n_rows() != 0)
6643 point[0] = projected_point[0];
6644 point[1] = projected_point[1];
6645 point[2] = first_patch.data.n_rows() != 0 ?
6649 projection_decomposition = svg_project_point(point,
6655 x_min_perspective = projection_decomposition[0];
6656 x_max_perspective = projection_decomposition[0];
6657 y_min_perspective = projection_decomposition[1];
6658 y_max_perspective = projection_decomposition[1];
6661 for (
const auto &patch : patches)
6664 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
6666 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
6669 {get_equispaced_location(patch, {i1, i2}, n_subdivisions),
6670 get_equispaced_location(patch, {i1 + 1, i2}, n_subdivisions),
6671 get_equispaced_location(patch, {i1, i2 + 1}, n_subdivisions),
6672 get_equispaced_location(patch,
6677 patch.
data.n_rows() == 0,
6680 patch.
data.n_rows()));
6682 const std::array<Point<3>, 4>
vertices = {
6685 patch.
data.n_rows() != 0 ?
6686 patch.
data(0, i1 * d1 + i2 * d2) :
6690 patch.
data.n_rows() != 0 ?
6691 patch.
data(0, (i1 + 1) * d1 + i2 * d2) :
6695 patch.
data.n_rows() != 0 ?
6696 patch.
data(0, i1 * d1 + (i2 + 1) * d2) :
6700 patch.
data.n_rows() != 0 ?
6701 patch.
data(0, (i1 + 1) * d1 + (i2 + 1) * d2) :
6704 projection_decompositions = {
6728 static_cast<double>(
6729 projection_decompositions[0][0]));
6732 static_cast<double>(
6733 projection_decompositions[1][0]));
6736 static_cast<double>(
6737 projection_decompositions[2][0]));
6740 static_cast<double>(
6741 projection_decompositions[3][0]));
6745 static_cast<double>(
6746 projection_decompositions[0][0]));
6749 static_cast<double>(
6750 projection_decompositions[1][0]));
6753 static_cast<double>(
6754 projection_decompositions[2][0]));
6757 static_cast<double>(
6758 projection_decompositions[3][0]));
6762 static_cast<double>(
6763 projection_decompositions[0][1]));
6766 static_cast<double>(
6767 projection_decompositions[1][1]));
6770 static_cast<double>(
6771 projection_decompositions[2][1]));
6774 static_cast<double>(
6775 projection_decompositions[3][1]));
6779 static_cast<double>(
6780 projection_decompositions[0][1]));
6783 static_cast<double>(
6784 projection_decompositions[1][1]));
6787 static_cast<double>(
6788 projection_decompositions[2][1]));
6791 static_cast<double>(
6792 projection_decompositions[3][1]));
6797 x_dimension_perspective = x_max_perspective - x_min_perspective;
6798 y_dimension_perspective = y_max_perspective - y_min_perspective;
6800 std::multiset<SvgCell> cells;
6803 for (
const auto &patch : patches)
6807 for (
unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
6809 for (
unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
6812 {get_equispaced_location(patch, {i1, i2}, n_subdivisions),
6813 get_equispaced_location(patch, {i1 + 1, i2}, n_subdivisions),
6814 get_equispaced_location(patch, {i1, i2 + 1}, n_subdivisions),
6815 get_equispaced_location(patch,
6820 patch.
data.n_rows() == 0,
6823 patch.
data.n_rows()));
6829 cell.vertices[0][2] = patch.
data.n_rows() != 0 ?
6830 patch.
data(0, i1 * d1 + i2 * d2) :
6835 cell.vertices[1][2] = patch.
data.n_rows() != 0 ?
6836 patch.
data(0, (i1 + 1) * d1 + i2 * d2) :
6841 cell.vertices[2][2] = patch.
data.n_rows() != 0 ?
6842 patch.
data(0, i1 * d1 + (i2 + 1) * d2) :
6847 cell.vertices[3][2] =
6848 patch.
data.n_rows() != 0 ?
6849 patch.
data(0, (i1 + 1) * d1 + (i2 + 1) * d2) :
6852 cell.projected_vertices[0] =
6853 svg_project_point(cell.vertices[0],
6858 cell.projected_vertices[1] =
6859 svg_project_point(cell.vertices[1],
6864 cell.projected_vertices[2] =
6865 svg_project_point(cell.vertices[2],
6870 cell.projected_vertices[3] =
6871 svg_project_point(cell.vertices[3],
6877 cell.center = .25 * (cell.vertices[0] + cell.vertices[1] +
6878 cell.vertices[2] + cell.vertices[3]);
6879 cell.projected_center = svg_project_point(cell.center,
6885 cell.depth = cell.center.distance(camera_position);
6895 width =
static_cast<unsigned int>(
6896 .5 + height * (x_dimension_perspective / y_dimension_perspective));
6897 unsigned int additional_width = 0;
6900 additional_width =
static_cast<unsigned int>(
6904 out <<
"<svg width=\"" << width + additional_width <<
"\" height=\""
6905 << height <<
"\" xmlns=\"http://www.w3.org/2000/svg\" version=\"1.1\">"
6907 <<
" <rect width=\"" << width + additional_width <<
"\" height=\""
6908 << height <<
"\" style=\"fill:white\"/>" <<
'\n'
6911 unsigned int triangle_counter = 0;
6914 for (
const auto &cell : cells)
6918 for (
unsigned int triangle_index = 0; triangle_index < 4;
6921 switch (triangle_index)
6924 points3d_triangle[0] = cell.vertices[0],
6925 points3d_triangle[1] = cell.vertices[1],
6926 points3d_triangle[2] = cell.center;
6929 points3d_triangle[0] = cell.vertices[1],
6930 points3d_triangle[1] = cell.vertices[3],
6931 points3d_triangle[2] = cell.center;
6934 points3d_triangle[0] = cell.vertices[3],
6935 points3d_triangle[1] = cell.vertices[2],
6936 points3d_triangle[2] = cell.center;
6939 points3d_triangle[0] = cell.vertices[2],
6940 points3d_triangle[1] = cell.vertices[0],
6941 points3d_triangle[2] = cell.center;
6948 svg_get_gradient_parameters(points3d_triangle);
6951 .667 - ((gradient_param[4] - z_min) / z_dimension) * .667;
6953 .667 - ((gradient_param[5] - z_min) / z_dimension) * .667;
6955 unsigned int start_r = 0;
6956 unsigned int start_g = 0;
6957 unsigned int start_b = 0;
6959 unsigned int stop_r = 0;
6960 unsigned int stop_g = 0;
6961 unsigned int stop_b = 0;
6963 unsigned int start_i =
static_cast<unsigned int>(start_h * 6.);
6964 unsigned int stop_i =
static_cast<unsigned int>(stop_h * 6.);
6966 double start_f = start_h * 6. - start_i;
6967 double start_q = 1. - start_f;
6969 double stop_f = stop_h * 6. - stop_i;
6970 double stop_q = 1. - stop_f;
6972 switch (start_i % 6)
6976 start_g =
static_cast<unsigned int>(.5 + 255. * start_f);
6979 start_r =
static_cast<unsigned int>(.5 + 255. * start_q),
6984 start_b =
static_cast<unsigned int>(.5 + 255. * start_f);
6987 start_g =
static_cast<unsigned int>(.5 + 255. * start_q),
6991 start_r =
static_cast<unsigned int>(.5 + 255. * start_f),
6996 start_b =
static_cast<unsigned int>(.5 + 255. * start_q);
7006 stop_g =
static_cast<unsigned int>(.5 + 255. * stop_f);
7009 stop_r =
static_cast<unsigned int>(.5 + 255. * stop_q),
7014 stop_b =
static_cast<unsigned int>(.5 + 255. * stop_f);
7017 stop_g =
static_cast<unsigned int>(.5 + 255. * stop_q),
7021 stop_r =
static_cast<unsigned int>(.5 + 255. * stop_f),
7026 stop_b =
static_cast<unsigned int>(.5 + 255. * stop_q);
7032 Point<3> gradient_start_point_3d, gradient_stop_point_3d;
7034 gradient_start_point_3d[0] = gradient_param[0];
7035 gradient_start_point_3d[1] = gradient_param[1];
7036 gradient_start_point_3d[2] = gradient_param[4];
7038 gradient_stop_point_3d[0] = gradient_param[2];
7039 gradient_stop_point_3d[1] = gradient_param[3];
7040 gradient_stop_point_3d[2] = gradient_param[5];
7043 svg_project_point(gradient_start_point_3d,
7049 svg_project_point(gradient_stop_point_3d,
7056 out <<
" <linearGradient id=\"" << triangle_counter
7057 <<
"\" gradientUnits=\"userSpaceOnUse\" "
7059 <<
static_cast<unsigned int>(
7061 ((gradient_start_point[0] - x_min_perspective) /
7062 x_dimension_perspective) *
7063 (width - (width / 100.) * 2. * margin_in_percent) +
7064 ((width / 100.) * margin_in_percent))
7067 <<
static_cast<unsigned int>(
7068 .5 + height - (height / 100.) * margin_in_percent -
7069 ((gradient_start_point[1] - y_min_perspective) /
7070 y_dimension_perspective) *
7071 (height - (height / 100.) * 2. * margin_in_percent))
7074 <<
static_cast<unsigned int>(
7076 ((gradient_stop_point[0] - x_min_perspective) /
7077 x_dimension_perspective) *
7078 (width - (width / 100.) * 2. * margin_in_percent) +
7079 ((width / 100.) * margin_in_percent))
7082 <<
static_cast<unsigned int>(
7083 .5 + height - (height / 100.) * margin_in_percent -
7084 ((gradient_stop_point[1] - y_min_perspective) /
7085 y_dimension_perspective) *
7086 (height - (height / 100.) * 2. * margin_in_percent))
7089 <<
" <stop offset=\"0\" style=\"stop-color:rgb(" << start_r
7090 <<
"," << start_g <<
"," << start_b <<
")\"/>" <<
'\n'
7091 <<
" <stop offset=\"1\" style=\"stop-color:rgb(" << stop_r
7092 <<
"," << stop_g <<
"," << stop_b <<
")\"/>" <<
'\n'
7093 <<
" </linearGradient>" <<
'\n';
7096 double x1 = 0, y1 = 0, x2 = 0, y2 = 0;
7097 double x3 = cell.projected_center[0];
7098 double y3 = cell.projected_center[1];
7100 switch (triangle_index)
7103 x1 = cell.projected_vertices[0][0],
7104 y1 = cell.projected_vertices[0][1],
7105 x2 = cell.projected_vertices[1][0],
7106 y2 = cell.projected_vertices[1][1];
7109 x1 = cell.projected_vertices[1][0],
7110 y1 = cell.projected_vertices[1][1],
7111 x2 = cell.projected_vertices[3][0],
7112 y2 = cell.projected_vertices[3][1];
7115 x1 = cell.projected_vertices[3][0],
7116 y1 = cell.projected_vertices[3][1],
7117 x2 = cell.projected_vertices[2][0],
7118 y2 = cell.projected_vertices[2][1];
7121 x1 = cell.projected_vertices[2][0],
7122 y1 = cell.projected_vertices[2][1],
7123 x2 = cell.projected_vertices[0][0],
7124 y2 = cell.projected_vertices[0][1];
7130 out <<
" <path d=\"M "
7131 <<
static_cast<unsigned int>(
7133 ((x1 - x_min_perspective) / x_dimension_perspective) *
7134 (width - (width / 100.) * 2. * margin_in_percent) +
7135 ((width / 100.) * margin_in_percent))
7137 <<
static_cast<unsigned int>(
7138 .5 + height - (height / 100.) * margin_in_percent -
7139 ((y1 - y_min_perspective) / y_dimension_perspective) *
7140 (height - (height / 100.) * 2. * margin_in_percent))
7142 <<
static_cast<unsigned int>(
7144 ((x2 - x_min_perspective) / x_dimension_perspective) *
7145 (width - (width / 100.) * 2. * margin_in_percent) +
7146 ((width / 100.) * margin_in_percent))
7148 <<
static_cast<unsigned int>(
7149 .5 + height - (height / 100.) * margin_in_percent -
7150 ((y2 - y_min_perspective) / y_dimension_perspective) *
7151 (height - (height / 100.) * 2. * margin_in_percent))
7153 <<
static_cast<unsigned int>(
7155 ((x3 - x_min_perspective) / x_dimension_perspective) *
7156 (width - (width / 100.) * 2. * margin_in_percent) +
7157 ((width / 100.) * margin_in_percent))
7159 <<
static_cast<unsigned int>(
7160 .5 + height - (height / 100.) * margin_in_percent -
7161 ((y3 - y_min_perspective) / y_dimension_perspective) *
7162 (height - (height / 100.) * 2. * margin_in_percent))
7164 <<
static_cast<unsigned int>(
7166 ((x1 - x_min_perspective) / x_dimension_perspective) *
7167 (width - (width / 100.) * 2. * margin_in_percent) +
7168 ((width / 100.) * margin_in_percent))
7170 <<
static_cast<unsigned int>(
7171 .5 + height - (height / 100.) * margin_in_percent -
7172 ((y1 - y_min_perspective) / y_dimension_perspective) *
7173 (height - (height / 100.) * 2. * margin_in_percent))
7174 <<
"\" style=\"stroke:black; fill:url(#" << triangle_counter
7185 out <<
'\n' <<
" <!-- colorbar -->" <<
'\n';
7187 unsigned int element_height =
static_cast<unsigned int>(
7188 ((height / 100.) * (71. - 2. * margin_in_percent)) / 4);
7189 unsigned int element_width =
7190 static_cast<unsigned int>(.5 + (height / 100.) * 2.5);
7192 additional_width = 0;
7195 static_cast<unsigned int>(.5 + (height / 100.) * 2.5);
7197 for (
unsigned int index = 0; index < 4; ++index)
7199 double start_h = .667 - ((index + 1) / 4.) * .667;
7200 double stop_h = .667 - (index / 4.) * .667;
7202 unsigned int start_r = 0;
7203 unsigned int start_g = 0;
7204 unsigned int start_b = 0;
7206 unsigned int stop_r = 0;
7207 unsigned int stop_g = 0;
7208 unsigned int stop_b = 0;
7210 unsigned int start_i =
static_cast<unsigned int>(start_h * 6.);
7211 unsigned int stop_i =
static_cast<unsigned int>(stop_h * 6.);
7213 double start_f = start_h * 6. - start_i;
7214 double start_q = 1. - start_f;
7216 double stop_f = stop_h * 6. - stop_i;
7217 double stop_q = 1. - stop_f;
7219 switch (start_i % 6)
7223 start_g =
static_cast<unsigned int>(.5 + 255. * start_f);
7226 start_r =
static_cast<unsigned int>(.5 + 255. * start_q),
7231 start_b =
static_cast<unsigned int>(.5 + 255. * start_f);
7234 start_g =
static_cast<unsigned int>(.5 + 255. * start_q),
7238 start_r =
static_cast<unsigned int>(.5 + 255. * start_f),
7243 start_b =
static_cast<unsigned int>(.5 + 255. * start_q);
7253 stop_g =
static_cast<unsigned int>(.5 + 255. * stop_f);
7256 stop_r =
static_cast<unsigned int>(.5 + 255. * stop_q),
7261 stop_b =
static_cast<unsigned int>(.5 + 255. * stop_f);
7264 stop_g =
static_cast<unsigned int>(.5 + 255. * stop_q),
7268 stop_r =
static_cast<unsigned int>(.5 + 255. * stop_f),
7273 stop_b =
static_cast<unsigned int>(.5 + 255. * stop_q);
7280 out <<
" <linearGradient id=\"colorbar_" << index
7281 <<
"\" gradientUnits=\"userSpaceOnUse\" "
7282 <<
"x1=\"" << width + additional_width <<
"\" "
7284 <<
static_cast<unsigned int>(.5 + (height / 100.) *
7285 (margin_in_percent + 29)) +
7286 (3 - index) * element_height
7288 <<
"x2=\"" << width + additional_width <<
"\" "
7290 <<
static_cast<unsigned int>(.5 + (height / 100.) *
7291 (margin_in_percent + 29)) +
7292 (4 - index) * element_height
7295 <<
" <stop offset=\"0\" style=\"stop-color:rgb(" << start_r
7296 <<
"," << start_g <<
"," << start_b <<
")\"/>" <<
'\n'
7297 <<
" <stop offset=\"1\" style=\"stop-color:rgb(" << stop_r
7298 <<
"," << stop_g <<
"," << stop_b <<
")\"/>" <<
'\n'
7299 <<
" </linearGradient>" <<
'\n';
7304 <<
" x=\"" << width + additional_width <<
"\" y=\""
7305 <<
static_cast<unsigned int>(.5 + (height / 100.) *
7306 (margin_in_percent + 29)) +
7307 (3 - index) * element_height
7308 <<
"\" width=\"" << element_width <<
"\" height=\""
7310 <<
"\" style=\"stroke:black; stroke-width:2; fill:url(#colorbar_"
7311 << index <<
")\"/>" <<
'\n';
7314 for (
unsigned int index = 0; index < 5; ++index)
7318 << width + additional_width +
7319 static_cast<unsigned int>(1.5 * element_width)
7321 <<
static_cast<unsigned int>(
7322 .5 + (height / 100.) * (margin_in_percent + 29) +
7323 (4. - index) * element_height + 30.)
7325 <<
" style=\"text-anchor:start; font-size:80; font-family:Helvetica";
7327 if (index == 0 || index == 4)
7328 out <<
"; font-weight:bold";
7331 <<
static_cast<float>(
7332 (
static_cast<int>((z_min + index * (z_dimension / 4.)) *
7341 out <<
"</text>" <<
'\n';
7346 out <<
'\n' <<
"</svg>";
7352 template <
int dim,
int spacedim>
7356 const std::vector<std::string> & data_names,
7358 std::tuple<
unsigned int,
7362 &nonscalar_data_ranges,
7371 out << dim <<
' ' << spacedim <<
'\n';
7374 out <<
"[deal.II intermediate format graphics data]" <<
'\n'
7380 out << data_names.size() <<
'\n';
7381 for (
const auto &data_name : data_names)
7382 out << data_name <<
'\n';
7384 out << patches.size() <<
'\n';
7385 for (
unsigned int i = 0; i < patches.size(); ++i)
7386 out << patches[i] <<
'\n';
7388 out << nonscalar_data_ranges.size() <<
'\n';
7389 for (
const auto &nonscalar_data_range : nonscalar_data_ranges)
7390 out << std::get<0>(nonscalar_data_range) <<
' '
7391 << std::get<1>(nonscalar_data_range) <<
'\n'
7392 << std::get<2>(nonscalar_data_range) <<
'\n';
7400 template <
int dim,
int spacedim>
7404 const std::vector<std::string> & data_names,
7406 std::tuple<
unsigned int,
7410 & nonscalar_data_ranges,
7412 const std::string & filename,
7416#ifndef DEAL_II_WITH_MPI
7419 (void)nonscalar_data_ranges;
7426 ExcMessage(
"This functionality requires MPI to be enabled."));
7445 std::vector<char> my_buffer;
7447 boost::iostreams::filtering_ostream f;
7453# ifdef DEAL_II_WITH_ZLIB
7454 f.push(boost::iostreams::zlib_compressor(
7455 get_boost_zlib_compression_level(compression)));
7460 "Compression requires deal.II to be configured with ZLIB support."));
7463 boost::iostreams::back_insert_device<std::vector<char>> inserter(
7467 write_deal_II_intermediate<dim, spacedim>(
7468 patches, data_names, nonscalar_data_ranges, flags, f);
7470 const std::uint64_t my_size = my_buffer.size();
7476 const ParallelIntermediateHeader header{
7479 static_cast<std::uint64_t
>(compression),
7489 std::vector<std::uint64_t> chunk_sizes(n_ranks);
7490 int ierr = MPI_Gather(&my_size,
7493 static_cast<std::uint64_t *
>(chunk_sizes.data()),
7501 MPI_Info_create(&info);
7504 ierr = MPI_File_open(
7505 comm, filename.c_str(), MPI_MODE_CREATE | MPI_MODE_WRONLY, info, &fh);
7507 ierr = MPI_Info_free(&info);
7511 ierr = MPI_File_set_size(fh, 0);
7515 ierr = MPI_Barrier(
comm);
7522 fh, 0, &header,
sizeof(header), MPI_CHAR, MPI_STATUS_IGNORE);
7537 std::uint64_t prefix_sum = 0;
7538 ierr = MPI_Exscan(&my_size, &prefix_sum, 1, MPI_UINT64_T, MPI_SUM,
comm);
7542 const MPI_Offset offset =
static_cast<MPI_Offset
>(
sizeof(header)) +
7543 n_ranks *
sizeof(std::uint64_t) + prefix_sum;
7546 fh, offset, my_buffer.data(), my_size, MPI_CHAR, MPI_STATUS_IGNORE);
7555 ierr = MPI_File_sync(fh);
7558 ierr = MPI_File_close(&fh);
7565 std::pair<unsigned int, unsigned int>
7570 unsigned int dim, spacedim;
7571 input >> dim >> spacedim;
7573 return std::make_pair(dim, spacedim);
7582template <
int dim,
int spacedim>
7584 : default_subdivisions(1)
7590template <
int dim,
int spacedim>
7595 get_dataset_names(),
7596 get_nonscalar_data_ranges(),
7603template <
int dim,
int spacedim>
7608 get_dataset_names(),
7609 get_nonscalar_data_ranges(),
7616template <
int dim,
int spacedim>
7621 get_dataset_names(),
7622 get_nonscalar_data_ranges(),
7629template <
int dim,
int spacedim>
7634 get_dataset_names(),
7635 get_nonscalar_data_ranges(),
7642template <
int dim,
int spacedim>
7647 get_dataset_names(),
7648 get_nonscalar_data_ranges(),
7655template <
int dim,
int spacedim>
7660 get_dataset_names(),
7661 get_nonscalar_data_ranges(),
7668template <
int dim,
int spacedim>
7673 get_dataset_names(),
7674 get_nonscalar_data_ranges(),
7681template <
int dim,
int spacedim>
7686 get_dataset_names(),
7687 get_nonscalar_data_ranges(),
7692template <
int dim,
int spacedim>
7697 get_dataset_names(),
7698 get_nonscalar_data_ranges(),
7703template <
int dim,
int spacedim>
7708 get_dataset_names(),
7709 get_nonscalar_data_ranges(),
7715template <
int dim,
int spacedim>
7718 const std::string &filename,
7721#ifndef DEAL_II_WITH_MPI
7725 std::ofstream f(filename);
7733 int ierr = MPI_Info_create(&info);
7736 ierr = MPI_File_open(
7737 comm, filename.c_str(), MPI_MODE_CREATE | MPI_MODE_WRONLY, info, &fh);
7740 ierr = MPI_File_set_size(fh, 0);
7744 ierr = MPI_Barrier(
comm);
7746 ierr = MPI_Info_free(&info);
7750 unsigned int header_size;
7751 std::uint64_t footer_offset;
7756 std::stringstream ss;
7758 header_size = ss.str().size();
7761 fh, 0, ss.str().c_str(), header_size, MPI_CHAR, MPI_STATUS_IGNORE);
7765 ierr = MPI_Bcast(&header_size, 1, MPI_UNSIGNED, 0,
comm);
7769 const auto & patches = get_patches();
7778 std::stringstream ss;
7779 if (my_n_patches > 0 || (global_n_patches == 0 && myrank == 0))
7781 get_dataset_names(),
7782 get_nonscalar_data_ranges(),
7787 const std::uint64_t size_on_proc = ss.str().size();
7788 std::uint64_t prefix_sum = 0;
7790 MPI_Exscan(&size_on_proc, &prefix_sum, 1, MPI_UINT64_T, MPI_SUM,
comm);
7794 const MPI_Offset offset =
static_cast<MPI_Offset
>(header_size) + prefix_sum;
7804 if (myrank == n_ranks - 1)
7807 footer_offset = size_on_proc + offset;
7809 std::stringstream ss;
7811 const unsigned int footer_size = ss.str().size();
7829 ierr = MPI_File_sync(fh);
7832 ierr = MPI_File_close(&fh);
7839template <
int dim,
int spacedim>
7843 const std::vector<std::string> &piece_names)
const
7847 get_dataset_names(),
7848 get_nonscalar_data_ranges(),
7854template <
int dim,
int spacedim>
7857 const std::string &directory,
7858 const std::string &filename_without_extension,
7859 const unsigned int counter,
7861 const unsigned int n_digits_for_counter,
7862 const unsigned int n_groups)
const
7865 const unsigned int n_ranks =
7867 const unsigned int n_files_written =
7868 (n_groups == 0 || n_groups > n_ranks) ? n_ranks : n_groups;
7873 const unsigned int n_digits =
7876 const unsigned int color = rank % n_files_written;
7877 const std::string filename =
7878 directory + filename_without_extension +
"_" +
7882 if (n_groups == 0 || n_groups > n_ranks)
7885 std::ofstream output(filename.c_str());
7887 this->write_vtu(output);
7889 else if (n_groups == 1)
7892 this->write_vtu_in_parallel(filename.c_str(), mpi_communicator);
7896#ifdef DEAL_II_WITH_MPI
7899 int ierr = MPI_Comm_split(mpi_communicator, color, rank, &comm_group);
7901 this->write_vtu_in_parallel(filename.c_str(), comm_group);
7909 const std::string pvtu_filename =
7910 filename_without_extension +
"_" +
7915 std::vector<std::string> filename_vector;
7916 for (
unsigned int i = 0; i < n_files_written; ++i)
7918 const std::string filename =
7919 filename_without_extension +
"_" +
7923 filename_vector.emplace_back(filename);
7926 std::ofstream pvtu_output((directory + pvtu_filename).c_str());
7927 this->write_pvtu_record(pvtu_output, filename_vector);
7930 return pvtu_filename;
7935template <
int dim,
int spacedim>
7938 std::ostream &out)
const
7941 get_dataset_names(),
7942 get_nonscalar_data_ranges(),
7943 deal_II_intermediate_flags,
7949template <
int dim,
int spacedim>
7952 const std::string & filename,
7958 get_dataset_names(),
7959 get_nonscalar_data_ranges(),
7960 deal_II_intermediate_flags,
7968template <
int dim,
int spacedim>
7972 const std::string & h5_filename,
7973 const double cur_time,
7976 return create_xdmf_entry(
7977 data_filter, h5_filename, h5_filename, cur_time,
comm);
7982template <
int dim,
int spacedim>
7986 const std::string & h5_mesh_filename,
7987 const std::string & h5_solution_filename,
7988 const double cur_time,
7992 ExcMessage(
"XDMF only supports 2 or 3 space dimensions."));
7994#ifndef DEAL_II_WITH_HDF5
7998 (void)h5_mesh_filename;
7999 (void)h5_solution_filename;
8008 std::uint64_t local_node_cell_count[2], global_node_cell_count[2];
8010 local_node_cell_count[0] = data_filter.
n_nodes();
8011 local_node_cell_count[1] = data_filter.
n_cells();
8015 int ierr = MPI_Allreduce(local_node_cell_count,
8016 global_node_cell_count,
8030 const bool have_data = (data_filter.
n_nodes() > 0);
8033 const int key = myrank;
8034 const int color = (have_data ? 1 : 0);
8035 const int ierr = MPI_Comm_split(
comm, color, key, &split_comm);
8039 const bool am_i_first_rank_with_data =
8042 ierr = MPI_Comm_free(&split_comm);
8045 const int tag = 47381;
8048 if (am_i_first_rank_with_data)
8050 const auto &patches = get_patches();
8055 for (
const auto &patch : patches)
8061 h5_solution_filename,
8063 global_node_cell_count[0],
8064 global_node_cell_count[1],
8067 patches[0].reference_cell);
8068 const unsigned int n_data_sets = data_filter.
n_data_sets();
8072 for (
unsigned int i = 0; i < n_data_sets; ++i)
8082 ierr = MPI_Send(buffer.data(), buffer.size(), MPI_BYTE, 0, tag,
comm);
8091 if (myrank == 0 && !am_i_first_rank_with_data)
8096 int ierr = MPI_Probe(MPI_ANY_SOURCE, tag,
comm, &status);
8100 ierr = MPI_Get_count(&status, MPI_BYTE, &len);
8103 std::vector<char> buffer(len);
8104 ierr = MPI_Recv(buffer.data(),
8113 return Utilities::unpack<XDMFEntry>(buffer,
false);
8121template <
int dim,
int spacedim>
8124 const std::vector<XDMFEntry> &entries,
8125 const std::string & filename,
8128#ifdef DEAL_II_WITH_MPI
8132 const int myrank = 0;
8138 std::ofstream xdmf_file(filename);
8140 xdmf_file <<
"<?xml version=\"1.0\" ?>\n";
8141 xdmf_file <<
"<!DOCTYPE Xdmf SYSTEM \"Xdmf.dtd\" []>\n";
8142 xdmf_file <<
"<Xdmf Version=\"2.0\">\n";
8143 xdmf_file <<
" <Domain>\n";
8145 <<
" <Grid Name=\"CellTime\" GridType=\"Collection\" CollectionType=\"Temporal\">\n";
8147 for (
const auto &entry : entries)
8149 xdmf_file << entry.get_xdmf_content(3);
8152 xdmf_file <<
" </Grid>\n";
8153 xdmf_file <<
" </Domain>\n";
8154 xdmf_file <<
"</Xdmf>\n";
8166template <
int dim,
int spacedim>
8172 get_dataset_names(),
8173 get_nonscalar_data_ranges(),
8180#ifdef DEAL_II_WITH_HDF5
8184 template <
int dim,
int spacedim>
8189 const bool write_mesh_file,
8190 const std::string & mesh_filename,
8191 const std::string & solution_filename,
8194 hid_t h5_mesh_file_id = -1, h5_solution_file_id, file_plist_id, plist_id;
8195 hid_t node_dataspace, node_dataset, node_file_dataspace,
8196 node_memory_dataspace, node_dataset_id;
8197 hid_t cell_dataspace, cell_dataset, cell_file_dataspace,
8198 cell_memory_dataspace;
8199 hid_t pt_data_dataspace, pt_data_dataset, pt_data_file_dataspace,
8200 pt_data_memory_dataspace;
8202 std::uint64_t local_node_cell_count[2];
8203 hsize_t count[2], offset[2], node_ds_dim[2], cell_ds_dim[2];
8204 std::vector<double> node_data_vec;
8205 std::vector<unsigned int> cell_data_vec;
8209 local_node_cell_count[0] = data_filter.
n_nodes();
8210 local_node_cell_count[1] = data_filter.
n_cells();
8213 file_plist_id = H5Pcreate(H5P_FILE_ACCESS);
8216# ifdef DEAL_II_WITH_MPI
8217# ifdef H5_HAVE_PARALLEL
8219 status = H5Pset_fapl_mpio(file_plist_id,
comm, MPI_INFO_NULL);
8224# ifndef DEAL_II_WITH_ZLIB
8231 std::uint64_t global_node_cell_count[2] = {0, 0};
8232 std::uint64_t global_node_cell_offsets[2] = {0, 0};
8234# ifdef DEAL_II_WITH_MPI
8235 int ierr = MPI_Allreduce(local_node_cell_count,
8236 global_node_cell_count,
8242 ierr = MPI_Exscan(local_node_cell_count,
8243 global_node_cell_offsets,
8250 global_node_cell_count[0] = local_node_cell_count[0];
8251 global_node_cell_count[1] = local_node_cell_count[1];
8252 global_node_cell_offsets[0] = global_node_cell_offsets[1] = 0;
8256 plist_id = H5Pcreate(H5P_DATASET_XFER);
8258# ifdef DEAL_II_WITH_MPI
8259# ifdef H5_HAVE_PARALLEL
8260 status = H5Pset_dxpl_mpio(plist_id, H5FD_MPIO_COLLECTIVE);
8265 if (write_mesh_file)
8268 h5_mesh_file_id = H5Fcreate(mesh_filename.c_str(),
8276 node_ds_dim[0] = global_node_cell_count[0];
8277 node_ds_dim[1] = (spacedim < 2) ? 2 : spacedim;
8278 node_dataspace = H5Screate_simple(2, node_ds_dim,
nullptr);
8281 cell_ds_dim[0] = global_node_cell_count[1];
8282 cell_ds_dim[1] = patches[0].reference_cell.n_vertices();
8283 cell_dataspace = H5Screate_simple(2, cell_ds_dim,
nullptr);
8287# if H5Gcreate_vers == 1
8288 node_dataset = H5Dcreate(h5_mesh_file_id,
8294 node_dataset_id = H5Pcreate(H5P_DATASET_CREATE);
8295# ifdef DEAL_II_WITH_ZLIB
8296 H5Pset_deflate(node_dataset_id,
8298 H5Pset_chunk(node_dataset_id, 2, node_ds_dim);
8300 node_dataset = H5Dcreate(h5_mesh_file_id,
8307 H5Pclose(node_dataset_id);
8310# if H5Gcreate_vers == 1
8311 cell_dataset = H5Dcreate(h5_mesh_file_id,
8317 node_dataset_id = H5Pcreate(H5P_DATASET_CREATE);
8318# ifdef DEAL_II_WITH_ZLIB
8319 H5Pset_deflate(node_dataset_id,
8321 H5Pset_chunk(node_dataset_id, 2, cell_ds_dim);
8323 cell_dataset = H5Dcreate(h5_mesh_file_id,
8330 H5Pclose(node_dataset_id);
8335 status = H5Sclose(node_dataspace);
8337 status = H5Sclose(cell_dataspace);
8342 count[0] = local_node_cell_count[0];
8343 count[1] = (spacedim < 2) ? 2 : spacedim;
8345 offset[0] = global_node_cell_offsets[0];
8348 node_memory_dataspace = H5Screate_simple(2, count,
nullptr);
8352 node_file_dataspace = H5Dget_space(node_dataset);
8354 status = H5Sselect_hyperslab(
8355 node_file_dataspace, H5S_SELECT_SET, offset,
nullptr, count,
nullptr);
8359 count[0] = local_node_cell_count[1];
8360 count[1] = patches[0].reference_cell.n_vertices();
8361 offset[0] = global_node_cell_offsets[1];
8363 cell_memory_dataspace = H5Screate_simple(2, count,
nullptr);
8366 cell_file_dataspace = H5Dget_space(cell_dataset);
8368 status = H5Sselect_hyperslab(
8369 cell_file_dataspace, H5S_SELECT_SET, offset,
nullptr, count,
nullptr);
8374 status = H5Dwrite(node_dataset,
8376 node_memory_dataspace,
8377 node_file_dataspace,
8379 node_data_vec.data());
8381 node_data_vec.clear();
8384 data_filter.
fill_cell_data(global_node_cell_offsets[0], cell_data_vec);
8385 status = H5Dwrite(cell_dataset,
8387 cell_memory_dataspace,
8388 cell_file_dataspace,
8390 cell_data_vec.data());
8392 cell_data_vec.clear();
8395 status = H5Sclose(node_file_dataspace);
8397 status = H5Sclose(cell_file_dataspace);
8401 status = H5Sclose(node_memory_dataspace);
8403 status = H5Sclose(cell_memory_dataspace);
8407 status = H5Dclose(node_dataset);
8409 status = H5Dclose(cell_dataset);
8413 if (mesh_filename != solution_filename)
8415 status = H5Fclose(h5_mesh_file_id);
8421 if (mesh_filename == solution_filename && write_mesh_file)
8423 h5_solution_file_id = h5_mesh_file_id;
8428 h5_solution_file_id = H5Fcreate(solution_filename.c_str(),
8438 std::string vector_name;
8447 node_ds_dim[0] = global_node_cell_count[0];
8448 node_ds_dim[1] = pt_data_vector_dim;
8449 pt_data_dataspace = H5Screate_simple(2, node_ds_dim,
nullptr);
8452# if H5Gcreate_vers == 1
8453 pt_data_dataset = H5Dcreate(h5_solution_file_id,
8454 vector_name.c_str(),
8459 node_dataset_id = H5Pcreate(H5P_DATASET_CREATE);
8460# ifdef DEAL_II_WITH_ZLIB
8461 H5Pset_deflate(node_dataset_id,
8463 H5Pset_chunk(node_dataset_id, 2, node_ds_dim);
8465 pt_data_dataset = H5Dcreate(h5_solution_file_id,
8466 vector_name.c_str(),
8472 H5Pclose(node_dataset_id);
8477 count[0] = local_node_cell_count[0];
8478 count[1] = pt_data_vector_dim;
8479 offset[0] = global_node_cell_offsets[0];
8481 pt_data_memory_dataspace = H5Screate_simple(2, count,
nullptr);
8485 pt_data_file_dataspace = H5Dget_space(pt_data_dataset);
8487 status = H5Sselect_hyperslab(pt_data_file_dataspace,
8496 status = H5Dwrite(pt_data_dataset,
8498 pt_data_memory_dataspace,
8499 pt_data_file_dataspace,
8505 status = H5Sclose(pt_data_dataspace);
8507 status = H5Sclose(pt_data_memory_dataspace);
8509 status = H5Sclose(pt_data_file_dataspace);
8512 status = H5Dclose(pt_data_dataset);
8517 status = H5Pclose(file_plist_id);
8521 status = H5Pclose(plist_id);
8525 status = H5Fclose(h5_solution_file_id);
8533template <
int dim,
int spacedim>
8537 const std::vector<std::string> & data_names,
8539 std::tuple<
unsigned int,
8543 & nonscalar_data_ranges,
8546 const unsigned int n_data_sets = data_names.size();
8548#ifndef DEAL_II_WITH_MPI
8555 Assert(patches.size() > 0, ExcNoPatches());
8557 if (patches.size() == 0)
8561 unsigned int n_nodes;
8562 std::tie(n_nodes, std::ignore) = count_nodes_and_cells(patches);
8577 [&patches]() {
return create_global_data_table(patches); });
8585 std::move(*create_global_data_table_task.
return_value());
8589 unsigned int i, n_th_vector, data_set, pt_data_vector_dim;
8590 std::string vector_name;
8591 for (n_th_vector = 0, data_set = 0; data_set < n_data_sets;)
8594 while (n_th_vector < nonscalar_data_ranges.size() &&
8595 std::get<0>(nonscalar_data_ranges[n_th_vector]) < data_set)
8599 if (n_th_vector < nonscalar_data_ranges.size() &&
8600 std::get<0>(nonscalar_data_ranges[n_th_vector]) == data_set)
8603 pt_data_vector_dim = std::get<1>(nonscalar_data_ranges[n_th_vector]) -
8604 std::get<0>(nonscalar_data_ranges[n_th_vector]) +
8609 std::get<1>(nonscalar_data_ranges[n_th_vector]) >=
8610 std::get<0>(nonscalar_data_ranges[n_th_vector]),
8611 ExcLowerRange(std::get<1>(nonscalar_data_ranges[n_th_vector]),
8612 std::get<0>(nonscalar_data_ranges[n_th_vector])));
8614 std::get<1>(nonscalar_data_ranges[n_th_vector]) < n_data_sets,
8615 ExcIndexRange(std::get<1>(nonscalar_data_ranges[n_th_vector]),
8621 if (!std::get<2>(nonscalar_data_ranges[n_th_vector]).empty())
8623 vector_name = std::get<2>(nonscalar_data_ranges[n_th_vector]);
8628 for (i = std::get<0>(nonscalar_data_ranges[n_th_vector]);
8629 i < std::get<1>(nonscalar_data_ranges[n_th_vector]);
8631 vector_name += data_names[i] +
"__";
8633 data_names[std::get<1>(nonscalar_data_ranges[n_th_vector])];
8639 pt_data_vector_dim = 1;
8640 vector_name = data_names[data_set];
8650 data_set += pt_data_vector_dim;
8656template <
int dim,
int spacedim>
8660 const std::string & filename,
8664 get_patches(), data_filter, hdf5_flags, filename,
comm);
8669template <
int dim,
int spacedim>
8673 const bool write_mesh_file,
8674 const std::string & mesh_filename,
8675 const std::string & solution_filename,
8689template <
int dim,
int spacedim>
8695 const std::string & filename,
8699 patches, data_filter, flags,
true, filename, filename,
comm);
8704template <
int dim,
int spacedim>
8710 const bool write_mesh_file,
8711 const std::string & mesh_filename,
8712 const std::string & solution_filename,
8718 "DataOutBase was asked to write HDF5 output for a space dimension of 1. "
8719 "HDF5 only supports datasets that live in 2 or 3 dimensions."));
8721#ifndef DEAL_II_WITH_HDF5
8727 (void)write_mesh_file;
8728 (void)mesh_filename;
8729 (void)solution_filename;
8738# ifndef H5_HAVE_PARALLEL
8742 "Serial HDF5 output on multiple processes is not yet supported."));
8752 Assert((patches.size() > 0) || (n_ranks > 1), ExcNoPatches());
8759 const bool have_patches = (patches.size() > 0);
8763 const int color = (have_patches ? 1 : 0);
8764 const int ierr = MPI_Comm_split(
comm, color, key, &split_comm);
8770 do_write_hdf5<dim, spacedim>(patches,
8779 const int ierr = MPI_Comm_free(&split_comm);
8787template <
int dim,
int spacedim>
8795 output_format = default_fmt;
8797 switch (output_format)
8843 write_deal_II_intermediate(out);
8853template <
int dim,
int spacedim>
8862template <
int dim,
int spacedim>
8863template <
typename FlagType>
8869 if (
typeid(flags) ==
typeid(dx_flags))
8871 else if (
typeid(flags) ==
typeid(ucd_flags))
8873 else if (
typeid(flags) ==
typeid(povray_flags))
8875 else if (
typeid(flags) ==
typeid(eps_flags))
8877 else if (
typeid(flags) ==
typeid(gmv_flags))
8879 else if (
typeid(flags) ==
typeid(hdf5_flags))
8881 else if (
typeid(flags) ==
typeid(tecplot_flags))
8884 else if (
typeid(flags) ==
typeid(vtk_flags))
8886 else if (
typeid(flags) ==
typeid(svg_flags))
8888 else if (
typeid(flags) ==
typeid(gnuplot_flags))
8891 else if (
typeid(flags) ==
typeid(deal_II_intermediate_flags))
8892 deal_II_intermediate_flags =
8900template <
int dim,
int spacedim>
8913template <
int dim,
int spacedim>
8920 "A name for the output format to be used");
8924 "Number of subdivisions of each mesh cell");
8970template <
int dim,
int spacedim>
8974 const std::string &output_name = prm.
get(
"Output format");
8976 default_subdivisions = prm.
get_integer(
"Subdivisions");
8979 dx_flags.parse_parameters(prm);
8983 ucd_flags.parse_parameters(prm);
8987 gnuplot_flags.parse_parameters(prm);
8991 povray_flags.parse_parameters(prm);
8995 eps_flags.parse_parameters(prm);
8999 gmv_flags.parse_parameters(prm);
9003 hdf5_flags.parse_parameters(prm);
9007 tecplot_flags.parse_parameters(prm);
9011 vtk_flags.parse_parameters(prm);
9015 deal_II_intermediate_flags.parse_parameters(prm);
9021template <
int dim,
int spacedim>
9025 return (
sizeof(default_fmt) +
9041template <
int dim,
int spacedim>
9043 std::tuple<
unsigned int,
9050 std::tuple<
unsigned int,
9057template <
int dim,
int spacedim>
9065 std::set<std::string> all_names;
9068 std::tuple<
unsigned int,
9072 ranges = this->get_nonscalar_data_ranges();
9073 const std::vector<std::string> data_names = this->get_dataset_names();
9074 const unsigned int n_data_sets = data_names.size();
9075 std::vector<bool> data_set_written(n_data_sets,
false);
9077 for (
const auto &range : ranges)
9079 const std::string &name = std::get<2>(range);
9082 Assert(all_names.find(name) == all_names.end(),
9084 "Error: names of fields in DataOut need to be unique, "
9086 name +
"' is used more than once."));
9087 all_names.insert(name);
9088 for (
unsigned int i = std::get<0>(range); i <= std::get<1>(range);
9090 data_set_written[i] =
true;
9094 for (
unsigned int data_set = 0; data_set < n_data_sets; ++data_set)
9095 if (data_set_written[data_set] ==
false)
9097 const std::string &name = data_names[data_set];
9098 Assert(all_names.find(name) == all_names.end(),
9100 "Error: names of fields in DataOut need to be unique, "
9102 name +
"' is used more than once."));
9103 all_names.insert(name);
9113template <
int dim,
int spacedim>
9121 std::vector<typename ::DataOutBase::Patch<dim, spacedim>> tmp;
9125 std::vector<std::string> tmp;
9126 tmp.swap(dataset_names);
9130 std::tuple<
unsigned int,
9135 tmp.swap(nonscalar_data_ranges);
9142 std::pair<unsigned int, unsigned int> dimension_info =
9145 (dimension_info.second == spacedim),
9146 ExcIncompatibleDimensions(
9147 dimension_info.first, dim, dimension_info.second, spacedim));
9156 getline(in, header);
9158 std::ostringstream s;
9159 s <<
"[deal.II intermediate format graphics data]";
9161 Assert(header == s.str(), ExcUnexpectedInput(s.str(), header));
9165 getline(in, header);
9167 std::ostringstream s;
9171 Assert(header == s.str(), ExcUnexpectedInput(s.str(), header));
9175 getline(in, header);
9177 std::ostringstream s;
9181 Assert(header == s.str(),
9183 "Invalid or incompatible file format. Intermediate format "
9184 "files can only be read by the same deal.II version as they "
9185 "are written by."));
9189 unsigned int n_datasets;
9191 dataset_names.resize(n_datasets);
9192 for (
unsigned int i = 0; i < n_datasets; ++i)
9193 in >> dataset_names[i];
9195 unsigned int n_patches;
9197 patches.resize(n_patches);
9198 for (
unsigned int i = 0; i < n_patches; ++i)
9201 unsigned int n_nonscalar_data_ranges;
9202 in >> n_nonscalar_data_ranges;
9203 nonscalar_data_ranges.resize(n_nonscalar_data_ranges);
9204 for (
unsigned int i = 0; i < n_nonscalar_data_ranges; ++i)
9206 in >> std::get<0>(nonscalar_data_ranges[i]) >>
9207 std::get<1>(nonscalar_data_ranges[i]);
9216 std::get<2>(nonscalar_data_ranges[i]) = name;
9224template <
int dim,
int spacedim>
9230 ParallelIntermediateHeader header;
9231 in.read(
reinterpret_cast<char *
>(&header),
sizeof(header));
9233 header.magic == 0x00dea111,
9235 "Invalid header of parallel deal.II intermediate format encountered."));
9239 "Incorrect header version of parallel deal.II intermediate format."));
9241 std::vector<std::uint64_t> chunk_sizes(header.n_ranks);
9242 in.read(
reinterpret_cast<char *
>(chunk_sizes.data()),
9243 header.n_ranks *
sizeof(std::uint64_t));
9245 for (
unsigned int n = 0; n < header.n_ranks; ++n)
9249 std::vector<char> temp_buffer(chunk_sizes[n]);
9250 in.read(temp_buffer.data(), chunk_sizes[n]);
9253 header.compression) !=
9257 boost::iostreams::filtering_istreambuf f;
9260#ifdef DEAL_II_WITH_ZLIB
9261 f.push(boost::iostreams::zlib_decompressor());
9266 "Decompression requires deal.II to be configured with ZLIB support."));
9269 boost::iostreams::basic_array_source<char> source(temp_buffer.data(),
9270 temp_buffer.size());
9273 std::stringstream datastream;
9274 boost::iostreams::copy(f, datastream);
9284 temp_reader.
read(datastream);
9292template <
int dim,
int spacedim>
9296 using Patch = typename ::DataOutBase::Patch<dim, spacedim>;
9299 const std::vector<Patch> &source_patches = source.
get_patches();
9304 ExcIncompatibleDatasetNames());
9307 Assert(get_nonscalar_data_ranges().size() ==
9309 ExcMessage(
"Both sources need to declare the same components "
9311 for (
unsigned int i = 0; i < get_nonscalar_data_ranges().size(); ++i)
9313 Assert(std::get<0>(get_nonscalar_data_ranges()[i]) ==
9315 ExcMessage(
"Both sources need to declare the same components "
9317 Assert(std::get<1>(get_nonscalar_data_ranges()[i]) ==
9319 ExcMessage(
"Both sources need to declare the same components "
9321 Assert(std::get<2>(get_nonscalar_data_ranges()[i]) ==
9323 ExcMessage(
"Both sources need to declare the same components "
9328 Assert(patches[0].n_subdivisions == source_patches[0].n_subdivisions,
9329 ExcIncompatiblePatchLists());
9330 Assert(patches[0].data.n_rows() == source_patches[0].data.n_rows(),
9331 ExcIncompatiblePatchLists());
9332 Assert(patches[0].data.n_cols() == source_patches[0].data.n_cols(),
9333 ExcIncompatiblePatchLists());
9337 const unsigned int old_n_patches = patches.size();
9338 patches.insert(patches.end(), source_patches.begin(), source_patches.end());
9341 for (
unsigned int i = old_n_patches; i < patches.size(); ++i)
9342 patches[i].patch_index += old_n_patches;
9345 for (
unsigned int i = old_n_patches; i < patches.size(); ++i)
9347 if (patches[i].neighbors[n] !=
9349 patches[i].neighbors[n] += old_n_patches;
9354template <
int dim,
int spacedim>
9355const std::vector<typename ::DataOutBase::Patch<dim, spacedim>> &
9363template <
int dim,
int spacedim>
9364std::vector<std::string>
9367 return dataset_names;
9372template <
int dim,
int spacedim>
9374 std::tuple<
unsigned int,
9380 return nonscalar_data_ranges;
9389 , h5_sol_filename(
"")
9390 , h5_mesh_filename(
"")
9392 , num_nodes(
numbers::invalid_unsigned_int)
9393 , num_cells(
numbers::invalid_unsigned_int)
9394 , dimension(
numbers::invalid_unsigned_int)
9395 , space_dimension(
numbers::invalid_unsigned_int)
9403 const std::uint64_t nodes,
9404 const std::uint64_t cells,
9405 const unsigned int dim)
9411 const std::uint64_t nodes,
9412 const std::uint64_t cells,
9413 const unsigned int dim,
9415 :
XDMFEntry(filename, filename, time, nodes, cells, dim, dim, cell_type)
9421 const std::string & solution_filename,
9423 const std::uint64_t nodes,
9424 const std::uint64_t cells,
9425 const unsigned int dim)
9439 const std::string & solution_filename,
9441 const std::uint64_t nodes,
9442 const std::uint64_t cells,
9443 const unsigned int dim,
9458 const std::string & solution_filename,
9460 const std::uint64_t nodes,
9461 const std::uint64_t cells,
9462 const unsigned int dim,
9463 const unsigned int spacedim)
9484 const unsigned int dimension)
9491 return ReferenceCells::get_hypercube<0>();
9493 return ReferenceCells::get_hypercube<1>();
9495 return ReferenceCells::get_hypercube<2>();
9497 return ReferenceCells::get_hypercube<3>();
9510 const std::string & solution_filename,
9512 const std::uint64_t nodes,
9513 const std::uint64_t cells,
9514 const unsigned int dim,
9515 const unsigned int spacedim,
9518 , h5_sol_filename(solution_filename)
9519 , h5_mesh_filename(mesh_filename)
9524 , space_dimension(spacedim)
9525 , cell_type(cell_type_hex_if_invalid(cell_type_, dim))
9532 const unsigned int dimension)
9545 indent(
const unsigned int indent_level)
9547 std::string res =
"";
9548 for (
unsigned int i = 0; i < indent_level; ++i)
9561 (void)reference_cell;
9563 ExcMessage(
"Incorrect ReferenceCell type passed in."));
9575 std::stringstream ss;
9577 ss << indent(indent_level + 0)
9578 <<
"<Grid Name=\"mesh\" GridType=\"Uniform\">\n";
9579 ss << indent(indent_level + 1) <<
"<Time Value=\"" <<
entry_time <<
"\"/>\n";
9580 ss << indent(indent_level + 1) <<
"<Geometry GeometryType=\""
9582 ss << indent(indent_level + 2) <<
"<DataItem Dimensions=\"" <<
num_nodes
9584 <<
"\" NumberType=\"Float\" Precision=\"8\" Format=\"HDF\">\n";
9586 ss << indent(indent_level + 2) <<
"</DataItem>\n";
9587 ss << indent(indent_level + 1) <<
"</Geometry>\n";
9592 ss << indent(indent_level + 1) <<
"<Topology TopologyType=\"";
9610 ss <<
"Quadrilateral";
9629 ss <<
"Tetrahedron";
9633 ss <<
"\" NumberOfElements=\"" <<
num_cells;
9635 ss <<
"\" NodesPerElement=\"1\">\n";
9637 ss <<
"\" NodesPerElement=\"2\">\n";
9642 ss << indent(indent_level + 2) <<
"<DataItem Dimensions=\"" <<
num_cells
9644 <<
"\" NumberType=\"UInt\" Format=\"HDF\">\n";
9647 ss << indent(indent_level + 2) <<
"</DataItem>\n";
9648 ss << indent(indent_level + 1) <<
"</Topology>\n";
9654 ss << indent(indent_level + 1)
9655 <<
"<Topology TopologyType=\"Polyvertex\" NumberOfElements=\""
9657 ss << indent(indent_level + 1) <<
"</Topology>\n";
9662 ss << indent(indent_level + 1) <<
"<Attribute Name=\""
9663 << attribute_dim.first <<
"\" AttributeType=\""
9664 << (attribute_dim.second > 1 ?
"Vector" :
"Scalar")
9665 <<
"\" Center=\"Node\">\n";
9667 ss << indent(indent_level + 2) <<
"<DataItem Dimensions=\"" <<
num_nodes
9668 <<
" " << (attribute_dim.second > 1 ? 3 : 1)
9669 <<
"\" NumberType=\"Float\" Precision=\"8\" Format=\"HDF\">\n";
9671 << attribute_dim.first <<
'\n';
9672 ss << indent(indent_level + 2) <<
"</DataItem>\n";
9673 ss << indent(indent_level + 1) <<
"</Attribute>\n";
9676 ss << indent(indent_level + 0) <<
"</Grid>\n";
9685 template <
int dim,
int spacedim>
9690 out <<
"[deal.II intermediate Patch<" << dim <<
',' << spacedim <<
">]"
9709 out << patch.
data.n_rows() <<
' ' << patch.
data.n_cols() <<
'\n';
9710 for (
unsigned int i = 0; i < patch.
data.n_rows(); ++i)
9711 for (
unsigned int j = 0; j < patch.
data.n_cols(); ++j)
9712 out << patch.
data[i][j] <<
' ';
9721 template <
int dim,
int spacedim>
9733 getline(in, header);
9734 while ((header.size() != 0) && (header.back() ==
' '))
9735 header.erase(header.size() - 1);
9737 while ((header.empty()) && in);
9739 std::ostringstream s;
9740 s <<
"[deal.II intermediate Patch<" << dim <<
',' << spacedim <<
">]";
9742 Assert(header == s.str(), ExcUnexpectedInput(s.str(), header));
9746#ifdef DEAL_II_HAVE_CXX17
9747 if constexpr (dim > 0)
9768 unsigned int n_subdivisions;
9769 in >> n_subdivisions;
9770#ifdef DEAL_II_HAVE_CXX17
9771 if constexpr (dim > 1)
9778 const_cast<unsigned int &
>(patch.
n_subdivisions) = n_subdivisions;
9783 unsigned int n_rows, n_cols;
9784 in >> n_rows >> n_cols;
9785 patch.
data.reinit(n_rows, n_cols);
9786 for (
unsigned int i = 0; i < patch.
data.n_rows(); ++i)
9787 for (
unsigned int j = 0; j < patch.
data.n_cols(); ++j)
9788 in >> patch.
data[i][j];
9799#include "data_out_base.inst"
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
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
Map3DPoint existing_points
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
void parse_parameters(ParameterHandler &prm)
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_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
std::size_t memory_consumption() const
void set_default_format(const DataOutBase::OutputFormat default_format)
void write(std::ostream &out, const DataOutBase::OutputFormat output_format=DataOutBase::default_format) const
void write_gnuplot(std::ostream &out) const
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
void validate_dataset_names() const
void set_flags(const FlagType &flags)
void write_vtk(std::ostream &out) const
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
void write_eps(std::ostream &out) const
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)
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
long int get_integer(const std::string &entry_string) const
bool get_bool(const std::string &entry_name) const
void declare_entry(const std::string &entry, const std::string &default_value, const Patterns::PatternBase &pattern=Patterns::Anything(), const std::string &documentation="", const bool has_to_be_set=false)
std::string get(const std::string &entry_string) const
double get_double(const std::string &entry_name) const
void enter_subsection(const std::string &subsection)
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
unsigned int n_vertices() const
bool is_hyper_cube() const
unsigned int vtk_linear_type() const
unsigned int vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const
unsigned int vtk_quadratic_type() const
unsigned int vtk_lagrange_type() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
static constexpr TableIndices< rank_ > unrolled_to_component_indices(const unsigned int i)
std::vector< RT > return_values()
internal::return_value< RT >::reference_type return_value()
std::string h5_sol_filename
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::map< std::string, unsigned int > attribute_dims
unsigned int space_dimension
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_PACKAGE_VERSION
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_FALLTHROUGH
#define DEAL_II_PACKAGE_NAME
Point< 2 > projected_vertices[4]
Point< 2 > projected_center
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcFileNotOpen(std::string arg1)
static ::ExceptionBase & ExcNotEnoughSpaceDimensionLabels()
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNoPatches()
#define DeclException2(Exception2, type1, type2, outsequence)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcLowerRange(int arg1, int arg2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNeedsHDF5()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcNotInitialized()
static ::ExceptionBase & ExcInvalidDatasetSize(int arg1, int arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Task< RT > new_task(const std::function< RT()> &function)
LogStream & operator<<(LogStream &log, const T &t)
DataComponentInterpretation
@ component_is_part_of_tensor
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)
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_nodes(const std::vector< Patch< dim, spacedim > > &patches, StreamType &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_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_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_vtu_header(std::ostream &out, const VtkFlags &flags)
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_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 write_data(const std::vector< Patch< dim, spacedim > > &patches, unsigned int n_data_sets, const bool double_precision, StreamType &out)
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_pvd_record(std::ostream &out, const std::vector< std::pair< double, std::string > > ×_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_vtu_footer(std::ostream &out)
void write_cells(const std::vector< Patch< dim, spacedim > > &patches, StreamType &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)
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)
OutputFormat parse_output_format(const std::string &format_name)
std::vector< Point< spacedim > > get_node_positions(const std::vector< Patch< dim, spacedim > > &patches)
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)
std::istream & operator>>(std::istream &in, Patch< dim, spacedim > &patch)
std::string get_output_format_names()
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)
void write_visit_record(std::ostream &out, const std::vector< std::string > &piece_names)
void write_high_order_cells(const std::vector< Patch< dim, spacedim > > &patches, StreamType &out, const bool legacy_format)
std::string default_suffix(const OutputFormat output_format)
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)
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_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_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)
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &Du)
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
constexpr const ReferenceCell Tetrahedron
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Pyramid
constexpr const ReferenceCell Invalid
constexpr const ReferenceCell Triangle
constexpr const ReferenceCell Hexahedron
constexpr const ReferenceCell Vertex
constexpr const ReferenceCell Line
int File_write_at_c(MPI_File fh, MPI_Offset offset, const void *buf, MPI_Count count, MPI_Datatype datatype, MPI_Status *status)
int File_write_at_all_c(MPI_File fh, MPI_Offset offset, const void *buf, MPI_Count count, MPI_Datatype datatype, MPI_Status *status)
T sum(const T &t, const MPI_Comm mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm mpi_communicator)
unsigned int this_mpi_process(const MPI_Comm mpi_communicator)
void free_communicator(MPI_Comm mpi_communicator)
constexpr T pow(const T base, const int iexp)
size_t pack(const T &object, std::vector< char > &dest_buffer, const bool allow_compression=true)
std::string encode_base64(const std::vector< unsigned char > &binary_input)
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
unsigned int needed_digits(const unsigned int max_number)
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
static constexpr double PI
static const unsigned int invalid_unsigned_int
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > cos(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sin(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)
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)
void parse_parameters(const ParameterHandler &prm)
bool filter_duplicate_vertices
DataOutFilterFlags(const bool filter_duplicate_vertices=false, const bool xdmf_hdf5_output=false)
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)
static void declare_parameters(ParameterHandler &prm)
std::size_t memory_consumption() const
ReferenceCell reference_cell
static const unsigned int no_neighbor
bool operator==(const Patch &patch) const
static const unsigned int space_dim
unsigned int n_subdivisions
std::array< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > vertices
void swap(Patch< dim, spacedim > &other_patch)
bool points_are_available
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
bool write_higher_order_cells
DataOutBase::CompressionLevel compression_level
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices()
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices()
bool operator<(const SynchronousIterators< Iterators > &a, const SynchronousIterators< Iterators > &b)