[Cleanup] Move VTK writing to a class

Now also writes out the velocity

src/one-body-sample.cc

@@ -320,6 +320,10 @@ int main(int argc, char *argv[]) {
     FrictionWriter<Dune::BitSetVector<1>> writer(frictionData, frictionalNodes);
     writer.writeInfo(alpha_initial, u_initial, v_initial);
 
+    MyVTKWriter<typename MyAssembler::VertexBasis,
+                typename MyAssembler::CellBasis>
+    vtkWriter(myAssembler.cellBasis, myAssembler.vertexBasis, "obs");
+
     // Set up TNNMG solver
     using NonlinearFactory =
         SolverFactory<dims, MyBlockProblem<ConvexProblem<
@@ -458,11 +462,10 @@ int main(int argc, char *argv[]) {
       relaxationWriter << std::endl;
 
       if (parset.get<bool>("io.writeVTK")) {
-        ScalarVector vonMisesStress;
+        ScalarVector stress;
         myAssembler.assembleVonMisesStress(youngModulus, poissonRatio, u,
-                                           vonMisesStress);
-        writeVtk(myAssembler.vertexBasis, u, alpha, myAssembler.cellBasis,
-                 vonMisesStress, (boost::format("obs%d") % run).str());
+                                           stress);
+        vtkWriter.write(u, v, alpha, stress);
       }
     }
     iterationWriter.close();

src/vtk.cc

@@ -8,28 +8,45 @@
 
 #include "vtk.hh"
 
-template <class VertexBasis, class CellBasis, class Vector, class ScalarVector>
-void writeVtk(VertexBasis const &vertexBasis, Vector const &displacement,
-              ScalarVector const &logState, CellBasis const &cellBasis,
-              ScalarVector const &stress, std::string const &filename) {
+template <class VertexBasis, class CellBasis>
+MyVTKWriter<VertexBasis, CellBasis>::MyVTKWriter(
+    CellBasis const &_cellBasis, VertexBasis const &_vertexBasis,
+    std::string _prefix)
+    : cellBasis(_cellBasis),
+      vertexBasis(_vertexBasis),
+      prefix(_prefix),
+      counter(0) {}
+
+template <class VertexBasis, class CellBasis>
+template <class Vector, class ScalarVector>
+void MyVTKWriter<VertexBasis, CellBasis>::write(Vector const &u,
+                                                Vector const &v,
+                                                ScalarVector const &alpha,
+                                                ScalarVector const &stress) {
   Dune::VTKWriter<typename VertexBasis::GridView> writer(
       vertexBasis.getGridView());
 
   auto const displacementPointer =
       std::make_shared<VTKBasisGridFunction<VertexBasis, Vector> const>(
-          vertexBasis, displacement, "displacement");
+          vertexBasis, u, "displacement");
   writer.addVertexData(displacementPointer);
 
+  auto const velocityPointer =
+      std::make_shared<VTKBasisGridFunction<VertexBasis, Vector> const>(
+          vertexBasis, v, "velocity");
+  writer.addVertexData(velocityPointer);
+
   auto const logStatePointer =
       std::make_shared<VTKBasisGridFunction<VertexBasis, ScalarVector> const>(
-          vertexBasis, logState, "logState");
+          vertexBasis, alpha, "logState");
   writer.addVertexData(logStatePointer);
 
-  auto const vonmisesPointer =
+  auto const stressPointer =
       std::make_shared<VTKBasisGridFunction<CellBasis, ScalarVector> const>(
           cellBasis, stress, "stress");
-  writer.addCellData(vonmisesPointer);
+  writer.addCellData(stressPointer);
 
+  std::string const filename = prefix + std::to_string(counter++);
   writer.write(filename.c_str());
 }
 

src/vtk.hh

 #ifndef VTK_HH
 #define VTK_HH
 
-template <class VertexBasis, class CellBasis, class Vector, class ScalarVector>
-void writeVtk(VertexBasis const &vertexBasis, Vector const &displacement,
-              ScalarVector const &logState, CellBasis const &cellBasis,
-              ScalarVector const &stress, std::string const &filename);
+template <class VertexBasis, class CellBasis> class MyVTKWriter {
+  CellBasis const &cellBasis;
+  VertexBasis const &vertexBasis;
+  std::string const prefix;
+
+  size_t counter;
+
+public:
+  MyVTKWriter(CellBasis const &cellBasis, VertexBasis const &vertexBasis,
+              std::string prefix);
+
+  template <class Vector, class ScalarVector>
+  void write(Vector const &u, Vector const &v, ScalarVector const &alpha,
+             ScalarVector const &stress);
+};
 
 #endif

src/vtk_tmpl.cc

@@ -8,10 +8,11 @@
 #include <dune/fufem/functionspacebases/p0basis.hh>
 #include <dune/fufem/functionspacebases/p1nodalbasis.hh>
 
-using P1Basis = P1NodalBasis<GridView, double>;
 using MyP0Basis = P0Basis<GridView, double>;
+using P1Basis = P1NodalBasis<GridView, double>;
+
+template class MyVTKWriter<P1Basis, MyP0Basis>;
 
-template void writeVtk<P1Basis, MyP0Basis, Vector, ScalarVector>(
-    P1Basis const &vertexBasis, Vector const &displacement,
-    ScalarVector const &logState, MyP0Basis const &cellBasis,
-    ScalarVector const &stress, std::string const &filename);
+template void MyVTKWriter<P1Basis, MyP0Basis>::write<Vector, ScalarVector>(
+    Vector const &u, Vector const &v, ScalarVector const &alpha,
+    ScalarVector const &stress);