]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Inline most of the remaining functions in VtuStream into write_vtu_main().
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 21 Nov 2022 18:46:01 +0000 (11:46 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Sat, 26 Nov 2022 03:57:50 +0000 (20:57 -0700)
source/base/data_out_base.cc

index c333d04335d604b2095a02035e846ceb7afc6c8f..09ee6e76fc3e6d61a0e72da52e4191bf41a5caeb 100644 (file)
@@ -1551,29 +1551,6 @@ namespace
   public:
     VtuStream(std::ostream &stream, const DataOutBase::VtkFlags &flags);
 
-    /**
-     * The order of vertices for these cells in different dimensions is
-     * <ol>
-     * <li> [0,1]
-     * <li> []
-     * <li> []
-     * </ol>
-     */
-    template <int dim>
-    void
-    write_cell(const unsigned int                   index,
-               const unsigned int                   start,
-               const std::array<unsigned int, dim> &offsets);
-
-    /**
-     * Print vertices [start, start+n_points[
-     */
-    void
-    write_cell_single(const unsigned int   index,
-                      const unsigned int   start,
-                      const unsigned int   n_points,
-                      const ReferenceCell &reference_cell);
-
     /**
      * Write a high-order cell type, i.e., a Lagrange cell
      * in the VTK terminology.
@@ -2062,137 +2039,6 @@ namespace
 
 
 
-  template <int dim>
-  void
-  VtuStream::write_cell(const unsigned int,
-                        const unsigned int                   start,
-                        const std::array<unsigned int, dim> &offsets)
-  {
-#ifdef DEAL_II_WITH_ZLIB
-    if (flags.compression_level != DataOutBase::CompressionLevel::plain_text)
-      {
-        switch (dim)
-          {
-            case 0:
-              {
-                cells.push_back(start);
-                break;
-              }
-
-            case 1:
-              {
-                const unsigned int d1 = offsets[0];
-                cells.push_back(start);
-                cells.push_back(start + d1);
-                break;
-              }
-
-            case 2:
-              {
-                const unsigned int d1 = offsets[0];
-                const unsigned int d2 = offsets[1];
-                cells.push_back(start);
-                cells.push_back(start + d1);
-                cells.push_back(start + d2 + d1);
-                cells.push_back(start + d2);
-                break;
-              }
-
-            case 3:
-              {
-                const unsigned int d1 = offsets[0];
-                const unsigned int d2 = offsets[1];
-                const unsigned int d3 = offsets[2];
-                cells.push_back(start);
-                cells.push_back(start + d1);
-                cells.push_back(start + d2 + d1);
-                cells.push_back(start + d2);
-                cells.push_back(start + d3);
-                cells.push_back(start + d3 + d1);
-                cells.push_back(start + d3 + d2 + d1);
-                cells.push_back(start + d3 + d2);
-                break;
-              }
-
-            default:
-              Assert(false, ExcNotImplemented());
-          }
-      }
-    else
-#endif
-      {
-        switch (dim)
-          {
-            case 0:
-              {
-                stream << start;
-                break;
-              }
-
-            case 1:
-              {
-                const unsigned int d1 = offsets[0];
-                stream << start << '\t' << start + d1;
-                break;
-              }
-
-            case 2:
-              {
-                const unsigned int d1 = offsets[0];
-                const unsigned int d2 = offsets[1];
-                stream << start << '\t' << start + d1 << '\t' << start + d2 + d1
-                       << '\t' << start + d2;
-                break;
-              }
-
-            case 3:
-              {
-                const unsigned int d1 = offsets[0];
-                const unsigned int d2 = offsets[1];
-                const unsigned int d3 = offsets[2];
-                stream << start << '\t' << start + d1 << '\t' << start + d2 + d1
-                       << '\t' << start + d2 << '\t' << start + d3 << '\t'
-                       << start + d3 + d1 << '\t' << start + d3 + d2 + d1
-                       << '\t' << start + d3 + d2;
-                break;
-              }
-
-            default:
-              Assert(false, ExcNotImplemented());
-          }
-        stream << '\n';
-      }
-  }
-
-  void
-  VtuStream::write_cell_single(const unsigned int   index,
-                               const unsigned int   start,
-                               const unsigned int   n_points,
-                               const ReferenceCell &reference_cell)
-  {
-    (void)index;
-
-    static const std::array<unsigned int, 5> table = {{0, 1, 3, 2, 4}};
-
-#ifdef DEAL_II_WITH_ZLIB
-    if (flags.compression_level != DataOutBase::CompressionLevel::plain_text)
-      {
-        for (unsigned int i = 0; i < n_points; ++i)
-          cells.push_back(
-            start + (reference_cell == ReferenceCells::Pyramid ? table[i] : i));
-      }
-    else
-#endif
-      {
-        for (unsigned int i = 0; i < n_points; ++i)
-          stream << '\t'
-                 << start + (reference_cell == ReferenceCells::Pyramid ?
-                               table[i] :
-                               i);
-        stream << '\n';
-      }
-  }
-
   template <int dim>
   void
   VtuStream::write_high_order_cell(const unsigned int,
@@ -2224,8 +2070,7 @@ namespace
         // then release the data
         *this << vtu_stringize_array(cells,
                                      flags.compression_level,
-                                     stream.precision())
-              << '\n';
+                                     stream.precision());
         cells.clear();
       }
 #endif
@@ -6279,28 +6124,61 @@ namespace DataOutBase
     out << "  </Points>\n\n";
     //-------------------------------
     // now for the cells
-    VtuStream vtu_out(out, flags);
-
     out << "  <Cells>\n";
     out << "    <DataArray type=\"Int32\" Name=\"connectivity\" format=\""
         << ascii_or_binary << "\">\n";
     if (flags.write_higher_order_cells)
-      write_high_order_cells(patches, vtu_out, /* legacy_format = */ false);
+      {
+        std::ostringstream o;
+        {
+          VtuStream vtu_out(o, flags);
+
+          write_high_order_cells(patches, vtu_out, /* legacy_format = */ false);
+          vtu_out.flush_cells();
+        }
+        out << o.str() << '\n';
+      }
     else
       {
         Assert(dim <= 3, ExcNotImplemented());
-        unsigned int count                 = 0;
-        unsigned int first_vertex_of_patch = 0;
+
+        std::vector<int32_t> cells;
+        unsigned int         first_vertex_of_patch = 0;
+
         for (const auto &patch : patches)
           {
             // special treatment of simplices since they are not subdivided
             if (patch.reference_cell != ReferenceCells::get_hypercube<dim>())
               {
-                vtu_out.write_cell_single(count++,
-                                          first_vertex_of_patch,
-                                          patch.data.n_cols(),
-                                          patch.reference_cell);
-                first_vertex_of_patch += patch.data.n_cols();
+                const unsigned int n_points = patch.data.n_cols();
+                static const std::array<unsigned int, 5>
+                  pyramid_index_translation_table = {{0, 1, 3, 2, 4}};
+
+#ifdef DEAL_II_WITH_ZLIB
+                if (flags.compression_level !=
+                    DataOutBase::CompressionLevel::plain_text)
+                  {
+                    for (unsigned int i = 0; i < n_points; ++i)
+                      cells.push_back(
+                        first_vertex_of_patch +
+                        (patch.reference_cell == ReferenceCells::Pyramid ?
+                           pyramid_index_translation_table[i] :
+                           i));
+                  }
+                else
+#endif
+                  {
+                    for (unsigned int i = 0; i < n_points; ++i)
+                      out << '\t'
+                          << first_vertex_of_patch +
+                               (patch.reference_cell ==
+                                    ReferenceCells::Pyramid ?
+                                  pyramid_index_translation_table[i] :
+                                  i);
+                    out << '\n';
+                  }
+
+                first_vertex_of_patch += n_points;
               }
             else
               {
@@ -6311,26 +6189,54 @@ namespace DataOutBase
                   {
                     case 0:
                       {
+                        auto write_cell =
+                          [&flags, &out, &cells](const unsigned int start) {
+#ifdef DEAL_II_WITH_ZLIB
+                            if (flags.compression_level !=
+                                DataOutBase::CompressionLevel::plain_text)
+                              {
+                                cells.push_back(start);
+                              }
+                            else
+#endif
+                              {
+                                out << start;
+                                out << '\n';
+                              }
+                          };
+
                         const unsigned int starting_offset =
                           first_vertex_of_patch;
-                        // First write line in x direction
-                        vtu_out.template write_cell<0>(count++,
-                                                       starting_offset,
-                                                       {});
+                        write_cell(starting_offset);
                         break;
                       }
 
                     case 1:
                       {
                         const unsigned int d1 = 1;
+
+                        auto write_cell =
+                          [&flags, &out, &cells, d1](const unsigned int start) {
+#ifdef DEAL_II_WITH_ZLIB
+                            if (flags.compression_level !=
+                                DataOutBase::CompressionLevel::plain_text)
+                              {
+                                cells.push_back(start);
+                                cells.push_back(start + d1);
+                              }
+                            else
+#endif
+                              {
+                                out << start << '\t' << start + d1;
+                                out << '\n';
+                              }
+                          };
+
                         for (unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
                           {
                             const unsigned int starting_offset =
                               first_vertex_of_patch + i1 * d1;
-                            // First write line in x direction
-                            vtu_out.template write_cell<1>(count++,
-                                                           starting_offset,
-                                                           {{d1}});
+                            write_cell(starting_offset);
                           }
                         break;
                       }
@@ -6339,15 +6245,33 @@ namespace DataOutBase
                       {
                         const unsigned int d1 = 1;
                         const unsigned int d2 = n;
+
+                        auto write_cell = [&flags, &out, &cells, d1, d2](
+                                            const unsigned int start) {
+#ifdef DEAL_II_WITH_ZLIB
+                          if (flags.compression_level !=
+                              DataOutBase::CompressionLevel::plain_text)
+                            {
+                              cells.push_back(start);
+                              cells.push_back(start + d1);
+                              cells.push_back(start + d2 + d1);
+                              cells.push_back(start + d2);
+                            }
+                          else
+#endif
+                            {
+                              out << start << '\t' << start + d1 << '\t'
+                                  << start + d2 + d1 << '\t' << start + d2;
+                              out << '\n';
+                            }
+                        };
+
                         for (unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
                           for (unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
                             {
                               const unsigned int starting_offset =
                                 first_vertex_of_patch + i2 * d2 + i1 * d1;
-                              // First write line in x direction
-                              vtu_out.template write_cell<2>(count++,
-                                                             starting_offset,
-                                                             {{d1, d2}});
+                              write_cell(starting_offset);
                             }
                         break;
                       }
@@ -6357,6 +6281,35 @@ namespace DataOutBase
                         const unsigned int d1 = 1;
                         const unsigned int d2 = n;
                         const unsigned int d3 = n * n;
+
+                        auto write_cell = [&flags, &out, &cells, d1, d2, d3](
+                                            const unsigned int start) {
+#ifdef DEAL_II_WITH_ZLIB
+                          if (flags.compression_level !=
+                              DataOutBase::CompressionLevel::plain_text)
+                            {
+                              cells.push_back(start);
+                              cells.push_back(start + d1);
+                              cells.push_back(start + d2 + d1);
+                              cells.push_back(start + d2);
+                              cells.push_back(start + d3);
+                              cells.push_back(start + d3 + d1);
+                              cells.push_back(start + d3 + d2 + d1);
+                              cells.push_back(start + d3 + d2);
+                            }
+                          else
+#endif
+                            {
+                              out << start << '\t' << start + d1 << '\t'
+                                  << start + d2 + d1 << '\t' << start + d2
+                                  << '\t' << start + d3 << '\t'
+                                  << start + d3 + d1 << '\t'
+                                  << start + d3 + d2 + d1 << '\t'
+                                  << start + d3 + d2;
+                              out << '\n';
+                            }
+                        };
+
                         for (unsigned int i3 = 0; i3 < n_subdivisions; ++i3)
                           for (unsigned int i2 = 0; i2 < n_subdivisions; ++i2)
                             for (unsigned int i1 = 0; i1 < n_subdivisions; ++i1)
@@ -6364,10 +6317,7 @@ namespace DataOutBase
                                 const unsigned int starting_offset =
                                   first_vertex_of_patch + i3 * d3 + i2 * d2 +
                                   i1 * d1;
-                                // First write line in x direction
-                                vtu_out.template write_cell<3>(count++,
-                                                               starting_offset,
-                                                               {{d1, d2, d3}});
+                                write_cell(starting_offset);
                               }
                         break;
                       }
@@ -6382,7 +6332,17 @@ namespace DataOutBase
               }
           }
 
-        vtu_out.flush_cells();
+          // Flush the 'cells' object we created herein.
+#ifdef DEAL_II_WITH_ZLIB
+        if (flags.compression_level !=
+            DataOutBase::CompressionLevel::plain_text)
+          {
+            out << vtu_stringize_array(cells,
+                                       flags.compression_level,
+                                       out.precision())
+                << '\n';
+          }
+#endif
       }
     out << "    </DataArray>\n";
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.