From ef91a25e2ef96bcb8e3d2cfb5fed012a6910c7e4 Mon Sep 17 00:00:00 2001 From: Alexander Grayver Date: Thu, 26 Jul 2018 17:48:11 +0200 Subject: [PATCH] Add changelog --- .../changes/minor/20180727AlexanderGrayver | 4 + include/deal.II/base/data_out_base.h | 6 +- source/base/data_out_base.cc | 169 ++++++++++-------- 3 files changed, 102 insertions(+), 77 deletions(-) create mode 100644 doc/news/changes/minor/20180727AlexanderGrayver diff --git a/doc/news/changes/minor/20180727AlexanderGrayver b/doc/news/changes/minor/20180727AlexanderGrayver new file mode 100644 index 0000000000..41e5f75e84 --- /dev/null +++ b/doc/news/changes/minor/20180727AlexanderGrayver @@ -0,0 +1,4 @@ +Added support for high-order VTK output by using newly +introduced Lagrange VTK cells. +
+(Alexander Grayver, 2018/07/27) diff --git a/include/deal.II/base/data_out_base.h b/include/deal.II/base/data_out_base.h index 4177a1a1f4..c448e9185e 100644 --- a/include/deal.II/base/data_out_base.h +++ b/include/deal.II/base/data_out_base.h @@ -1162,9 +1162,9 @@ namespace DataOutBase VtkFlags( const double time = std::numeric_limits::min(), const unsigned int cycle = std::numeric_limits::min(), - const bool print_date_and_time = true, - const ZlibCompressionLevel compression_level = best_compression, - const bool write_lagrange_cells = false); + const bool print_date_and_time = true, + const ZlibCompressionLevel compression_level = best_compression, + const bool write_lagrange_cells = false); }; diff --git a/source/base/data_out_base.cc b/source/base/data_out_base.cc index 0ae8dc3554..daa7f29e5c 100644 --- a/source/base/data_out_base.cc +++ b/source/base/data_out_base.cc @@ -692,19 +692,20 @@ namespace // DataOutBase in general Wolfgang // supposed that we don't need it in general, but however this choice avoids a // -Warray-bounds check warning - const unsigned int vtk_cell_type[5] = {1, // VTK_VERTEX - 3, // VTK_LINE - 9, // VTK_QUAD - 12, // VTK_HEXAHEDRON + const unsigned int vtk_cell_type[5] = {1, // VTK_VERTEX + 3, // VTK_LINE + 9, // VTK_QUAD + 12, // VTK_HEXAHEDRON static_cast(-1)}; // VTK cell ids defined in vtk_cell_type are used for linear cells, // the ones defined below are used when Lagrange cells are written. - const unsigned int vtk_lagrange_cell_type[5] = {1, // VTK_VERTEX - 68, // VTK_LAGRANGE_CURVE - 70, // VTK_LAGRANGE_QUADRILATERAL - 72, // VTK_LAGRANGE_HEXAHEDRON - static_cast(-1)}; + const unsigned int vtk_lagrange_cell_type[5] = { + 1, // VTK_VERTEX + 68, // VTK_LAGRANGE_CURVE + 70, // VTK_LAGRANGE_QUADRILATERAL + 72, // VTK_LAGRANGE_HEXAHEDRON + static_cast(-1)}; //----------------------------------------------------------------------// // Auxiliary functions @@ -797,8 +798,11 @@ namespace * Modified from * https://github.com/Kitware/VTK/blob/265ca48a79a36538c95622c237da11133608bbe5/Common/DataModel/vtkLagrangeQuadrilateral.cxx#L558 */ - int vtk_point_index_from_ijk(int i, int j, int, - const std::array &order) + int + vtk_point_index_from_ijk(int i, + int j, + int, + const std::array &order) { bool ibdy = (i == 0 || i == order[0]); bool jbdy = (j == 0 || j == order[1]); @@ -806,23 +810,25 @@ namespace int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0); if (nbdy == 2) // Vertex DOF - { // ijk is a corner node. Return the proper index (somewhere in [0,7]): - return (i ? (j ? 2 : 1) : (j ? 3 : 0)); - } + { // ijk is a corner node. Return the proper index (somewhere in [0,7]): + return (i ? (j ? 2 : 1) : (j ? 3 : 0)); + } int offset = 4; if (nbdy == 1) // Edge DOF - { - if (!ibdy) - { // On i axis - return (i - 1) + (j ? order[0] - 1 + order[1] - 1 : 0) + offset; - } + { + if (!ibdy) + { // On i axis + return (i - 1) + (j ? order[0] - 1 + order[1] - 1 : 0) + offset; + } if (!jbdy) - { // On j axis - return (j - 1) + (i ? order[0] - 1 : 2 * (order[0] - 1) + order[1] - 1) + offset; + { // On j axis + return (j - 1) + + (i ? order[0] - 1 : 2 * (order[0] - 1) + order[1] - 1) + + offset; + } } - } offset += 2 * (order[0] - 1 + order[1] - 1); // nbdy == 0: Face DOF @@ -835,8 +841,11 @@ namespace * Modified from * https://github.com/Kitware/VTK/blob/265ca48a79a36538c95622c237da11133608bbe5/Common/DataModel/vtkLagrangeHexahedron.cxx#L734 */ - int vtk_point_index_from_ijk(int i, int j, int k, - const std::array &order) + int + vtk_point_index_from_ijk(int i, + int j, + int k, + const std::array &order) { bool ibdy = (i == 0 || i == order[0]); bool jbdy = (j == 0 || j == order[1]); @@ -845,57 +854,67 @@ namespace int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0) + (kbdy ? 1 : 0); if (nbdy == 3) // Vertex DOF - { // ijk is a corner node. Return the proper index (somewhere in [0,7]): - return (i ? (j ? 2 : 1) : (j ? 3 : 0)) + (k ? 4 : 0); - } + { // ijk is a corner node. Return the proper index (somewhere in [0,7]): + return (i ? (j ? 2 : 1) : (j ? 3 : 0)) + (k ? 4 : 0); + } int offset = 8; if (nbdy == 2) // Edge DOF - { - if (!ibdy) - { // On i axis - return (i - 1) + (j ? order[0] - 1 + order[1] - 1 : 0) + (k ? 2 * (order[0] - 1 + order[1] - 1) : 0) + offset; - } - if (!jbdy) - { // On j axis - return (j - 1) + (i ? order[0] - 1 : 2 * (order[0] - 1) + order[1] - 1) + (k ? 2 * (order[0] - 1 + order[1] - 1) : 0) + offset; + { + if (!ibdy) + { // On i axis + return (i - 1) + (j ? order[0] - 1 + order[1] - 1 : 0) + + (k ? 2 * (order[0] - 1 + order[1] - 1) : 0) + offset; + } + if (!jbdy) + { // On j axis + return (j - 1) + + (i ? order[0] - 1 : 2 * (order[0] - 1) + order[1] - 1) + + (k ? 2 * (order[0] - 1 + order[1] - 1) : 0) + offset; + } + // !kbdy, On k axis + offset += 4 * (order[0] - 1) + 4 * (order[1] - 1); + return (k - 1) + (order[2] - 1) * (i ? (j ? 3 : 1) : (j ? 2 : 0)) + + offset; } - // !kbdy, On k axis - offset += 4 * (order[0] - 1) + 4 * (order[1] - 1); - return (k - 1) + (order[2] - 1) * (i ? (j ? 3 : 1) : (j ? 2 : 0)) + offset; - } offset += 4 * (order[0] - 1 + order[1] - 1 + order[2] - 1); if (nbdy == 1) // Face DOF - { - if (ibdy) // On i-normal face - { - return (j - 1) + ((order[1] - 1) * (k - 1)) + (i ? (order[1] - 1) * (order[2] - 1) : 0) + offset; - } - offset += 2 * (order[1] - 1) * (order[2] - 1); - if (jbdy) // On j-normal face { - return (i - 1) + ((order[0] - 1) * (k - 1)) + (j ? (order[2] - 1) * (order[0] - 1) : 0) + offset; + if (ibdy) // On i-normal face + { + return (j - 1) + ((order[1] - 1) * (k - 1)) + + (i ? (order[1] - 1) * (order[2] - 1) : 0) + offset; + } + offset += 2 * (order[1] - 1) * (order[2] - 1); + if (jbdy) // On j-normal face + { + return (i - 1) + ((order[0] - 1) * (k - 1)) + + (j ? (order[2] - 1) * (order[0] - 1) : 0) + offset; + } + offset += 2 * (order[2] - 1) * (order[0] - 1); + // kbdy, On k-normal face + return (i - 1) + ((order[0] - 1) * (j - 1)) + + (k ? (order[0] - 1) * (order[1] - 1) : 0) + offset; } - offset += 2 * (order[2] - 1) * (order[0] - 1); - // kbdy, On k-normal face - return (i - 1) + ((order[0] - 1) * (j - 1)) + (k ? (order[0] - 1) * (order[1] - 1) : 0) + offset; - } // nbdy == 0: Body DOF - offset += 2 * ((order[1] - 1) * (order[2] - 1) + (order[2] - 1) * (order[0] - 1) + (order[0] - 1) * (order[1] - 1)); - return offset + (i - 1) + (order[0] - 1) * ((j - 1) + (order[1] - 1) * ((k - 1))); + offset += + 2 * ((order[1] - 1) * (order[2] - 1) + (order[2] - 1) * (order[0] - 1) + + (order[0] - 1) * (order[1] - 1)); + return offset + (i - 1) + + (order[0] - 1) * ((j - 1) + (order[1] - 1) * ((k - 1))); } - int vtk_point_index_from_ijk(int, int, int, - const std::array &) + int + vtk_point_index_from_ijk(int, int, int, const std::array &) { ExcNotImplemented(); return 0; } - int vtk_point_index_from_ijk(int, int, int, - const std::array &) + int + vtk_point_index_from_ijk(int, int, int, const std::array &) { ExcNotImplemented(); return 0; @@ -1196,8 +1215,8 @@ namespace */ template void - write_lagrange_cell(const unsigned int index, - const unsigned int start, + write_lagrange_cell(const unsigned int index, + const unsigned int start, const std::vector &connectivity); }; @@ -1569,11 +1588,11 @@ namespace template void VtkStream::write_lagrange_cell(const unsigned int, - const unsigned int start, + const unsigned int start, const std::vector &connectivity) { stream << connectivity.size(); - for(const auto& c: connectivity) + for (const auto &c : connectivity) stream << '\t' << start + c; stream << '\n'; } @@ -2550,7 +2569,8 @@ namespace DataOutBase template void - write_lagrange_cells(const std::vector> &patches, StreamType &out) + write_lagrange_cells(const std::vector> &patches, + StreamType & out) { Assert(dim <= 3 && dim > 1, ExcNotImplemented()); unsigned int first_vertex_of_patch = 0; @@ -2584,7 +2604,8 @@ namespace DataOutBase for (unsigned int i1 = 0; i1 <= n1; ++i1) { const unsigned int local_offset = i3 * d3 + i2 * d2 + i1 * d1; - const unsigned int connectivity_idx = vtk_point_index_from_ijk(i1, i2, i3, cell_order); + const unsigned int connectivity_idx = + vtk_point_index_from_ijk(i1, i2, i3, cell_order); connectivity[connectivity_idx] = local_offset; } @@ -5257,11 +5278,11 @@ namespace DataOutBase // as cell order and adjust variables accordingly, otherwise // each patch is written as a linear cell. unsigned int n_points_per_cell = GeometryInfo::vertices_per_cell; - if(flags.write_lagrange_cells) - { - n_cells = patches.size(); - n_points_per_cell = n_nodes / n_cells; - } + if (flags.write_lagrange_cells) + { + n_cells = patches.size(); + n_points_per_cell = n_nodes / n_cells; + } // in gmv format the vertex coordinates and the data have an order that is a // bit unpleasant (first all x coordinates, then all y coordinate, ...; @@ -5292,9 +5313,9 @@ namespace DataOutBase out << '\n'; ///////////////////////////////// // now for the cells - out << "CELLS " << n_cells << ' ' - << n_cells * (n_points_per_cell + 1) << '\n'; - if(flags.write_lagrange_cells) + out << "CELLS " << n_cells << ' ' << n_cells * (n_points_per_cell + 1) + << '\n'; + if (flags.write_lagrange_cells) write_lagrange_cells(patches, vtk_out); else write_cells(patches, vtk_out); @@ -5304,9 +5325,9 @@ namespace DataOutBase out << "CELL_TYPES " << n_cells << '\n'; // need to distinguish between linear and Lagrange cells - const unsigned int vtk_cell_id = flags.write_lagrange_cells - ? vtk_lagrange_cell_type[dim] - : vtk_cell_type[dim]; + const unsigned int vtk_cell_id = flags.write_lagrange_cells ? + vtk_lagrange_cell_type[dim] : + vtk_cell_type[dim]; for (unsigned int i = 0; i < n_cells; ++i) out << ' ' << vtk_cell_id; out << '\n'; -- 2.39.5