32 #include <boost/algorithm/string.hpp> 33 #include <boost/archive/binary_oarchive.hpp> 35 #ifdef DEAL_II_GMSH_WITH_API 54 const bool write_faces,
55 const bool write_diameter,
56 const bool write_measure,
57 const bool write_all_faces)
58 : write_cells(write_cells)
59 , write_faces(write_faces)
60 , write_diameter(write_diameter)
61 , write_measure(write_measure)
62 , write_all_faces(write_all_faces)
71 "Write the mesh connectivity as DX grid cells");
75 "Write faces of cells. These may be boundary faces " 76 "or all faces between mesh cells, according to " 77 "\"Write all faces\"");
81 "If cells are written, additionally write their" 82 " diameter as data for visualization");
86 "Write the volume of each cell as data");
90 "Write all faces, not only boundary");
105 : write_faces(write_faces)
106 , write_lines(write_lines)
128 : write_preamble(write_preamble)
129 , write_faces(write_faces)
130 , write_lines(write_lines)
155 const unsigned int n_extra_curved_line_points,
156 const bool curved_inner_cells,
157 const bool write_additional_boundary_lines)
158 : write_cell_numbers(write_cell_numbers)
159 , n_extra_curved_line_points(n_extra_curved_line_points)
160 , curved_inner_cells(curved_inner_cells)
161 , write_additional_boundary_lines(write_additional_boundary_lines)
183 const unsigned int size,
184 const double line_width,
185 const bool color_lines_on_user_flag,
186 const unsigned int n_boundary_face_points,
187 const bool color_lines_level)
188 : size_type(size_type)
190 , line_width(line_width)
191 , color_lines_on_user_flag(color_lines_on_user_flag)
192 , n_boundary_face_points(n_boundary_face_points)
193 , color_lines_level(color_lines_level)
203 "Depending on this parameter, either the " 205 "of the eps is scaled to \"Size\"");
209 "Size of the output in points");
213 "Width of the lines drawn in points");
217 "Draw lines with user flag set in different color");
221 "Number of points on boundary edges. " 222 "Increase this beyond 2 to see curved boundaries.");
226 "Draw different colors according to grid level.");
233 if (param.
get(
"Size by") == std::string(
"width"))
235 else if (param.
get(
"Size by") == std::string(
"height"))
247 const unsigned int size,
254 color_lines_on_user_flag,
255 n_boundary_face_points)
273 const unsigned int size,
274 const double line_width,
275 const bool color_lines_on_user_flag,
276 const unsigned int n_boundary_face_points,
277 const bool write_cell_numbers,
278 const bool write_cell_number_level,
279 const bool write_vertex_numbers,
280 const bool color_lines_level)
284 color_lines_on_user_flag,
285 n_boundary_face_points,
287 , write_cell_numbers(write_cell_numbers)
288 , write_cell_number_level(write_cell_number_level)
289 , write_vertex_numbers(write_vertex_numbers)
299 "(2D only) Write cell numbers" 300 " into the centers of cells");
304 "(2D only) if \"Cell number\" is true, write " 305 "numbers in the form level.number");
309 "Write numbers for each vertex");
317 write_cell_numbers = param.
get_bool(
"Cell number");
318 write_cell_number_level = param.
get_bool(
"Level number");
319 write_vertex_numbers = param.
get_bool(
"Vertex number");
325 const unsigned int size,
326 const double line_width,
327 const bool color_lines_on_user_flag,
328 const unsigned int n_boundary_face_points,
329 const double azimut_angle,
330 const double turn_angle)
334 color_lines_on_user_flag,
335 n_boundary_face_points)
336 , azimut_angle(azimut_angle)
337 , turn_angle(turn_angle)
347 "Azimuth of the viw point, that is, the angle " 348 "in the plane from the x-axis.");
352 "Elevation of the view point above the xy-plane.");
360 azimut_angle = 90 - param.
get_double(
"Elevation");
367 : draw_boundary(true)
370 , n_boundary_face_points(0)
376 , boundary_thickness(3)
410 const unsigned int boundary_line_thickness,
413 const int azimuth_angle,
414 const int polar_angle,
416 const bool convert_level_number_to_height,
417 const bool label_level_number,
418 const bool label_cell_index,
419 const bool label_material_id,
420 const bool label_subdomain_id,
421 const bool draw_colorbar,
422 const bool draw_legend,
423 const bool label_boundary_id)
426 , line_thickness(line_thickness)
427 , boundary_line_thickness(boundary_line_thickness)
429 , background(background)
430 , azimuth_angle(azimuth_angle)
431 , polar_angle(polar_angle)
433 , convert_level_number_to_height(convert_level_number_to_height)
434 , level_height_factor(0.3f)
435 , cell_font_scaling(1.f)
436 , label_level_number(label_level_number)
437 , label_cell_index(label_cell_index)
438 , label_material_id(label_material_id)
439 , label_subdomain_id(label_subdomain_id)
440 , label_level_subdomain_id(false)
441 , label_boundary_id(label_boundary_id)
442 , draw_colorbar(draw_colorbar)
443 , draw_legend(draw_legend)
447 : draw_bounding_box(false)
560 switch (output_format)
603 if (format_name ==
"none" || format_name ==
"false")
606 if (format_name ==
"dx")
609 if (format_name ==
"ucd")
612 if (format_name ==
"gnuplot")
615 if (format_name ==
"eps")
618 if (format_name ==
"xfig")
621 if (format_name ==
"msh")
624 if (format_name ==
"svg")
627 if (format_name ==
"mathgl")
630 if (format_name ==
"vtk")
633 if (format_name ==
"vtu")
646 return "none|dx|gnuplot|eps|ucd|xfig|msh|svg|mathgl|vtk|vtu";
779 template <
int dim,
int spacedim>
782 std::ostream & out)
const 798 std::vector<unsigned int> renumber(vertices.size());
801 unsigned int new_number = 0;
802 for (
unsigned int i = 0; i < vertices.size(); ++i)
804 renumber[i] = new_number++;
808 out <<
"object \"vertices\" class array type float rank 1 shape " << dim
809 <<
" items " << n_vertices <<
" data follows" <<
'\n';
811 for (
unsigned int i = 0; i < vertices.size(); ++i)
813 out <<
'\t' << vertices[i] <<
'\n';
820 const unsigned int n_faces =
828 out <<
"object \"cells\" class array type int rank 1 shape " 829 << n_vertices_per_cell <<
" items " <<
n_cells <<
" data follows" 840 out <<
"attribute \"element type\" string \"";
848 <<
"attribute \"ref\" string \"positions\"" <<
'\n' 853 out <<
"object \"material\" class array type int rank 0 items " << n_cells
854 <<
" data follows" <<
'\n';
856 out <<
' ' << cell->material_id();
857 out <<
'\n' <<
"attribute \"dep\" string \"connections\"" <<
'\n' <<
'\n';
859 out <<
"object \"level\" class array type int rank 0 items " << n_cells
860 <<
" data follows" <<
'\n';
862 out <<
' ' << cell->level();
863 out <<
'\n' <<
"attribute \"dep\" string \"connections\"" <<
'\n' <<
'\n';
867 out <<
"object \"measure\" class array type float rank 0 items " 868 << n_cells <<
" data follows" <<
'\n';
870 out <<
'\t' << cell->measure();
872 <<
"attribute \"dep\" string \"connections\"" <<
'\n' 878 out <<
"object \"diameter\" class array type float rank 0 items " 879 << n_cells <<
" data follows" <<
'\n';
881 out <<
'\t' << cell->diameter();
883 <<
"attribute \"dep\" string \"connections\"" <<
'\n' 890 out <<
"object \"faces\" class array type int rank 1 shape " 891 << n_vertices_per_face <<
" items " << n_faces <<
" data follows" 896 for (
const unsigned int f : cell->face_indices())
901 for (
unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_face;
904 << renumber[face->vertex_index(
909 out <<
"attribute \"element type\" string \"";
915 <<
"attribute \"ref\" string \"positions\"" <<
'\n' 921 out <<
"object \"boundary\" class array type int rank 0 items " << n_faces
922 <<
" data follows" <<
'\n';
929 <<
static_cast<std::make_signed<types::boundary_id>::type
>(
930 cell->face(f)->boundary_id());
934 out <<
"attribute \"dep\" string \"connections\"" <<
'\n' <<
'\n';
938 out <<
"object \"face measure\" class array type float rank 0 items " 939 << n_faces <<
" data follows" <<
'\n';
943 out <<
' ' << cell->face(f)->measure();
946 out <<
"attribute \"dep\" string \"connections\"" <<
'\n' <<
'\n';
951 out <<
"object \"face diameter\" class array type float rank 0 items " 952 << n_faces <<
" data follows" <<
'\n';
956 out <<
' ' << cell->face(f)->diameter();
959 out <<
"attribute \"dep\" string \"connections\"" <<
'\n' <<
'\n';
973 out <<
"object \"deal data\" class field" <<
'\n' 974 <<
"component \"positions\" value \"vertices\"" <<
'\n' 975 <<
"component \"connections\" value \"cells\"" <<
'\n';
979 out <<
"object \"cell data\" class field" <<
'\n' 980 <<
"component \"positions\" value \"vertices\"" <<
'\n' 981 <<
"component \"connections\" value \"cells\"" <<
'\n';
982 out <<
"component \"material\" value \"material\"" <<
'\n';
983 out <<
"component \"level\" value \"level\"" <<
'\n';
985 out <<
"component \"measure\" value \"measure\"" <<
'\n';
987 out <<
"component \"diameter\" value \"diameter\"" <<
'\n';
992 out <<
"object \"face data\" class field" <<
'\n' 993 <<
"component \"positions\" value \"vertices\"" <<
'\n' 994 <<
"component \"connections\" value \"faces\"" <<
'\n';
995 out <<
"component \"boundary\" value \"boundary\"" <<
'\n';
997 out <<
"component \"measure\" value \"face measure\"" <<
'\n';
999 out <<
"component \"diameter\" value \"face diameter\"" <<
'\n';
1002 out <<
'\n' <<
"object \"grid data\" class group" <<
'\n';
1004 out <<
"member \"cells\" value \"cell data\"" <<
'\n';
1006 out <<
"member \"faces\" value \"face data\"" <<
'\n';
1007 out <<
"end" <<
'\n';
1018 template <
int dim,
int spacedim>
1021 std::ostream & out)
const 1049 out <<
"@f$NOD" <<
'\n' << n_vertices <<
'\n';
1054 for (
unsigned int i = 0; i < vertices.size(); ++i)
1058 <<
" " << vertices[i];
1059 for (
unsigned int d = spacedim + 1;
d <= 3; ++
d)
1065 out <<
"@f$ENDNOD" <<
'\n' 1104 unsigned int elm_type;
1124 out << cell->active_cell_index() + 1 <<
' ' << elm_type <<
' ' 1125 << cell->material_id() <<
' ' << cell->subdomain_id() <<
' ' 1126 << cell->n_vertices() <<
' ';
1147 out <<
"@f$ENDELM\n";
1157 template <
int dim,
int spacedim>
1160 std::ostream & out)
const 1178 std::time_t time1 = std::time(
nullptr);
1179 std::tm * time = std::localtime(&time1);
1181 <<
"# This file was generated by the deal.II library." <<
'\n' 1182 <<
"# Date = " << time->tm_year + 1900 <<
"/" << time->tm_mon + 1
1183 <<
"/" << time->tm_mday <<
'\n' 1184 <<
"# Time = " << time->tm_hour <<
":" << std::setw(2) << time->tm_min
1185 <<
":" << std::setw(2) << time->tm_sec <<
'\n' 1187 <<
"# For a description of the UCD format see the AVS Developer's guide." 1193 out << n_vertices <<
' ' 1203 for (
unsigned int i = 0; i < vertices.size(); ++i)
1207 <<
" " << vertices[i];
1208 for (
unsigned int d = spacedim + 1;
d <= 3; ++
d)
1217 out << cell->active_cell_index() + 1 <<
' ' << cell->material_id() <<
' ';
1274 template <
int dim,
int spacedim>
1292 const int spacedim = 2;
1298 out <<
"#FIG 3.2\nLandscape\nCenter\nInches" << std::endl
1299 <<
"A4\n100.00\nSingle" 1302 <<
"-3" << std::endl
1303 <<
"# generated by deal.II GridOut class" << std::endl
1304 <<
"# reduce first number to scale up image" << std::endl
1305 <<
"1200 2" << std::endl;
1308 unsigned int colno = 32;
1309 out <<
"0 " << colno++ <<
" #ff0000" << std::endl;
1310 out <<
"0 " << colno++ <<
" #ff8000" << std::endl;
1311 out <<
"0 " << colno++ <<
" #ffd000" << std::endl;
1312 out <<
"0 " << colno++ <<
" #ffff00" << std::endl;
1313 out <<
"0 " << colno++ <<
" #c0ff00" << std::endl;
1314 out <<
"0 " << colno++ <<
" #80ff00" << std::endl;
1315 out <<
"0 " << colno++ <<
" #00f000" << std::endl;
1316 out <<
"0 " << colno++ <<
" #00f0c0" << std::endl;
1317 out <<
"0 " << colno++ <<
" #00f0ff" << std::endl;
1318 out <<
"0 " << colno++ <<
" #00c0ff" << std::endl;
1319 out <<
"0 " << colno++ <<
" #0080ff" << std::endl;
1320 out <<
"0 " << colno++ <<
" #0040ff" << std::endl;
1321 out <<
"0 " << colno++ <<
" #0000c0" << std::endl;
1322 out <<
"0 " << colno++ <<
" #5000ff" << std::endl;
1323 out <<
"0 " << colno++ <<
" #8000ff" << std::endl;
1324 out <<
"0 " << colno++ <<
" #b000ff" << std::endl;
1325 out <<
"0 " << colno++ <<
" #ff00ff" << std::endl;
1326 out <<
"0 " << colno++ <<
" #ff80ff" << std::endl;
1328 for (
unsigned int i = 0; i < 8; ++i)
1329 out <<
"0 " << colno++ <<
" #" << std::hex << 32 * i + 31 << 32 * i + 31
1330 << 32 * i + 31 << std::dec << std::endl;
1332 for (
unsigned int i = 1; i < 16; ++i)
1333 out <<
"0 " << colno++ <<
" #00" << std::hex << 16 * i + 15 << std::dec
1334 <<
"00" << std::endl;
1336 for (
unsigned int i = 1; i < 16; ++i)
1337 out <<
"0 " << colno++ <<
" #" << std::hex << 16 * i + 15 << 16 * i + 15
1338 << std::dec <<
"00" << std::endl;
1340 for (
unsigned int i = 1; i < 16; ++i)
1341 out <<
"0 " << colno++ <<
" #" << std::hex << 16 * i + 15 << std::dec
1342 <<
"0000" << std::endl;
1344 for (
unsigned int i = 1; i < 16; ++i)
1345 out <<
"0 " << colno++ <<
" #" << std::hex << 16 * i + 15 <<
"00" 1346 << 16 * i + 15 << std::dec << std::endl;
1348 for (
unsigned int i = 1; i < 16; ++i)
1349 out <<
"0 " << colno++ <<
" #0000" << std::hex << 16 * i + 15 << std::dec
1352 for (
unsigned int i = 1; i < 16; ++i)
1353 out <<
"0 " << colno++ <<
" #00" << std::hex << 16 * i + 15 << 16 * i + 15
1354 << std::dec << std::endl;
1376 out << cell->material_id() + 32;
1379 out << cell->level() + 8;
1382 out << cell->subdomain_id() + 32;
1385 out << cell->level_subdomain_id() + 32;
1394 (900 + cell->material_id()))
1400 << nv + 1 << std::endl;
1406 for (
unsigned int k = 0; k <= nv; ++k)
1410 for (
unsigned int d = 0; d < static_cast<unsigned int>(dim); ++
d)
1414 out <<
'\t' << ((
d == 0) ? val : -val);
1419 static const unsigned int face_reorder[4] = {2, 1, 3, 0};
1421 for (
const unsigned int f : face_reorder)
1449 for (
unsigned int k = 0;
1450 k < GeometryInfo<dim>::vertices_per_face;
1454 for (
unsigned int d = 0; d < static_cast<unsigned int>(dim);
1460 out <<
'\t' << ((
d == 0) ? val : -val);
1477 #ifdef DEAL_II_GMSH_WITH_API 1478 template <
int dim,
int spacedim>
1481 const std::string & filename)
const 1484 const std::array<int, 8> dealii_to_gmsh_type = {{15, 1, 2, 3, 4, 7, 6, 5}};
1487 const std::array<std::vector<unsigned int>, 8> dealii_to_gmsh = {
1494 {{0, 1, 2, 3, 4, 5}},
1495 {{0, 1, 3, 2, 4, 5, 7, 6}}}};
1500 std::vector<double> coords(3 *
vertices.size());
1501 std::vector<std::size_t> nodes(
vertices.size());
1507 for (
unsigned int d = 0;
d < spacedim; ++
d)
1508 coords[i * 3 +
d] = p[
d];
1520 using IdPair = std::pair<types::material_id, types::manifold_id>;
1521 std::map<IdPair, int> id_pair_to_entity_tag;
1522 std::vector<IdPair> all_pairs;
1524 std::set<IdPair> set_of_pairs;
1527 set_of_pairs.insert({cell->material_id(), cell->manifold_id()});
1528 for (
const auto &f : cell->face_iterators())
1530 (f->boundary_id() != 0 &&
1532 set_of_pairs.insert({f->boundary_id(), f->manifold_id()});
1534 for (
const auto l : cell->line_indices())
1536 const auto &f = cell->line(
l);
1538 (f->boundary_id() != 0 &&
1540 set_of_pairs.insert({f->boundary_id(), f->manifold_id()});
1543 all_pairs = {set_of_pairs.begin(), set_of_pairs.end()};
1546 for (
const auto &p : set_of_pairs)
1547 id_pair_to_entity_tag[p] = entity++;
1550 const auto n_entity_tags = id_pair_to_entity_tag.size();
1553 std::vector<std::vector<std::vector<std::size_t>>> element_ids(
1554 n_entity_tags, std::vector<std::vector<std::size_t>>(8));
1555 std::vector<std::vector<std::vector<std::size_t>>> element_nodes(
1556 n_entity_tags, std::vector<std::vector<std::size_t>>(8));
1559 std::size_t element_id = 1;
1561 const auto add_element = [&](
const auto &element,
const int &entity_tag) {
1562 const auto type = element->reference_cell();
1567 for (
const auto v : element->vertex_indices())
1568 element_nodes[entity_tag - 1][type].emplace_back(
1569 element->vertex_index(dealii_to_gmsh[type][v]) + 1);
1572 element_ids[entity_tag - 1][type].emplace_back(element_id);
1580 std::set<std::pair<int, int>> dim_entity_tag;
1582 auto maybe_add_element =
1583 [&](
const auto & element,
1585 const auto struct_dim = element->structure_dimension;
1589 const bool non_default_boundary_or_material_id =
1590 (boundary_or_material_id != 0 &&
1592 const bool non_default_manifold =
1594 if (struct_dim == dim || non_default_boundary_or_material_id ||
1595 non_default_manifold)
1597 const auto entity_tag =
1598 id_pair_to_entity_tag[{boundary_or_material_id,
manifold_id}];
1599 add_element(element, entity_tag);
1600 dim_entity_tag.insert({struct_dim, entity_tag});
1607 maybe_add_element(cell, cell->material_id());
1608 for (
const auto &face : cell->face_iterators())
1609 maybe_add_element(face, face->boundary_id());
1611 for (
const auto l : cell->line_indices())
1612 maybe_add_element(cell->line(
l), cell->line(
l)->boundary_id());
1617 gmsh::option::setNumber(
"General.Verbosity", 0);
1618 gmsh::model::add(
"Grid generated in deal.II");
1619 for (
const auto &p : dim_entity_tag)
1621 gmsh::model::addDiscreteEntity(p.first, p.second);
1622 gmsh::model::mesh::addNodes(p.first, p.second, nodes, coords);
1625 for (
unsigned int entity_tag = 0; entity_tag < n_entity_tags; ++entity_tag)
1626 for (
unsigned int t = 1; t < 8; ++t)
1628 const auto all_element_ids = element_ids[entity_tag][t];
1629 const auto all_element_nodes = element_nodes[entity_tag][t];
1630 const auto gmsh_t = dealii_to_gmsh_type[t];
1631 if (all_element_ids.size() > 0)
1632 gmsh::model::mesh::addElementsByType(entity_tag + 1,
1639 gmsh::model::mesh::reclassifyNodes();
1640 gmsh::model::mesh::removeDuplicateNodes();
1644 for (
const auto &it : dim_entity_tag)
1646 const auto &
d = it.first;
1647 const auto &entity_tag = it.second;
1648 const auto &
boundary_id = all_pairs[entity_tag - 1].first;
1649 const auto &
manifold_id = all_pairs[entity_tag - 1].second;
1651 std::string physical_name;
1662 std::string sep = physical_name !=
"" ?
", " :
"";
1665 sep +
"ManifoldID:" +
1667 const auto physical_tag =
1668 gmsh::model::addPhysicalGroup(
d, {entity_tag}, -1);
1669 if (physical_name !=
"")
1670 gmsh::model::setPhysicalName(
d, physical_tag, physical_name);
1674 gmsh::write(filename);
1699 const float camera_focus)
1705 camera_focus / ((point - camera_position) * camera_direction);
1708 camera_position + phi * (point - camera_position);
1710 return {(projection - camera_position - camera_focus * camera_direction) *
1712 (projection - camera_position - camera_focus * camera_direction) *
1719 template <
int dim,
int spacedim>
1722 std::ostream & )
const 1725 ExcMessage(
"Mesh output in SVG format is not implemented for anything " 1726 "other than two-dimensional meshes in two-dimensional " 1727 "space. That's because three-dimensional meshes are best " 1728 "viewed in programs that allow changing the viewpoint, " 1729 "but SVG format does not allow this: It is an inherently " 1730 "2d format, and for three-dimensional meshes would " 1731 "require choosing one, fixed viewpoint." 1733 "You probably want to output your mesh in a format such " 1734 "as VTK, VTU, or gnuplot."));
1743 unsigned int min_level, max_level;
1751 Assert(height != 0 || width != 0,
1752 ExcMessage(
"You have to set at least one of width and height"));
1754 unsigned int margin_in_percent = 0;
1756 margin_in_percent = 8;
1759 unsigned int cell_label_font_size;
1772 float x_max_perspective, x_min_perspective;
1773 float y_max_perspective, y_min_perspective;
1775 float x_dimension_perspective, y_dimension_perspective;
1779 double x_min = tria.
begin()->vertex(0)[0];
1780 double x_max = x_min;
1781 double y_min = tria.
begin()->vertex(0)[1];
1782 double y_max = y_min;
1784 double x_dimension, y_dimension;
1786 min_level = max_level = tria.
begin()->level();
1789 std::set<unsigned int> materials;
1792 std::set<unsigned int> levels;
1795 std::set<unsigned int> subdomains;
1798 std::set<int> level_subdomains;
1806 for (
unsigned int vertex_index = 0; vertex_index < cell->n_vertices();
1809 if (cell->vertex(vertex_index)[0] < x_min)
1810 x_min = cell->vertex(vertex_index)[0];
1811 if (cell->vertex(vertex_index)[0] > x_max)
1812 x_max = cell->vertex(vertex_index)[0];
1814 if (cell->vertex(vertex_index)[1] < y_min)
1815 y_min = cell->vertex(vertex_index)[1];
1816 if (cell->vertex(vertex_index)[1] > y_max)
1817 y_max = cell->vertex(vertex_index)[1];
1820 if (static_cast<unsigned int>(cell->level()) < min_level)
1821 min_level = cell->level();
1822 if (static_cast<unsigned int>(cell->level()) > max_level)
1823 max_level = cell->level();
1825 materials.insert(cell->material_id());
1826 levels.insert(cell->level());
1827 if (cell->is_active())
1828 subdomains.insert(cell->subdomain_id() + 2);
1829 level_subdomains.insert(cell->level_subdomain_id() + 2);
1832 x_dimension = x_max - x_min;
1833 y_dimension = y_max - y_min;
1836 const unsigned int n_materials = materials.size();
1839 const unsigned int n_levels = levels.size();
1842 const unsigned int n_subdomains = subdomains.size();
1845 const unsigned int n_level_subdomains = level_subdomains.size();
1859 n = n_level_subdomains;
1868 camera_position[0] = 0;
1869 camera_position[1] = 0;
1870 camera_position[2] = 2. *
std::max(x_dimension, y_dimension);
1873 camera_direction[0] = 0;
1874 camera_direction[1] = 0;
1875 camera_direction[2] = -1;
1878 camera_horizontal[0] = 1;
1879 camera_horizontal[1] = 0;
1880 camera_horizontal[2] = 0;
1882 camera_focus = .5 *
std::max(x_dimension, y_dimension);
1888 const double angle_factor = 3.14159265 / 180.;
1891 camera_position_temp[1] =
1894 camera_position_temp[2] =
1898 camera_direction_temp[1] =
1901 camera_direction_temp[2] =
1905 camera_horizontal_temp[1] =
1908 camera_horizontal_temp[2] =
1912 camera_position[1] = camera_position_temp[1];
1913 camera_position[2] = camera_position_temp[2];
1915 camera_direction[1] = camera_direction_temp[1];
1916 camera_direction[2] = camera_direction_temp[2];
1918 camera_horizontal[1] = camera_horizontal_temp[1];
1919 camera_horizontal[2] = camera_horizontal_temp[2];
1922 camera_position_temp[0] =
1925 camera_position_temp[1] =
1929 camera_direction_temp[0] =
1932 camera_direction_temp[1] =
1936 camera_horizontal_temp[0] =
1939 camera_horizontal_temp[1] =
1943 camera_position[0] = camera_position_temp[0];
1944 camera_position[1] = camera_position_temp[1];
1946 camera_direction[0] = camera_direction_temp[0];
1947 camera_direction[1] = camera_direction_temp[1];
1949 camera_horizontal[0] = camera_horizontal_temp[0];
1950 camera_horizontal[1] = camera_horizontal_temp[1];
1953 camera_position[0] = x_min + .5 * x_dimension;
1954 camera_position[1] = y_min + .5 * y_dimension;
1956 camera_position[0] += 2. *
std::max(x_dimension, y_dimension) *
1959 camera_position[1] -= 2. *
std::max(x_dimension, y_dimension) *
1966 point[0] = tria.
begin()->vertex(0)[0];
1967 point[1] = tria.
begin()->vertex(0)[1];
1970 float min_level_min_vertex_distance = 0;
1975 (
static_cast<float>(tria.
begin()->level()) /
1976 static_cast<float>(n_levels)) *
1977 std::max(x_dimension, y_dimension);
1980 projection_decomposition = svg_project_point(
1981 point, camera_position, camera_direction, camera_horizontal, camera_focus);
1983 x_max_perspective = projection_decomposition[0];
1984 x_min_perspective = projection_decomposition[0];
1986 y_max_perspective = projection_decomposition[1];
1987 y_min_perspective = projection_decomposition[1];
1991 point[0] = cell->vertex(0)[0];
1992 point[1] = cell->vertex(0)[1];
1999 (
static_cast<float>(cell->level()) / static_cast<float>(n_levels)) *
2000 std::max(x_dimension, y_dimension);
2003 projection_decomposition = svg_project_point(point,
2009 if (x_max_perspective < projection_decomposition[0])
2010 x_max_perspective = projection_decomposition[0];
2011 if (x_min_perspective > projection_decomposition[0])
2012 x_min_perspective = projection_decomposition[0];
2014 if (y_max_perspective < projection_decomposition[1])
2015 y_max_perspective = projection_decomposition[1];
2016 if (y_min_perspective > projection_decomposition[1])
2017 y_min_perspective = projection_decomposition[1];
2019 point[0] = cell->vertex(1)[0];
2020 point[1] = cell->vertex(1)[1];
2022 projection_decomposition = svg_project_point(point,
2028 if (x_max_perspective < projection_decomposition[0])
2029 x_max_perspective = projection_decomposition[0];
2030 if (x_min_perspective > projection_decomposition[0])
2031 x_min_perspective = projection_decomposition[0];
2033 if (y_max_perspective < projection_decomposition[1])
2034 y_max_perspective = projection_decomposition[1];
2035 if (y_min_perspective > projection_decomposition[1])
2036 y_min_perspective = projection_decomposition[1];
2038 point[0] = cell->vertex(2)[0];
2039 point[1] = cell->vertex(2)[1];
2041 projection_decomposition = svg_project_point(point,
2047 if (x_max_perspective < projection_decomposition[0])
2048 x_max_perspective = projection_decomposition[0];
2049 if (x_min_perspective > projection_decomposition[0])
2050 x_min_perspective = projection_decomposition[0];
2052 if (y_max_perspective < projection_decomposition[1])
2053 y_max_perspective = projection_decomposition[1];
2054 if (y_min_perspective > projection_decomposition[1])
2055 y_min_perspective = projection_decomposition[1];
2057 if (cell->n_vertices() == 4)
2059 point[0] = cell->vertex(3)[0];
2060 point[1] = cell->vertex(3)[1];
2062 projection_decomposition = svg_project_point(point,
2068 if (x_max_perspective < projection_decomposition[0])
2069 x_max_perspective = projection_decomposition[0];
2070 if (x_min_perspective > projection_decomposition[0])
2071 x_min_perspective = projection_decomposition[0];
2073 if (y_max_perspective < projection_decomposition[1])
2074 y_max_perspective = projection_decomposition[1];
2075 if (y_min_perspective > projection_decomposition[1])
2076 y_min_perspective = projection_decomposition[1];
2079 if (static_cast<unsigned int>(cell->level()) == min_level)
2080 min_level_min_vertex_distance = cell->minimum_vertex_distance();
2083 x_dimension_perspective = x_max_perspective - x_min_perspective;
2084 y_dimension_perspective = y_max_perspective - y_min_perspective;
2088 width =
static_cast<unsigned int>(
2089 .5 + height * (x_dimension_perspective / y_dimension_perspective));
2090 else if (height == 0)
2091 height =
static_cast<unsigned int>(
2092 .5 + width * (y_dimension_perspective / x_dimension_perspective));
2093 unsigned int additional_width = 0;
2095 unsigned int font_size =
2096 static_cast<unsigned int>(.5 + (height / 100.) * 1.75);
2097 cell_label_font_size =
static_cast<unsigned int>(
2099 min_level_min_vertex_distance /
std::min(x_dimension, y_dimension)));
2106 additional_width =
static_cast<unsigned int>(
2111 additional_width =
static_cast<unsigned int>(
2112 .5 + height * .175);
2124 out <<
"<svg width=\"" << width + additional_width <<
"\" height=\"" << height
2125 <<
"\" xmlns=\"http://www.w3.org/2000/svg\" version=\"1.1\">" <<
'\n' 2132 <<
" <linearGradient id=\"background_gradient\" gradientUnits=\"userSpaceOnUse\" x1=\"0\" y1=\"0\" x2=\"0\" y2=\"" 2133 << height <<
"\">" <<
'\n' 2134 <<
" <stop offset=\"0\" style=\"stop-color:white\"/>" <<
'\n' 2135 <<
" <stop offset=\"1\" style=\"stop-color:lightsteelblue\"/>" <<
'\n' 2136 <<
" </linearGradient>" <<
'\n';
2142 out <<
"<!-- internal style sheet -->" <<
'\n' 2143 <<
"<style type=\"text/css\"><![CDATA[" <<
'\n';
2147 out <<
" rect.background{fill:url(#background_gradient)}" <<
'\n';
2149 out <<
" rect.background{fill:white}" <<
'\n';
2151 out <<
" rect.background{fill:none}" <<
'\n';
2154 out <<
" rect{fill:none; stroke:rgb(25,25,25); stroke-width:" 2156 <<
" text{font-family:Helvetica; text-anchor:middle; fill:rgb(25,25,25)}" 2158 <<
" line{stroke:rgb(25,25,25); stroke-width:" 2160 <<
" path{fill:none; stroke:rgb(25,25,25); stroke-width:" 2162 <<
" circle{fill:white; stroke:black; stroke-width:2}" <<
'\n' 2168 unsigned int labeling_index = 0;
2169 auto materials_it = materials.begin();
2170 auto levels_it = levels.begin();
2171 auto subdomains_it = subdomains.begin();
2172 auto level_subdomains_it = level_subdomains.begin();
2174 for (
unsigned int index = 0; index < n; index++)
2179 h = .6 - (index / (n - 1.)) * .6;
2187 unsigned int i =
static_cast<unsigned int>(h * 6);
2189 double f = h * 6 - i;
2196 r = 255, g =
static_cast<unsigned int>(.5 + 255 * t);
2199 r =
static_cast<unsigned int>(.5 + 255 * q), g = 255;
2202 g = 255, b =
static_cast<unsigned int>(.5 + 255 * t);
2205 g =
static_cast<unsigned int>(.5 + 255 * q), b = 255;
2208 r =
static_cast<unsigned int>(.5 + 255 * t), b = 255;
2211 r = 255, b =
static_cast<unsigned int>(.5 + 255 * q);
2220 labeling_index = *materials_it++;
2223 labeling_index = *levels_it++;
2226 labeling_index = *subdomains_it++;
2229 labeling_index = *level_subdomains_it++;
2235 out <<
" path.p" << labeling_index <<
"{fill:rgb(" << r <<
',' << g
2236 <<
',' << b <<
"); " 2237 <<
"stroke:rgb(25,25,25); stroke-width:" 2240 out <<
" path.ps" << labeling_index <<
"{fill:rgb(" 2241 <<
static_cast<unsigned int>(.5 + .75 * r) <<
',' 2242 << static_cast<unsigned int>(.5 + .75 * g) <<
',' 2243 <<
static_cast<unsigned int>(.5 + .75 *
b) <<
"); " 2244 <<
"stroke:rgb(20,20,20); stroke-width:" 2247 out <<
" rect.r" << labeling_index <<
"{fill:rgb(" << r <<
',' << g
2248 <<
',' << b <<
"); " 2249 <<
"stroke:rgb(25,25,25); stroke-width:" 2256 out <<
"]]></style>" <<
'\n' <<
'\n';
2259 out <<
" <rect class=\"background\" width=\"" << width <<
"\" height=\"" 2260 << height <<
"\"/>" <<
'\n';
2264 unsigned int x_offset = 0;
2267 x_offset =
static_cast<unsigned int>(.5 + (height / 100.) *
2268 (margin_in_percent / 2.));
2270 x_offset =
static_cast<unsigned int>(.5 + height * .025);
2273 <<
" <text x=\"" << x_offset <<
"\" y=\"" 2274 <<
static_cast<unsigned int>(.5 + height * .0525) <<
'\"' 2275 <<
" style=\"font-weight:100; fill:lightsteelblue; text-anchor:start; font-family:Courier; font-size:" 2276 <<
static_cast<unsigned int>(.5 + height * .045) <<
"px\">" 2278 <<
"</text>" <<
'\n';
2297 out <<
" <!-- cells -->" <<
'\n';
2299 for (
unsigned int level_index = min_level; level_index <= max_level;
2312 out <<
" class=\"p";
2314 if (!cell->is_active() &&
2321 out << cell->material_id();
2324 out << static_cast<unsigned int>(cell->level());
2327 if (cell->is_active())
2328 out << cell->subdomain_id() + 2;
2333 out << cell->level_subdomain_id() + 2;
2344 point[0] = cell->vertex(0)[0];
2345 point[1] = cell->vertex(0)[1];
2351 (
static_cast<float>(cell->level()) /
2352 static_cast<float>(n_levels)) *
2353 std::max(x_dimension, y_dimension);
2356 projection_decomposition = svg_project_point(point,
2362 out << static_cast<unsigned int>(
2364 ((projection_decomposition[0] - x_min_perspective) /
2365 x_dimension_perspective) *
2366 (width - (width / 100.) * 2. * margin_in_percent) +
2367 ((width / 100.) * margin_in_percent))
2369 <<
static_cast<unsigned int>(
2370 .5 + height - (height / 100.) * margin_in_percent -
2371 ((projection_decomposition[1] - y_min_perspective) /
2372 y_dimension_perspective) *
2373 (height - (height / 100.) * 2. * margin_in_percent));
2377 point[0] = cell->vertex(1)[0];
2378 point[1] = cell->vertex(1)[1];
2380 projection_decomposition = svg_project_point(point,
2386 out << static_cast<unsigned int>(
2388 ((projection_decomposition[0] - x_min_perspective) /
2389 x_dimension_perspective) *
2390 (width - (width / 100.) * 2. * margin_in_percent) +
2391 ((width / 100.) * margin_in_percent))
2393 <<
static_cast<unsigned int>(
2394 .5 + height - (height / 100.) * margin_in_percent -
2395 ((projection_decomposition[1] - y_min_perspective) /
2396 y_dimension_perspective) *
2397 (height - (height / 100.) * 2. * margin_in_percent));
2401 if (cell->n_vertices() == 4)
2403 point[0] = cell->vertex(3)[0];
2404 point[1] = cell->vertex(3)[1];
2406 projection_decomposition = svg_project_point(point,
2412 out << static_cast<unsigned int>(
2414 ((projection_decomposition[0] - x_min_perspective) /
2415 x_dimension_perspective) *
2416 (width - (width / 100.) * 2. * margin_in_percent) +
2417 ((width / 100.) * margin_in_percent))
2419 <<
static_cast<unsigned int>(
2420 .5 + height - (height / 100.) * margin_in_percent -
2421 ((projection_decomposition[1] - y_min_perspective) /
2422 y_dimension_perspective) *
2423 (height - (height / 100.) * 2. * margin_in_percent));
2428 point[0] = cell->vertex(2)[0];
2429 point[1] = cell->vertex(2)[1];
2431 projection_decomposition = svg_project_point(point,
2437 out << static_cast<unsigned int>(
2439 ((projection_decomposition[0] - x_min_perspective) /
2440 x_dimension_perspective) *
2441 (width - (width / 100.) * 2. * margin_in_percent) +
2442 ((width / 100.) * margin_in_percent))
2444 <<
static_cast<unsigned int>(
2445 .5 + height - (height / 100.) * margin_in_percent -
2446 ((projection_decomposition[1] - y_min_perspective) /
2447 y_dimension_perspective) *
2448 (height - (height / 100.) * 2. * margin_in_percent));
2452 point[0] = cell->vertex(0)[0];
2453 point[1] = cell->vertex(0)[1];
2455 projection_decomposition = svg_project_point(point,
2461 out << static_cast<unsigned int>(
2463 ((projection_decomposition[0] - x_min_perspective) /
2464 x_dimension_perspective) *
2465 (width - (width / 100.) * 2. * margin_in_percent) +
2466 ((width / 100.) * margin_in_percent))
2468 <<
static_cast<unsigned int>(
2469 .5 + height - (height / 100.) * margin_in_percent -
2470 ((projection_decomposition[1] - y_min_perspective) /
2471 y_dimension_perspective) *
2472 (height - (height / 100.) * 2. * margin_in_percent));
2474 out <<
"\"/>" <<
'\n';
2481 point[0] = cell->center()[0];
2482 point[1] = cell->center()[1];
2488 (
static_cast<float>(cell->level()) /
2489 static_cast<float>(n_levels)) *
2490 std::max(x_dimension, y_dimension);
2493 const double distance_to_camera =
2494 std::sqrt(std::pow(point[0] - camera_position[0], 2.) +
2495 std::pow(point[1] - camera_position[1], 2.) +
2496 std::pow(point[2] - camera_position[2], 2.));
2497 const double distance_factor =
2498 distance_to_camera / (2. *
std::max(x_dimension, y_dimension));
2500 projection_decomposition = svg_project_point(point,
2506 const unsigned int font_size_this_cell =
2507 static_cast<unsigned int>(
2509 cell_label_font_size *
2510 std::pow(.5, cell->level() - 4. + 3.5 * distance_factor));
2514 <<
static_cast<unsigned int>(
2516 ((projection_decomposition[0] - x_min_perspective) /
2517 x_dimension_perspective) *
2518 (width - (width / 100.) * 2. * margin_in_percent) +
2519 ((width / 100.) * margin_in_percent))
2521 <<
static_cast<unsigned int>(
2522 .5 + height - (height / 100.) * margin_in_percent -
2523 ((projection_decomposition[1] - y_min_perspective) /
2524 y_dimension_perspective) *
2525 (height - (height / 100.) * 2. * margin_in_percent) +
2526 0.5 * font_size_this_cell)
2527 <<
"\" style=\"font-size:" << font_size_this_cell <<
"px\">";
2531 out << cell->level();
2538 out << cell->index();
2547 <<
static_cast<std::make_signed<types::material_id>::type
>(
2548 cell->material_id());
2556 if (cell->is_active())
2558 std::make_signed<types::subdomain_id>::type
>(
2559 cell->subdomain_id());
2572 <<
static_cast<std::make_signed<types::subdomain_id>::type
>(
2573 cell->level_subdomain_id());
2576 out <<
"</text>" <<
'\n';
2583 for (
auto faceIndex : cell->face_indices())
2585 if (cell->at_boundary(faceIndex))
2587 point[0] = cell->face(faceIndex)->vertex(0)[0];
2588 point[1] = cell->face(faceIndex)->vertex(0)[1];
2594 (
static_cast<float>(cell->level()) /
2595 static_cast<float>(n_levels)) *
2596 std::max(x_dimension, y_dimension);
2599 projection_decomposition =
2600 svg_project_point(point,
2606 out <<
" <line x1=\"" 2607 <<
static_cast<unsigned int>(
2609 ((projection_decomposition[0] -
2610 x_min_perspective) /
2611 x_dimension_perspective) *
2613 (width / 100.) * 2. * margin_in_percent) +
2614 ((width / 100.) * margin_in_percent))
2616 <<
static_cast<unsigned int>(
2618 (height / 100.) * margin_in_percent -
2619 ((projection_decomposition[1] -
2620 y_min_perspective) /
2621 y_dimension_perspective) *
2623 (height / 100.) * 2. * margin_in_percent));
2625 point[0] = cell->face(faceIndex)->vertex(1)[0];
2626 point[1] = cell->face(faceIndex)->vertex(1)[1];
2632 (
static_cast<float>(cell->level()) /
2633 static_cast<float>(n_levels)) *
2634 std::max(x_dimension, y_dimension);
2637 projection_decomposition =
2638 svg_project_point(point,
2645 <<
static_cast<unsigned int>(
2647 ((projection_decomposition[0] -
2648 x_min_perspective) /
2649 x_dimension_perspective) *
2651 (width / 100.) * 2. * margin_in_percent) +
2652 ((width / 100.) * margin_in_percent))
2654 <<
static_cast<unsigned int>(
2656 (height / 100.) * margin_in_percent -
2657 ((projection_decomposition[1] -
2658 y_min_perspective) /
2659 y_dimension_perspective) *
2661 (height / 100.) * 2. * margin_in_percent))
2667 const double distance_to_camera = std::sqrt(
2668 std::pow(point[0] - camera_position[0], 2.) +
2669 std::pow(point[1] - camera_position[1], 2.) +
2670 std::pow(point[2] - camera_position[2], 2.));
2671 const double distance_factor =
2672 distance_to_camera /
2673 (2. *
std::max(x_dimension, y_dimension));
2675 const unsigned int font_size_this_edge =
2676 static_cast<unsigned int>(
2677 .5 + .5 * cell_label_font_size *
2679 cell->level() - 4. +
2680 3.5 * distance_factor));
2682 point[0] = cell->face(faceIndex)->center()[0];
2683 point[1] = cell->face(faceIndex)->center()[1];
2689 (
static_cast<float>(cell->level()) /
2690 static_cast<float>(n_levels)) *
2691 std::max(x_dimension, y_dimension);
2694 projection_decomposition =
2695 svg_project_point(point,
2701 const unsigned int xc =
static_cast<unsigned int>(
2703 ((projection_decomposition[0] - x_min_perspective) /
2704 x_dimension_perspective) *
2706 (width / 100.) * 2. * margin_in_percent) +
2707 ((width / 100.) * margin_in_percent));
2708 const unsigned int yc =
static_cast<unsigned int>(
2709 .5 + height - (height / 100.) * margin_in_percent -
2710 ((projection_decomposition[1] - y_min_perspective) /
2711 y_dimension_perspective) *
2713 (height / 100.) * 2. * margin_in_percent));
2715 out <<
" <circle cx=\"" << xc <<
"\" cy=\"" << yc
2716 <<
"\" r=\"" << font_size_this_edge <<
"\" />" 2719 out <<
" <text x=\"" << xc <<
"\" y=\"" << yc
2720 <<
"\" style=\"font-size:" << font_size_this_edge
2721 <<
"px\" dominant-baseline=\"middle\">" 2722 <<
static_cast<int>(
2723 cell->face(faceIndex)->boundary_id())
2724 <<
"</text>" <<
'\n';
2736 out <<
'\n' <<
" <!-- legend -->" <<
'\n';
2738 additional_width = 0;
2740 additional_width =
static_cast<unsigned int>(.5 + (height / 100.) * 2.5);
2748 unsigned int line_offset = 0;
2749 out <<
" <rect x=\"" << width + additional_width <<
"\" y=\"" 2750 <<
static_cast<unsigned int>(.5 + (height / 100.) * margin_in_percent)
2752 << static_cast<unsigned int>(.5 + (height / 100.) *
2753 (40. - margin_in_percent))
2754 <<
"\" height=\"" <<
static_cast<unsigned int>(.5 + height * .215)
2757 out <<
" <text x=\"" 2758 << width + additional_width +
2759 static_cast<unsigned int>(.5 + (height / 100.) * 1.25)
2761 <<
static_cast<unsigned int>(.5 +
2762 (height / 100.) * margin_in_percent +
2763 (++line_offset) * 1.5 * font_size)
2764 <<
"\" style=\"text-anchor:start; font-weight:bold; font-size:" 2765 << font_size <<
"px\">" 2767 <<
"</text>" <<
'\n';
2771 out <<
" <text x=\"" 2772 << width + additional_width +
2773 static_cast<unsigned int>(.5 + (height / 100.) * 2.)
2775 <<
static_cast<unsigned int>(.5 +
2776 (height / 100.) * margin_in_percent +
2777 (++line_offset) * 1.5 * font_size)
2778 <<
"\" style=\"text-anchor:start; font-style:oblique; font-size:" 2779 << font_size <<
"px\">" 2787 out <<
"</text>" <<
'\n';
2792 out <<
" <text x=\"" 2793 << width + additional_width +
2794 static_cast<unsigned int>(.5 + (height / 100.) * 2.)
2796 <<
static_cast<unsigned int>(.5 +
2797 (height / 100.) * margin_in_percent +
2798 (++line_offset) * 1.5 * font_size)
2799 <<
"\" style=\"text-anchor:start; font-style:oblique; font-size:" 2800 << font_size <<
"px\">" 2807 out <<
"</text>" <<
'\n';
2812 out <<
" <text x=\"" 2813 << width + additional_width +
2814 static_cast<unsigned int>(.5 + (height / 100.) * 2.)
2816 <<
static_cast<unsigned int>(.5 +
2817 (height / 100.) * margin_in_percent +
2818 (++line_offset) * 1.5 * font_size)
2819 <<
"\" style=\"text-anchor:start; font-style:oblique; font-size:" 2820 << font_size <<
"px\">" 2827 out <<
"</text>" <<
'\n';
2832 out <<
" <text x= \"" 2833 << width + additional_width +
2834 static_cast<unsigned int>(.5 + (height / 100.) * 2.)
2836 <<
static_cast<unsigned int>(.5 +
2837 (height / 100.) * margin_in_percent +
2838 (++line_offset) * 1.5 * font_size)
2839 <<
"\" style=\"text-anchor:start; font-style:oblique; font-size:" 2840 << font_size <<
"px\">" 2846 out <<
"</text>" <<
'\n';
2851 out <<
" <text x= \"" 2852 << width + additional_width +
2853 static_cast<unsigned int>(.5 + (height / 100.) * 2.)
2855 <<
static_cast<unsigned int>(.5 +
2856 (height / 100.) * margin_in_percent +
2857 (++line_offset) * 1.5 * font_size)
2858 <<
"\" style=\"text-anchor:start; font-style:oblique; font-size:" 2859 << font_size <<
"px\">" 2860 <<
"level_subdomain_id" 2861 <<
"</text>" <<
'\n';
2866 out <<
" <text x=\"" 2867 << width + additional_width +
2868 static_cast<unsigned int>(.5 + (height / 100.) * 1.25)
2870 <<
static_cast<unsigned int>(.5 +
2871 (height / 100.) * margin_in_percent +
2872 (++line_offset) * 1.5 * font_size)
2873 <<
"\" style=\"text-anchor:start; font-weight:bold; font-size:" 2874 << font_size <<
"px\">" 2876 <<
"</text>" <<
'\n';
2878 out <<
" <text x= \"" 2879 << width + additional_width +
2880 static_cast<unsigned int>(.5 + (height / 100.) * 2.)
2882 <<
static_cast<unsigned int>(.5 +
2883 (height / 100.) * margin_in_percent +
2884 (++line_offset) * 1.5 * font_size)
2885 <<
"\" style=\"text-anchor:start; font-style:oblique; font-size:" 2886 << font_size <<
"px\">" 2888 <<
"</text>" <<
'\n';
2896 out <<
" <text x=\"" << width + additional_width <<
"\" y=\"" 2897 <<
static_cast<unsigned int>(
2898 .5 + (height / 100.) * margin_in_percent + 13.75 * font_size)
2899 <<
"\" style=\"text-anchor:start; font-size:" << font_size <<
"px\">" 2908 out <<
'\n' <<
" <!-- colorbar -->" <<
'\n';
2910 out <<
" <text x=\"" << width + additional_width <<
"\" y=\"" 2911 <<
static_cast<unsigned int>(
2912 .5 + (height / 100.) * (margin_in_percent + 29.) -
2914 <<
"\" style=\"text-anchor:start; font-weight:bold; font-size:" 2915 << font_size <<
"px\">";
2920 out <<
"material_id";
2923 out <<
"level_number";
2926 out <<
"subdomain_id";
2929 out <<
"level_subdomain_id";
2935 out <<
"</text>" <<
'\n';
2937 unsigned int element_height =
static_cast<unsigned int>(
2938 ((height / 100.) * (71. - 2. * margin_in_percent)) / n);
2939 unsigned int element_width =
2940 static_cast<unsigned int>(.5 + (height / 100.) * 2.5);
2942 int labeling_index = 0;
2943 auto materials_it = materials.begin();
2944 auto levels_it = levels.begin();
2945 auto subdomains_it = subdomains.begin();
2946 auto level_subdomains_it = level_subdomains.begin();
2948 for (
unsigned int index = 0; index < n; index++)
2953 labeling_index = *materials_it++;
2956 labeling_index = *levels_it++;
2959 labeling_index = *subdomains_it++;
2962 labeling_index = *level_subdomains_it++;
2968 out <<
" <rect class=\"r" << labeling_index <<
"\" x=\"" 2969 << width + additional_width <<
"\" y=\"" 2970 <<
static_cast<unsigned int>(.5 + (height / 100.) *
2971 (margin_in_percent + 29)) +
2972 (n - index - 1) * element_height
2973 <<
"\" width=\"" << element_width <<
"\" height=\"" 2974 << element_height <<
"\"/>" <<
'\n';
2976 out <<
" <text x=\"" 2977 << width + additional_width + 1.5 * element_width <<
"\" y=\"" 2978 <<
static_cast<unsigned int>(.5 + (height / 100.) *
2979 (margin_in_percent + 29)) +
2980 (n - index - 1 + .5) * element_height +
2981 static_cast<unsigned int>(.5 + font_size * .35)
2983 <<
" style=\"text-anchor:start; font-size:" 2984 <<
static_cast<unsigned int>(.5 + font_size) <<
"px";
2986 if (index == 0 || index == n - 1)
2987 out <<
"; font-weight:bold";
2989 out <<
"\">" << labeling_index;
2996 out <<
"</text>" <<
'\n';
3004 out <<
'\n' <<
"</svg>";
3020 template <
int dim,
int spacedim>
3023 std::ostream & out)
const 3030 const std::time_t time1 = std::time(
nullptr);
3031 const std::tm * time = std::localtime(&time1);
3035 <<
"\n# This file was generated by the deal.II library." 3036 <<
"\n# Date = " << time->tm_year + 1900 <<
"/" << std::setfill(
'0')
3037 << std::setw(2) << time->tm_mon + 1 <<
"/" << std::setfill(
'0')
3038 << std::setw(2) << time->tm_mday <<
"\n# Time = " << std::setfill(
'0')
3039 << std::setw(2) << time->tm_hour <<
":" << std::setfill(
'0')
3040 << std::setw(2) << time->tm_min <<
":" << std::setfill(
'0')
3041 << std::setw(2) << time->tm_sec <<
"\n#" 3042 <<
"\n# For a description of the MathGL script format see the MathGL manual. " 3044 <<
"\n# Note: This file is understood by MathGL v2.1 and higher only, and can " 3045 <<
"\n# be quickly viewed in a graphical environment using \'mglview\'. " 3052 const std::string axes =
"xyz";
3068 out <<
"\nsetsize 800 800";
3069 out <<
"\nrotate 0 0";
3072 out <<
"\nsetsize 800 800";
3073 out <<
"\nrotate 60 40";
3082 <<
"\n# Vertex ordering." 3083 <<
"\n# list <vertex order> <vertex indices>" 3093 out <<
"\nlist f 0 1 2 3" 3098 <<
"\nlist f 0 2 4 6 | 1 3 5 7 | 0 4 1 5 | 2 6 3 7 | 0 1 2 3 | 4 5 6 7" 3107 <<
"\n# List of vertices." 3108 <<
"\n# list <id> <vertices>" 3117 for (
unsigned int i = 0; i < dim; ++i)
3124 out <<
"\nlist " << axes[i] << cell->active_cell_index() <<
" ";
3126 out << cell->vertex(j)[i] <<
" ";
3133 <<
"\n# List of cells to quadplot." 3134 <<
"\n# quadplot <vertex order> <id> <style>" 3139 out <<
"\nquadplot f ";
3140 for (
unsigned int j = 0; j < dim; ++j)
3141 out << axes[j] << i <<
" ";
3167 template <
int dim,
int spacedim,
typename ITERATOR,
typename END>
3169 generate_triangulation_patches(
3175 for (; cell !=
end; ++cell)
3182 for (
const unsigned int v : cell->vertex_indices())
3184 patch.
vertices[v] = cell->vertex(v);
3185 patch.
data(0, v) = cell->level();
3187 static_cast<std::make_signed<types::manifold_id>::type
>(
3188 cell->manifold_id());
3190 static_cast<std::make_signed<types::material_id>::type
>(
3191 cell->material_id());
3192 if (cell->is_active())
3194 static_cast<std::make_signed<types::subdomain_id>::type
>(
3195 cell->subdomain_id());
3197 patch.
data(3, v) = -1;
3199 static_cast<std::make_signed<types::subdomain_id>::type
>(
3200 cell->level_subdomain_id());
3202 patches.push_back(patch);
3208 std::vector<std::string>
3209 triangulation_patch_data_names()
3211 std::vector<std::string> v(5);
3216 v[4] =
"level_subdomain";
3223 std::vector<typename Triangulation<3, 3>::active_line_iterator>
3226 std::vector<typename Triangulation<3, 3>::active_line_iterator> res;
3228 std::vector<bool> flags;
3233 for (
const auto l : face->line_indices())
3235 const auto line = face->line(
l);
3236 if (line->user_flag_set() || line->has_children())
3239 line->set_user_flag();
3240 if (line->at_boundary())
3241 res.emplace_back(line);
3252 template <
int dim,
int spacedim>
3253 std::vector<typename Triangulation<dim, spacedim>::active_line_iterator>
3265 std::vector<typename Triangulation<3, 3>::active_line_iterator>
3268 std::vector<typename Triangulation<3, 3>::active_line_iterator> res;
3270 std::vector<bool> flags;
3275 for (
const auto l : face->line_indices())
3277 const auto line = face->line(
l);
3278 if (line->user_flag_set() || line->has_children())
3281 line->set_user_flag();
3283 (line->boundary_id() != 0 &&
3285 res.emplace_back(line);
3295 template <
int dim,
int spacedim>
3296 std::vector<typename Triangulation<dim, spacedim>::active_line_iterator>
3307 template <
int dim,
int spacedim>
3308 std::vector<typename Triangulation<dim, spacedim>::active_face_iterator>
3311 std::vector<typename Triangulation<dim, spacedim>::active_face_iterator>
3318 res.push_back(face);
3329 template <
int dim,
int spacedim>
3330 std::vector<typename Triangulation<dim, spacedim>::active_face_iterator>
3333 std::vector<typename Triangulation<dim, spacedim>::active_face_iterator>
3340 (face->boundary_id() != 0 &&
3342 res.push_back(face);
3350 template <
int dim,
int spacedim>
3353 std::ostream & out)
const 3360 const auto n_vertices = vertices.size();
3362 out <<
"# vtk DataFile Version 3.0\n" 3363 <<
"Triangulation generated with deal.II\n" 3365 <<
"DATASET UNSTRUCTURED_GRID\n" 3366 <<
"POINTS " << n_vertices <<
" double\n";
3369 for (
const auto &v : vertices)
3372 for (
unsigned int d = spacedim + 1;
d <= 3; ++
d)
3378 get_relevant_face_iterators(tria) :
3379 get_boundary_face_iterators(tria);
3381 get_relevant_edge_iterators(tria) :
3382 get_boundary_edge_iterators(tria);
3388 "At least one of the flags (output_cells, output_faces, output_edges) has to be enabled!"));
3405 cells_size += cell->n_vertices() + 1;
3408 for (
const auto &face : faces)
3409 cells_size += face->n_vertices() + 1;
3412 for (
const auto &edge : edges)
3413 cells_size += edge->n_vertices() + 1;
3417 out <<
"\nCELLS " << n_cells <<
' ' << cells_size <<
'\n';
3432 static const std::array<int, 8> table = {{1, 3, 5, 9, 10, 14, 13, 12}};
3438 out << cell->n_vertices();
3439 for (
const unsigned int i : cell->vertex_indices())
3452 out << cell->vertex_index(i);
3455 static const std::array<unsigned int, 5> permutation_table{
3457 out << cell->vertex_index(permutation_table[i]);
3465 for (
const auto &face : faces)
3467 out << face->n_vertices();
3468 constexpr
int face_dim = dim > 1 ? dim - 1 : 1;
3469 for (
const unsigned int i : face->vertex_indices())
3473 face->n_vertices() ?
3480 for (
const auto &edge : edges)
3483 for (
const unsigned int i : edge->vertex_indices())
3484 out <<
' ' << edge->vertex_index(i);
3489 out <<
"\nCELL_TYPES " << n_cells <<
'\n';
3493 out << table[static_cast<int>(cell->reference_cell())] <<
' ';
3498 for (
const auto &face : faces)
3499 out << table[static_cast<int>(face->reference_cell())] <<
' ';
3504 for (
const auto &edge : edges)
3505 out << table[static_cast<int>(edge->reference_cell())] <<
' ';
3507 out <<
"\n\nCELL_DATA " << n_cells <<
'\n' 3508 <<
"SCALARS MaterialID int 1\n" 3509 <<
"LOOKUP_TABLE default\n";
3516 out << static_cast<std::make_signed<types::material_id>::type>(
3517 cell->material_id())
3524 for (
const auto &face : faces)
3526 out << static_cast<std::make_signed<types::boundary_id>::type>(
3527 face->boundary_id())
3534 for (
const auto &edge : edges)
3536 out << static_cast<std::make_signed<types::boundary_id>::type>(
3537 edge->boundary_id())
3542 out <<
"\n\nSCALARS ManifoldID int 1\n" 3543 <<
"LOOKUP_TABLE default\n";
3550 out << static_cast<std::make_signed<types::manifold_id>::type>(
3551 cell->manifold_id())
3558 for (
const auto &face : faces)
3560 out << static_cast<std::make_signed<types::manifold_id>::type>(
3561 face->manifold_id())
3568 for (
const auto &edge : edges)
3570 out << static_cast<std::make_signed<types::manifold_id>::type>(
3571 edge->manifold_id())
3584 template <
int dim,
int spacedim>
3587 std::ostream & out)
const 3595 std::vector<DataOutBase::Patch<dim, spacedim>> patches;
3597 generate_triangulation_patches(patches, tria.
begin_active(), tria.
end());
3602 triangulation_patch_data_names(),
3604 std::tuple<
unsigned int,
3612 out <<
" </UnstructuredGrid>\n";
3613 out <<
"<dealiiData encoding=\"base64\">";
3614 std::stringstream outstring;
3615 boost::archive::binary_oarchive ia(outstring);
3619 out <<
"\n</dealiiData>\n";
3620 out <<
"</VTKFile>\n";
3631 template <
int dim,
int spacedim>
3635 const std::string & filename_without_extension,
3636 const bool view_levels,
3637 const bool include_artificial)
const 3639 std::vector<DataOutBase::Patch<dim, spacedim>> patches;
3640 const unsigned int n_datasets = 4;
3641 std::vector<std::string> data_names;
3642 data_names.emplace_back(
"level");
3643 data_names.emplace_back(
"subdomain");
3644 data_names.emplace_back(
"level_subdomain");
3645 data_names.emplace_back(
"proc_writing");
3653 if (cell->has_children())
3655 if (!include_artificial &&
3659 else if (!include_artificial)
3661 if (cell->has_children() &&
3664 else if (cell->is_active() &&
3665 cell->level_subdomain_id() ==
3675 for (
unsigned int vertex = 0; vertex < n_q_points; ++vertex)
3677 patch.
vertices[vertex] = cell->vertex(vertex);
3678 patch.
data(0, vertex) = cell->level();
3679 if (cell->is_active())
3680 patch.
data(1, vertex) =
static_cast<double>(
3681 static_cast<std::make_signed<types::subdomain_id>::type
>(
3682 cell->subdomain_id()));
3684 patch.
data(1, vertex) = -1.0;
3685 patch.
data(2, vertex) =
static_cast<double>(
3686 static_cast<std::make_signed<types::subdomain_id>::type
>(
3687 cell->level_subdomain_id()));
3693 patches.push_back(patch);
3699 std::string new_file = filename_without_extension +
".vtu";
3703 new_file = filename_without_extension +
".proc" +
3708 if (tr->locally_owned_subdomain() == 0)
3710 std::vector<std::string> filenames;
3715 std::size_t pos = filename_without_extension.find_last_of(
'/');
3716 if (pos == std::string::npos)
3720 const unsigned int n_procs =
3722 for (
unsigned int i = 0; i < n_procs; ++i)
3723 filenames.push_back(filename_without_extension.substr(pos) +
3727 const std::string pvtu_filename =
3728 (filename_without_extension +
".pvtu");
3729 std::ofstream pvtu_output(pvtu_filename.c_str());
3737 data_out.add_data_vector(dummy_vector,
"level");
3738 data_out.add_data_vector(dummy_vector,
"subdomain");
3739 data_out.add_data_vector(dummy_vector,
"level_subdomain");
3740 data_out.add_data_vector(dummy_vector,
"proc_writing");
3742 data_out.build_patches();
3744 data_out.write_pvtu_record(pvtu_output, filenames);
3748 std::ofstream out(new_file.c_str());
3750 std::tuple<
unsigned int,
3812 template <
int dim,
int spacedim>
3817 unsigned int n_faces = 0;
3820 if ((face->at_boundary()) && (face->boundary_id() != 0))
3828 template <
int dim,
int spacedim>
3835 std::vector<bool> line_flags;
3839 .clear_user_flags_line();
3841 unsigned int n_lines = 0;
3844 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++
l)
3845 if (cell->line(
l)->at_boundary() && (cell->line(
l)->boundary_id() != 0) &&
3846 (cell->line(
l)->user_flag_set() ==
false))
3849 cell->line(
l)->set_user_flag();
3864 const unsigned int next_element_index,
3865 std::ostream &)
const 3867 return next_element_index;
3873 const unsigned int next_element_index,
3874 std::ostream &)
const 3876 return next_element_index;
3881 const unsigned int next_element_index,
3882 std::ostream &)
const 3884 return next_element_index;
3890 const unsigned int next_element_index,
3891 std::ostream &)
const 3893 return next_element_index;
3898 const unsigned int next_element_index,
3899 std::ostream &)
const 3901 return next_element_index;
3907 const unsigned int next_element_index,
3908 std::ostream &)
const 3910 return next_element_index;
3916 const unsigned int next_element_index,
3917 std::ostream &)
const 3919 return next_element_index;
3924 const unsigned int next_element_index,
3925 std::ostream &)
const 3927 return next_element_index;
3932 template <
int dim,
int spacedim>
3935 const unsigned int next_element_index,
3936 std::ostream & out)
const 3938 unsigned int current_element_index = next_element_index;
3941 if (face->at_boundary() && (face->boundary_id() != 0))
3943 out << current_element_index <<
' ';
3955 out << static_cast<unsigned int>(face->boundary_id()) <<
' ' 3956 << static_cast<unsigned int>(face->boundary_id()) <<
' ' 3959 for (
unsigned int vertex = 0;
3960 vertex < GeometryInfo<dim>::vertices_per_face;
3963 << face->vertex_index(
3968 ++current_element_index;
3970 return current_element_index;
3974 template <
int dim,
int spacedim>
3977 const unsigned int next_element_index,
3978 std::ostream & out)
const 3980 unsigned int current_element_index = next_element_index;
3985 std::vector<bool> line_flags;
3989 .clear_user_flags_line();
3992 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++
l)
3993 if (cell->line(
l)->at_boundary() && (cell->line(
l)->boundary_id() != 0) &&
3994 (cell->line(
l)->user_flag_set() ==
false))
3996 out << next_element_index <<
" 1 ";
3997 out << static_cast<unsigned int>(cell->line(
l)->boundary_id()) <<
' ' 3998 << static_cast<unsigned int>(cell->line(
l)->boundary_id())
4001 for (
unsigned int vertex = 0; vertex < 2; ++vertex)
4003 << cell->line(
l)->vertex_index(
4011 ++current_element_index;
4012 cell->line(
l)->set_user_flag();
4020 return current_element_index;
4027 const unsigned int next_element_index,
4028 std::ostream &)
const 4030 return next_element_index;
4035 const unsigned int next_element_index,
4036 std::ostream &)
const 4038 return next_element_index;
4043 const unsigned int next_element_index,
4044 std::ostream &)
const 4046 return next_element_index;
4051 const unsigned int next_element_index,
4052 std::ostream &)
const 4054 return next_element_index;
4059 const unsigned int next_element_index,
4060 std::ostream &)
const 4062 return next_element_index;
4068 const unsigned int next_element_index,
4069 std::ostream &)
const 4071 return next_element_index;
4077 const unsigned int next_element_index,
4078 std::ostream &)
const 4080 return next_element_index;
4085 const unsigned int next_element_index,
4086 std::ostream &)
const 4088 return next_element_index;
4093 template <
int dim,
int spacedim>
4096 const unsigned int next_element_index,
4097 std::ostream & out)
const 4099 unsigned int current_element_index = next_element_index;
4103 if (face->at_boundary() && (face->boundary_id() != 0))
4105 out << current_element_index <<
" " 4106 <<
static_cast<unsigned int>(face->boundary_id()) <<
" ";
4119 for (
unsigned int vertex = 0;
4120 vertex < GeometryInfo<dim>::vertices_per_face;
4122 out << face->vertex_index(
4128 ++current_element_index;
4130 return current_element_index;
4135 template <
int dim,
int spacedim>
4138 const unsigned int next_element_index,
4139 std::ostream & out)
const 4141 unsigned int current_element_index = next_element_index;
4146 std::vector<bool> line_flags;
4150 .clear_user_flags_line();
4153 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++
l)
4154 if (cell->line(
l)->at_boundary() && (cell->line(
l)->boundary_id() != 0) &&
4155 (cell->line(
l)->user_flag_set() ==
false))
4157 out << current_element_index <<
" " 4158 <<
static_cast<unsigned int>(cell->line(
l)->boundary_id())
4161 for (
unsigned int vertex = 0; vertex < 2; ++vertex)
4162 out << cell->line(
l)->vertex_index(
4171 ++current_element_index;
4172 cell->line(
l)->set_user_flag();
4179 return current_element_index;
4196 template <
int spacedim>
4200 while (points.size() > 2)
4203 first_difference /= first_difference.
norm();
4205 second_difference /= second_difference.
norm();
4207 if ((first_difference - second_difference).
norm() < 1
e-10)
4208 points.erase(points.begin() + 1);
4216 template <
int spacedim>
4225 for (
const auto &cell : tria.active_cell_iterators())
4228 out <<
"# cell " << cell <<
'\n';
4230 out << cell->vertex(0) <<
' ' << cell->level() <<
' ' 4231 << cell->material_id() <<
'\n' 4232 << cell->vertex(1) <<
' ' << cell->level() <<
' ' 4233 << cell->material_id() <<
'\n' 4246 template <
int spacedim>
4257 const unsigned int n_additional_points =
4259 const unsigned int n_points = 2 + n_additional_points;
4264 std::vector<
Point<dim - 1>> boundary_points;
4265 if (mapping !=
nullptr)
4267 boundary_points.resize(n_points);
4268 boundary_points[0][0] = 0;
4269 boundary_points[n_points - 1][0] = 1;
4270 for (
unsigned int i = 1; i < n_points - 1; ++i)
4271 boundary_points[i](0) = 1. * i / (n_points - 1);
4273 std::vector<double> dummy_weights(n_points, 1. / n_points);
4274 Quadrature<dim - 1> quadrature(boundary_points, dummy_weights);
4280 for (
const auto &cell : tria.active_cell_iterators())
4283 out <<
"# cell " << cell <<
'\n';
4285 if (mapping ==
nullptr ||
4297 << cell->level() <<
' ' << cell->material_id() <<
'\n';
4298 out << cell->vertex(0) <<
' ' << cell->level() <<
' ' 4299 << cell->material_id() <<
'\n' 4307 for (
const unsigned int face_no :
4311 spacedim>::face_iterator
4312 face = cell->face(face_no);
4313 if (dim != spacedim || face->at_boundary() ||
4319 std::vector<Point<spacedim>> line_points;
4322 const unsigned int offset = face_no * n_points;
4323 for (
unsigned int i = 0; i < n_points; ++i)
4324 line_points.push_back(
4326 cell, q_projector.
point(offset + i)));
4327 internal::remove_colinear_points(line_points);
4330 out <<
point <<
' ' << cell->level() <<
' ' 4331 << cell->material_id() <<
'\n';
4333 out <<
'\n' <<
'\n';
4339 out << face->vertex(0) <<
' ' << cell->level() <<
' ' 4340 << cell->material_id() <<
'\n' 4341 << face->vertex(1) <<
' ' << cell->level() <<
' ' 4342 << cell->material_id() <<
'\n' 4358 template <
int spacedim>
4369 const unsigned int n_additional_points =
4371 const unsigned int n_points = 2 + n_additional_points;
4375 std::unique_ptr<Quadrature<dim>> q_projector;
4376 std::vector<Point<1>> boundary_points;
4377 if (mapping !=
nullptr)
4379 boundary_points.resize(n_points);
4380 boundary_points[0][0] = 0;
4381 boundary_points[n_points - 1][0] = 1;
4382 for (
unsigned int i = 1; i < n_points - 1; ++i)
4383 boundary_points[i](0) = 1. * i / (n_points - 1);
4385 std::vector<double> dummy_weights(n_points, 1. / n_points);
4389 QIterated<dim - 1> quadrature(quadrature1d, 1);
4390 q_projector = std::make_unique<Quadrature<dim>>(
4394 for (
const auto &cell : tria.active_cell_iterators())
4397 out <<
"# cell " << cell <<
'\n';
4399 if (mapping ==
nullptr || n_points == 2 ||
4400 (!cell->has_boundary_lines() &&
4404 out << cell->vertex(0) <<
' ' << cell->level() <<
' ' 4405 << cell->material_id() <<
'\n' 4406 << cell->vertex(1) <<
' ' << cell->level() <<
' ' 4407 << cell->material_id() <<
'\n' 4408 << cell->vertex(5) <<
' ' << cell->level() <<
' ' 4409 << cell->material_id() <<
'\n' 4410 << cell->vertex(4) <<
' ' << cell->level() <<
' ' 4411 << cell->material_id() <<
'\n' 4412 << cell->vertex(0) <<
' ' << cell->level() <<
' ' 4413 << cell->material_id() <<
'\n' 4416 out << cell->vertex(2) <<
' ' << cell->level() <<
' ' 4417 << cell->material_id() <<
'\n' 4418 << cell->vertex(3) <<
' ' << cell->level() <<
' ' 4419 << cell->material_id() <<
'\n' 4420 << cell->vertex(7) <<
' ' << cell->level() <<
' ' 4421 << cell->material_id() <<
'\n' 4422 << cell->vertex(6) <<
' ' << cell->level() <<
' ' 4423 << cell->material_id() <<
'\n' 4424 << cell->vertex(2) <<
' ' << cell->level() <<
' ' 4425 << cell->material_id() <<
'\n' 4429 out << cell->vertex(0) <<
' ' << cell->level() <<
' ' 4430 << cell->material_id() <<
'\n' 4431 << cell->vertex(2) <<
' ' << cell->level() <<
' ' 4432 << cell->material_id() <<
'\n' 4434 out << cell->vertex(1) <<
' ' << cell->level() <<
' ' 4435 << cell->material_id() <<
'\n' 4436 << cell->vertex(3) <<
' ' << cell->level() <<
' ' 4437 << cell->material_id() <<
'\n' 4439 out << cell->vertex(5) <<
' ' << cell->level() <<
' ' 4440 << cell->material_id() <<
'\n' 4441 << cell->vertex(7) <<
' ' << cell->level() <<
' ' 4442 << cell->material_id() <<
'\n' 4444 out << cell->vertex(4) <<
' ' << cell->level() <<
' ' 4445 << cell->material_id() <<
'\n' 4446 << cell->vertex(6) <<
' ' << cell->level() <<
' ' 4447 << cell->material_id() <<
'\n' 4452 for (
const unsigned int face_no :
4456 spacedim>::face_iterator
4457 face = cell->face(face_no);
4459 if (face->at_boundary() &&
4462 const unsigned int offset = face_no * n_points * n_points;
4463 for (
unsigned int i = 0; i < n_points - 1; ++i)
4464 for (
unsigned int j = 0; j < n_points - 1; ++j)
4469 q_projector->point(offset + i * n_points + j));
4470 out << p0 <<
' ' << cell->level() <<
' ' 4471 << cell->material_id() <<
'\n';
4475 offset + (i + 1) * n_points + j)))
4476 <<
' ' << cell->level() <<
' ' 4477 << cell->material_id() <<
'\n';
4481 offset + (i + 1) * n_points + j + 1)))
4482 <<
' ' << cell->level() <<
' ' 4483 << cell->material_id() <<
'\n';
4486 q_projector->point(offset + i * n_points +
4488 <<
' ' << cell->level() <<
' ' 4489 << cell->material_id() <<
'\n';
4491 out << p0 <<
' ' << cell->level() <<
' ' 4492 << cell->material_id() <<
'\n';
4493 out <<
'\n' <<
'\n';
4498 for (
unsigned int l = 0;
4499 l < GeometryInfo<dim>::lines_per_face;
4502 const typename ::Triangulation<dim, spacedim>::
4503 line_iterator line = face->line(
l);
4506 &v1 = line->vertex(1);
4507 if (line->at_boundary() ||
4513 std::vector<Point<spacedim>> line_points;
4522 for (
unsigned int i = 0; i < n_points; ++i)
4523 line_points.push_back(
4526 (1 - boundary_points[i][0]) * u0 +
4527 boundary_points[i][0] * u1));
4528 internal::remove_colinear_points(line_points);
4530 out <<
point <<
' ' << cell->level() <<
' ' 4531 <<
static_cast<unsigned int>(
4532 cell->material_id())
4536 out << v0 <<
' ' << cell->level() <<
' ' 4537 << cell->material_id() <<
'\n' 4538 << v1 <<
' ' << cell->level() <<
' ' 4539 << cell->material_id() <<
'\n';
4541 out <<
'\n' <<
'\n';
4558 template <
int dim,
int spacedim>
4582 const unsigned int l)
4602 write_eps(const ::Triangulation<1, 2> &,
4612 write_eps(const ::Triangulation<1, 3> &,
4622 write_eps(const ::Triangulation<2, 3> &,
4633 template <
int dim,
int spacedim>
4635 write_eps(const ::Triangulation<dim, spacedim> &tria,
4641 using LineList = std::list<LineEntry>;
4653 static_cast<const GridOutFlags::EpsFlagsBase &>(eps_flags_3);
4681 for (
const auto &cell : tria.active_cell_iterators())
4682 for (
const unsigned int line_no : cell->line_indices())
4684 typename ::Triangulation<dim, spacedim>::line_iterator
4685 line = cell->line(line_no);
4697 if (!line->has_children() &&
4698 (mapping ==
nullptr || !line->at_boundary()))
4718 line_list.emplace_back(
4719 Point<2>(line->vertex(0)(0), line->vertex(0)(1)),
4720 Point<2>(line->vertex(1)(0), line->vertex(1)(1)),
4721 line->user_flag_set(),
4731 if (mapping !=
nullptr)
4738 std::vector<
Point<dim - 1>> boundary_points(n_points);
4740 for (
unsigned int i = 0; i < n_points; ++i)
4741 boundary_points[i](0) = 1. * (i + 1) / (n_points + 1);
4743 Quadrature<dim - 1> quadrature(boundary_points);
4751 for (
const auto &cell : tria.active_cell_iterators())
4752 for (
const unsigned int face_no :
4755 const typename ::Triangulation<dim, spacedim>::
4756 face_iterator face = cell->face(face_no);
4758 if (face->at_boundary())
4768 const unsigned int offset = face_no * n_points;
4769 for (
unsigned int i = 0; i < n_points; ++i)
4773 cell, q_projector.
point(offset + i)));
4774 const Point<2> p1(p1_dim(0), p1_dim(1));
4776 line_list.emplace_back(p0,
4778 face->user_flag_set(),
4785 const Point<2> p1(p1_dim(0), p1_dim(1));
4786 line_list.emplace_back(p0,
4788 face->user_flag_set(),
4823 const double turn_angle = eps_flags_3.
turn_angle;
4825 -std::sin(z_angle * 2. * pi / 360.) *
4826 std::sin(turn_angle * 2. * pi / 360.),
4827 +std::sin(z_angle * 2. * pi / 360.) *
4828 std::cos(turn_angle * 2. * pi / 360.),
4829 -std::cos(z_angle * 2. * pi / 360.));
4837 ((
Point<dim>(0, 0, 1) * view_direction) * view_direction);
4846 ((
Point<dim>(1, 0, 0) * view_direction) * view_direction) -
4847 ((
Point<dim>(1, 0, 0) * unit_vector1) * unit_vector1));
4851 for (
const auto &cell : tria.active_cell_iterators())
4852 for (
const unsigned int line_no : cell->line_indices())
4854 typename ::Triangulation<dim, spacedim>::line_iterator
4855 line = cell->line(line_no);
4856 line_list.emplace_back(
4857 Point<2>(line->vertex(0) * unit_vector2,
4858 line->vertex(0) * unit_vector1),
4859 Point<2>(line->vertex(1) * unit_vector2,
4860 line->vertex(1) * unit_vector1),
4861 line->user_flag_set(),
4877 double x_min = tria.begin_active()->vertex(0)(0);
4878 double x_max = x_min;
4879 double y_min = tria.begin_active()->vertex(0)(1);
4880 double y_max = y_min;
4881 unsigned int max_level = line_list.begin()->level;
4883 for (LineList::const_iterator line = line_list.begin();
4884 line != line_list.end();
4887 x_min =
std::min(x_min, line->first(0));
4888 x_min =
std::min(x_min, line->second(0));
4890 x_max =
std::max(x_max, line->first(0));
4891 x_max =
std::max(x_max, line->second(0));
4893 y_min =
std::min(y_min, line->first(1));
4894 y_min =
std::min(y_min, line->second(1));
4896 y_max =
std::max(y_max, line->first(1));
4897 y_max =
std::max(y_max, line->second(1));
4899 max_level =
std::max(max_level, line->level);
4907 const double scale =
4908 (eps_flags_base.
size /
4919 std::time_t time1 = std::time(
nullptr);
4920 std::tm * time = std::localtime(&time1);
4921 out <<
"%!PS-Adobe-2.0 EPSF-1.2" <<
'\n' 4922 <<
"%%Title: deal.II Output" <<
'\n' 4923 <<
"%%Creator: the deal.II library" <<
'\n' 4924 <<
"%%Creation Date: " << time->tm_year + 1900 <<
"/" 4925 << time->tm_mon + 1 <<
"/" << time->tm_mday <<
" - " 4926 << time->tm_hour <<
":" << std::setw(2) << time->tm_min <<
":" 4927 << std::setw(2) << time->tm_sec <<
'\n' 4928 <<
"%%BoundingBox: " 4932 <<
static_cast<unsigned int>(
4935 << static_cast<unsigned int>(
4945 out <<
"/m {moveto} bind def" <<
'\n' 4946 <<
"/x {lineto stroke} bind def" <<
'\n' 4947 <<
"/b {0 0 0 setrgbcolor} def" <<
'\n' 4948 <<
"/r {1 0 0 setrgbcolor} def" <<
'\n';
4955 out <<
"/l { neg " << (max_level) <<
" add " 4956 << (0.66666 /
std::max(1
U, (max_level - 1)))
4957 <<
" mul 1 0.8 sethsbcolor} def" <<
'\n';
4971 << (
"/R {rmoveto} bind def\n" 4972 "/Symbol-Oblique /Symbol findfont [1 0 .167 1 0 0] makefont\n" 4973 "dup length dict begin {1 index /FID eq {pop pop} {def} ifelse} forall\n" 4974 "currentdict end definefont\n" 4975 "/MFshow {{dup dup 0 get findfont exch 1 get scalefont setfont\n" 4976 "[ currentpoint ] exch dup 2 get 0 exch rmoveto dup dup 5 get exch 4 get\n" 4977 "{show} {stringwidth pop 0 rmoveto}ifelse dup 3 get\n" 4978 "{2 get neg 0 exch rmoveto pop} {pop aload pop moveto}ifelse} forall} bind def\n" 4979 "/MFwidth {0 exch {dup 3 get{dup dup 0 get findfont exch 1 get scalefont setfont\n" 4980 "5 get stringwidth pop add}\n" 4981 "{pop} ifelse} forall} bind def\n" 4982 "/MCshow { currentpoint stroke m\n" 4983 "exch dup MFwidth -2 div 3 -1 roll R MFshow } def\n")
4987 out <<
"%%EndProlog" <<
'\n' <<
'\n';
4990 out << eps_flags_base.
line_width <<
" setlinewidth" <<
'\n';
4994 const Point<2> offset(x_min, y_min);
4996 for (LineList::const_iterator line = line_list.begin();
4997 line != line_list.end();
5000 out << line->level <<
" l " << (line->first - offset) * scale <<
" m " 5001 << (line->second - offset) * scale <<
" x" <<
'\n';
5006 << (line->first - offset) * scale <<
" m " 5007 << (line->second - offset) * scale <<
" x" <<
'\n';
5013 out <<
"(Helvetica) findfont 140 scalefont setfont" <<
'\n';
5015 for (
const auto &cell : tria.active_cell_iterators())
5017 out << (cell->center()(0) - offset(0)) * scale <<
' ' 5018 << (cell->center()(1) - offset(1)) * scale <<
" m" <<
'\n' 5019 <<
"[ [(Helvetica) 12.0 0.0 true true (";
5023 out << cell->index();
5026 <<
"] -6 MCshow" <<
'\n';
5033 out <<
"(Helvetica) findfont 140 scalefont setfont" <<
'\n';
5039 std::set<unsigned int> treated_vertices;
5040 for (
const auto &cell : tria.active_cell_iterators())
5041 for (
const unsigned int vertex_no : cell->vertex_indices())
5042 if (treated_vertices.find(cell->vertex_index(vertex_no)) ==
5043 treated_vertices.end())
5045 treated_vertices.insert(cell->vertex_index(vertex_no));
5047 out << (cell->vertex(vertex_no)(0) - offset(0)) * scale <<
' ' 5048 << (cell->vertex(vertex_no)(1) - offset(1)) * scale
5050 <<
"[ [(Helvetica) 10.0 0.0 true true (" 5051 << cell->vertex_index(vertex_no) <<
")] " 5052 <<
"] -6 MCshow" <<
'\n';
5056 out <<
"showpage" <<
'\n';
5068 template <
int dim,
int spacedim>
5078 template <
int dim,
int spacedim>
5085 switch (output_format)
5135 template <
int dim,
int spacedim>
5146 #include "grid_out.inst" void parse_parameters(ParameterHandler ¶m)
void parse_parameters(ParameterHandler ¶m)
std::string encode_base64(const std::vector< unsigned char > &binary_input)
ReferenceCell reference_cell
unsigned int n_boundary_faces(const Triangulation< dim, spacedim > &tria) const
unsigned int n_active_cells() const
long int get_integer(const std::string &entry_string) const
unsigned int n_used_vertices() const
constexpr const ReferenceCell Pyramid
const types::manifold_id flat_manifold_id
static const unsigned int invalid_unsigned_int
void write_vtu_footer(std::ostream &out)
OutputFormat default_format
static void declare_parameters(ParameterHandler &prm)
DX(const bool write_cells=true, const bool write_faces=false, const bool write_diameter=false, const bool write_measure=false, const bool write_all_faces=true)
static OutputFormat parse_output_format(const std::string &format_name)
bool serialize_triangulation
void parse_parameters(ParameterHandler ¶m)
types::global_dof_index size_type
void attach_triangulation(const Triangulation< DoFHandlerType::dimension, DoFHandlerType::space_dimension > &)
bool write_additional_boundary_lines
void write_xfig(const Triangulation< dim, spacedim > &tria, std::ostream &out, const Mapping< dim, spacedim > *mapping=nullptr) const
static ::ExceptionBase & ExcIO()
unsigned int n_boundary_face_points
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
static Quadrature< dim > project_to_all_faces(const Quadrature< dim - 1 > &quadrature)
static void declare_parameters(ParameterHandler ¶m)
GridOutFlags::Eps< 2 > eps_flags_2
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)
IteratorRange< active_cell_iterator > active_cell_iterators() const
void write_vtu(const Triangulation< dim, spacedim > &tria, std::ostream &out) const
std::string get(const std::string &entry_string) const
constexpr const ReferenceCell Triangle
void parse_parameters(ParameterHandler ¶m)
void write_mesh_per_processor_as_vtu(const Triangulation< dim, spacedim > &tria, const std::string &filename_without_extension, const bool view_levels=false, const bool include_artificial=false) const
void write_vtu_header(std::ostream &out, const VtkFlags &flags)
const Point< dim > & point(const unsigned int i) const
static void declare_parameters(ParameterHandler ¶m)
Convert the level number into the cell color.
write() calls write_eps()
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim >>> &Du)
active_cell_iterator begin_active(const unsigned int level=0) const
unsigned int line_thickness
bool convert_level_number_to_height
#define AssertThrow(cond, exc)
void write_mathgl(const Triangulation< dim, spacedim > &tria, std::ostream &out) const
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const =0
unsigned int boundary_line_thickness
Point< spacedim > vertices[GeometryInfo< dim >::vertices_per_cell]
cell_iterator begin(const unsigned int level=0) const
constexpr Tensor< 1, dim, typename ProductType< Number1, Number2 >::type > cross_product_3d(const Tensor< 1, dim, Number1 > &src1, const Tensor< 1, dim, Number2 > &src2)
unsigned int write_msh_faces(const Triangulation< dim, spacedim > &tria, const unsigned int next_element_index, std::ostream &out) const
void write_eps(const std::vector< Patch< 2, spacedim >> &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation >> &nonscalar_data_ranges, const EpsFlags &flags, std::ostream &out)
void parse_parameters(ParameterHandler ¶m)
enum GridOutFlags::XFig::Coloring color_by
void parse_parameters(ParameterHandler ¶m)
unsigned int n_extra_curved_line_points
constexpr const ReferenceCell Line
write() calls write_ucd()
constexpr const ReferenceCell Wedge
bool output_only_relevant
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_gnuplot(const Triangulation< dim, spacedim > &tria, std::ostream &out, const Mapping< dim, spacedim > *mapping=nullptr) const
cell_iterator end() const
unsigned int n_boundary_lines(const Triangulation< dim, spacedim > &tria) const
void enter_subsection(const std::string &subsection)
GridOutFlags::Vtk vtk_flags
GridOutFlags::Gnuplot gnuplot_flags
void write_svg(const Triangulation< 2, 2 > &tria, std::ostream &out) const
void write_cells(const std::vector< Patch< dim, spacedim >> &patches, StreamType &out)
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
Convert the global subdomain id into the cell color.
Convert the material id into the cell color.
static ::ExceptionBase & ExcMessage(std::string arg1)
static void declare_parameters(ParameterHandler ¶m)
std::size_t memory_consumption() const
std::string compress(const std::string &input)
double get_double(const std::string &entry_name) const
GridOutFlags::MathGL mathgl_flags
void save_user_flags_line(std::ostream &out) 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)
write() calls write_mathgl()
constexpr const ReferenceCell Tetrahedron
GridOutFlags::DX dx_flags
GridOutFlags::Msh msh_flags
GridOutFlags::Eps< 1 > eps_flags_1
void reinit(const TableIndices< N > &new_size, const bool omit_default_initialization=false)
write() calls write_gnuplot()
#define Assert(cond, exc)
bool write_cell_number_level
void parse_parameters(const ParameterHandler &prm)
EpsFlagsBase(const SizeType size_type=width, const unsigned int size=300, const double line_width=0.5, const bool color_lines_on_user_flag=false, const unsigned int n_boundary_face_points=2, const bool color_lines_level=false)
void reference_cell(const ReferenceCell &reference_cell, Triangulation< dim, spacedim > &tria)
Abstract base class for mapping classes.
const types::boundary_id invalid_boundary_id
GridOutFlags::XFig xfig_flags
Ucd(const bool write_preamble=false, const bool write_faces=false, const bool write_lines=false)
const std::vector< Point< spacedim > > & get_vertices() const
#define DEAL_II_NAMESPACE_CLOSE
static void declare_parameters(ParameterHandler ¶m)
VectorType::value_type * end(VectorType &V)
constexpr const ReferenceCell Hexahedron
bool label_level_subdomain_id
static std::string get_output_format_names()
Convert the subdomain id into the cell color.
Msh(const bool write_faces=false, const bool write_lines=false)
unsigned int n_subdivisions
void save(Archive &ar, const unsigned int version) const
bool get_bool(const std::string &entry_name) const
void write_vtk(const Triangulation< dim, spacedim > &tria, std::ostream &out) const
IteratorRange< active_face_iterator > active_face_iterators() const
GridOutFlags::Ucd ucd_flags
static void declare_parameters(ParameterHandler ¶m)
TriangulationBase< dim, spacedim > Triangulation
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static void declare_parameters(ParameterHandler ¶m)
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Convert the level subdomain id into the cell color.
write() calls write_xfig()
DataComponentInterpretation
bool write_vertex_numbers
Convert the level subdomain id into the cell color.
const types::subdomain_id artificial_subdomain_id
Convert the level into the cell color.
void parse_parameters(ParameterHandler ¶m)
unsigned int n_boundary_face_points
static void declare_parameters(ParameterHandler ¶m)
unsigned int write_ucd_faces(const Triangulation< dim, spacedim > &tria, const unsigned int next_element_index, std::ostream &out) const
write() calls write_msh()
static constexpr double PI
const std::vector< bool > & get_used_vertices() const
#define DEAL_II_NAMESPACE_OPEN
void write_ucd(const Triangulation< dim, spacedim > &tria, std::ostream &out) const
T min(const T &t, const MPI_Comm &mpi_communicator)
constexpr const ReferenceCell Quadrilateral
unsigned int write_ucd_lines(const Triangulation< dim, spacedim > &tria, const unsigned int next_element_index, std::ostream &out) const
write() calls write_svg()
unsigned int write_msh_lines(const Triangulation< dim, spacedim > &tria, const unsigned int next_element_index, std::ostream &out) const
GridOutFlags::Svg svg_flags
bool color_lines_on_user_flag
void parse_parameters(ParameterHandler ¶m)
constexpr const ReferenceCell Vertex
std::array< unsigned int, GeometryInfo< dim >::faces_per_cell > neighbors
void write_eps(const Triangulation< dim, spacedim > &tria, std::ostream &out, const Mapping< dim, spacedim > *mapping=nullptr) const
static ::ExceptionBase & ExcNotImplemented()
Convert the material id into the cell color (default)
IteratorRange< cell_iterator > cell_iterators() const
float level_height_factor
GridOutFlags::Eps< 3 > eps_flags_3
static ::ExceptionBase & ExcInvalidState()
const types::boundary_id internal_face_boundary_id
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)
write() calls write_vtu()
numbers::NumberTraits< Number >::real_type norm() const
bool points_are_available
static void declare_parameters(ParameterHandler ¶m)
void set_flags(const GridOutFlags::DX &flags)
Svg(const unsigned int line_thickness=2, const unsigned int boundary_line_thickness=4, const bool margin=true, const Background background=white, const int azimuth_angle=0, const int polar_angle=0, const Coloring coloring=level_number, const bool convert_level_number_to_height=false, const bool label_level_number=false, const bool label_cell_index=false, const bool label_material_id=false, const bool label_subdomain_id=false, const bool draw_colorbar=false, const bool draw_legend=false, const bool label_boundary_id=false)
void parse_parameters(ParameterHandler ¶m)
void write_dx(const Triangulation< dim, spacedim > &tria, std::ostream &out) const
void write(const Triangulation< dim, spacedim > &tria, std::ostream &out, const OutputFormat output_format, const Mapping< dim, spacedim > *mapping=nullptr) const
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const =0
GridOutFlags::Vtu vtu_flags
T max(const T &t, const MPI_Comm &mpi_communicator)
write() calls write_vtk()
virtual types::subdomain_id locally_owned_subdomain() const
void parse_parameters(ParameterHandler ¶m)
void parse_parameters(ParameterHandler ¶m)
Gnuplot(const bool write_cell_number=false, const unsigned int n_extra_curved_line_points=2, const bool curved_inner_cells=false, const bool write_additional_boundary_lines=true)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()
void write_msh(const Triangulation< dim, spacedim > &tria, std::ostream &out) const
Expression floor(const Expression &x)
std::string default_suffix() const