src/hdf5/iteration-writer.hh → dune/tectonic/io/hdf5/iteration-writer.hh

@@ -3,7 +3,7 @@
 
 #include <dune/fufem/hdf5/sequenceio.hh>
 
-#include "../time-stepping/adaptivetimestepper.hh"
+#include "../../time-stepping/adaptivetimestepper.hh"
 
 class IterationWriter {
 public:

src/hdf5/patchinfo-writer.cc → dune/tectonic/io/hdf5/patchinfo-writer.cc

@@ -8,19 +8,19 @@
 #include "patchinfo-writer.hh"
 
 template <class LocalVector, class GridView>
-GridEvaluator<LocalVector, GridView>::GridEvaluator(
+GridEvaluator<LocalVector, GridView>::GridEvaluator(const CuboidGeometry& cuboidGeometry,
     ConvexPolyhedron<LocalVector> const &weakPatch, GridView const &gridView) {
-  double const bufferWidth = 0.05 * MyGeometry::lengthScale;
+  double const bufferWidth = 0.05 * cuboidGeometry.lengthScale();
   auto const xminmax = std::minmax_element(
       weakPatch.vertices.begin(), weakPatch.vertices.end(),
       [](LocalVector const &a, LocalVector const &b) { return a[0] < b[0]; });
   double const xmin = (*xminmax.first)[0] - bufferWidth;
   double const xspan = (*xminmax.second)[0] + bufferWidth - xmin;
 
-  double const zmin = -MyGeometry::depth / 2.0;
-  double const zspan = MyGeometry::depth;
+  double const zmin = -cuboidGeometry.depth() / 2.0;
+  double const zspan = cuboidGeometry.depth();
 
-  double spatialResolution = 0.01 * MyGeometry::lengthScale;
+  double spatialResolution = 0.01 * cuboidGeometry.lengthScale();
   size_t const xsteps = std::round(xspan / spatialResolution);
   size_t const zsteps = std::round(zspan / spatialResolution);
 
@@ -32,7 +32,7 @@ GridEvaluator<LocalVector, GridView>::GridEvaluator(
     zCoordinates[zi] = zmin + zi * zspan / zsteps;
 
   HierarchicApproximation<typename GridView::Grid, GridView> const
-      hApproximation(gridView.grid(), gridView, 1e-6 * MyGeometry::lengthScale);
+      hApproximation(gridView.grid(), gridView, 1e-6 * cuboidGeometry.lengthScale());
 
   LocalVector global(0);
   localInfo.resize(xsteps + 1);
@@ -68,8 +68,10 @@ template <class ProgramState, class VertexBasis, class GridView>
 PatchInfoWriter<ProgramState, VertexBasis, GridView>::PatchInfoWriter(
     HDF5::Grouplike &file, VertexBasis const &vertexBasis,
     Patch const &frictionalBoundary,
-    ConvexPolyhedron<LocalVector> const &weakPatch)
-    : group_(file, "weakPatchGrid"),
+    ConvexPolyhedron<LocalVector> const &weakPatch,
+    size_t const id)
+    : id_(id),
+      group_(file, "weakPatchGrid"+std::to_string(id_)),
       vertexBasis_(vertexBasis),
       gridEvaluator_(weakPatch, frictionalBoundary.gridView()),
       weakPatchGridVelocityWriter_(
@@ -87,7 +89,8 @@ PatchInfoWriter<ProgramState, VertexBasis, GridView>::PatchInfoWriter(
 template <class ProgramState, class VertexBasis, class GridView>
 void PatchInfoWriter<ProgramState, VertexBasis, GridView>::write(
     ProgramState const &programState) {
-  BasisGridFunction<VertexBasis, Vector> velocity(vertexBasis_, programState.v);
+    // TODO
+  BasisGridFunction<VertexBasis, Vector> velocity(vertexBasis_, programState.v[id_]);
   auto const gridVelocity = gridEvaluator_.evaluate(velocity);
   addEntry(weakPatchGridVelocityWriter_, programState.timeStep, gridVelocity);
 }

src/hdf5/patchinfo-writer.hh → dune/tectonic/io/hdf5/patchinfo-writer.hh

@@ -8,14 +8,16 @@
 #include <dune/fufem/geometry/convexpolyhedron.hh>
 #include <dune/fufem/hdf5/sequenceio.hh>
 
-#include "../one-body-problem-data/mygeometry.hh"
+#include "../../problem-data/grid/cuboidgeometry.hh"
 
 template <class LocalVector, class GridView> class GridEvaluator {
   using Element = typename GridView::Grid::template Codim<0>::Entity;
+  using CuboidGeometry = CuboidGeometry<typename LocalVector::field_type>;
 
 public:
-  GridEvaluator(ConvexPolyhedron<LocalVector> const &weakPatch,
-                GridView const &gridView);
+  GridEvaluator(const CuboidGeometry& cuboidGeometry,
+                const ConvexPolyhedron<LocalVector>& weakPatch,
+                const GridView& gridView);
 
   template <class Function>
   Dune::Matrix<typename Function::RangeType> evaluate(Function const &f) const;
@@ -36,11 +38,14 @@ class PatchInfoWriter {
 public:
   PatchInfoWriter(HDF5::Grouplike &file, VertexBasis const &vertexBasis,
                   Patch const &frictionalBoundary,
-                  ConvexPolyhedron<LocalVector> const &weakPatch);
+                  ConvexPolyhedron<LocalVector> const &weakPatch,
+                  size_t const id);
 
   void write(ProgramState const &programState);
 
 private:
+  const size_t id_;
+
   HDF5::Group group_;
 
   VertexBasis const &vertexBasis_;

dune/tectonic/io/hdf5/patchinfo-writer_tmpl.cc

+#include "../../explicitvectors.hh"
+#include "../../explicitgrid.hh"
+
+#include "../../data-structures/program_state.hh"
+
+using MyProgramState = ProgramState<Vector, ScalarVector>;
+using P1Basis = P1NodalBasis<DeformedGrid::LeafGridView, double>;
+using MyFunction = BasisGridFunction<P1Basis, Vector>;
+
+template class GridEvaluator<LocalVector, DeformedGrid::LeafGridView>;
+
+template Dune::Matrix<typename MyFunction::RangeType>
+GridEvaluator<LocalVector, DeformedGrid::LeafGridView>::evaluate(MyFunction const &f) const;
+
+template class PatchInfoWriter<MyProgramState, P1Basis, DeformedGrid::LeafGridView>;

src/hdf5/restart-io.cc → dune/tectonic/io/hdf5/restart-io.cc

@@ -5,26 +5,66 @@
 #include "restart-io.hh"
 
 template <class ProgramState>
-RestartIO<ProgramState>::RestartIO(HDF5::Grouplike &group, size_t vertexCount)
-    : displacementWriter_(group, "displacement", vertexCount,
-                          Vector::block_type::dimension),
-      velocityWriter_(group, "velocity", vertexCount,
-                      Vector::block_type::dimension),
-      accelerationWriter_(group, "acceleration", vertexCount,
-                          Vector::block_type::dimension),
-      stateWriter_(group, "state", vertexCount),
-      weightedNormalStressWriter_(group, "weightedNormalStress", vertexCount),
+RestartBodyIO<ProgramState>::RestartBodyIO(HDF5::Grouplike &group, const size_t vertexCount, const size_t id)
+    : id_(id),
+      displacementWriter_(group, "displacement"+std::to_string(id_), vertexCount,
+                          ProgramState::Vector::block_type::dimension),
+      velocityWriter_(group, "velocity"+std::to_string(id_), vertexCount,
+                      ProgramState::Vector::block_type::dimension),
+      accelerationWriter_(group, "acceleration"+std::to_string(id_), vertexCount,
+                          ProgramState::Vector::block_type::dimension),
+      stateWriter_(group, "state"+std::to_string(id_), vertexCount),
+      weightedNormalStressWriter_(group, "weightedNormalStress"+std::to_string(id_), vertexCount) {}
+
+template <class ProgramState>
+void RestartBodyIO<ProgramState>::write(const ProgramState& programState) {
+
+  addEntry(displacementWriter_, programState.timeStep, programState.u[id_]);
+  addEntry(velocityWriter_, programState.timeStep, programState.v[id_]);
+  addEntry(accelerationWriter_, programState.timeStep, programState.a[id_]);
+  addEntry(stateWriter_, programState.timeStep, programState.alpha[id_]);
+  addEntry(weightedNormalStressWriter_, programState.timeStep,
+           programState.weightedNormalStress[id_]);
+}
+
+template <class ProgramState>
+void RestartBodyIO<ProgramState>::read(size_t timeStep,
+                                   ProgramState& programState) {
+
+  readEntry(displacementWriter_, timeStep, programState.u[id_]);
+  readEntry(velocityWriter_, timeStep, programState.v[id_]);
+  readEntry(accelerationWriter_, timeStep, programState.a[id_]);
+  readEntry(stateWriter_, timeStep, programState.alpha[id_]);
+  readEntry(weightedNormalStressWriter_, timeStep,
+            programState.weightedNormalStress[id_]);
+}
+
+template <class ProgramState>
+RestartIO<ProgramState>::RestartIO(HDF5::Grouplike& group, const std::vector<size_t>& vertexCounts)
+    : bodyIO_(vertexCounts.size()),
       relativeTimeWriter_(group, "relativeTime"),
-      relativeTimeIncrementWriter_(group, "relativeTimeIncrement") {}
+      relativeTimeIncrementWriter_(group, "relativeTimeIncrement") {
+
+  for (size_t i=0; i<bodyIO_.size(); i++) {
+    bodyIO_[i] = new RestartBodyIO<ProgramState>(group, vertexCounts[i], i);
+  }
+}
+
+template <class ProgramState>
+RestartIO<ProgramState>::~RestartIO() {
+  for (size_t i=0; i<bodyIO_.size(); i++) {
+    delete bodyIO_[i];
+  }
+}
 
 template <class ProgramState>
 void RestartIO<ProgramState>::write(ProgramState const &programState) {
-  addEntry(displacementWriter_, programState.timeStep, programState.u);
-  addEntry(velocityWriter_, programState.timeStep, programState.v);
-  addEntry(accelerationWriter_, programState.timeStep, programState.a);
-  addEntry(stateWriter_, programState.timeStep, programState.alpha);
-  addEntry(weightedNormalStressWriter_, programState.timeStep,
-           programState.weightedNormalStress);
+  assert(programState.size() == bodyIO_.size());
+
+  for (size_t i=0; i<bodyIO_.size(); i++) {
+    bodyIO_[i]->write(programState);
+  }
+
   addEntry(relativeTimeWriter_, programState.timeStep,
            programState.relativeTime);
   addEntry(relativeTimeIncrementWriter_, programState.timeStep,
@@ -33,14 +73,15 @@ void RestartIO<ProgramState>::write(ProgramState const &programState) {
 
 template <class ProgramState>
 void RestartIO<ProgramState>::read(size_t timeStep,
-                                   ProgramState &programState) {
+                                   ProgramState& programState) {
+  assert(programState.size() == bodyIO_.size());
+
   programState.timeStep = timeStep;
-  readEntry(displacementWriter_, timeStep, programState.u);
-  readEntry(velocityWriter_, timeStep, programState.v);
-  readEntry(accelerationWriter_, timeStep, programState.a);
-  readEntry(stateWriter_, timeStep, programState.alpha);
-  readEntry(weightedNormalStressWriter_, timeStep,
-            programState.weightedNormalStress);
+
+  for (size_t i=0; i<bodyIO_.size(); i++) {
+    bodyIO_[i]->read(timeStep, programState);
+  }
+
   readEntry(relativeTimeWriter_, timeStep, programState.relativeTime);
   readEntry(relativeTimeIncrementWriter_, timeStep, programState.relativeTau);
 }

src/hdf5/restart-io.hh → dune/tectonic/io/hdf5/restart-io.hh

@@ -4,23 +4,37 @@
 #include <dune/fufem/hdf5/file.hh>
 #include <dune/fufem/hdf5/sequenceio.hh>
 
-template <class ProgramState> class RestartIO {
-  using ScalarVector = typename ProgramState::ScalarVector;
-  using Vector = typename ProgramState::Vector;
-
+template <class ProgramState> class RestartBodyIO {
 public:
-  RestartIO(HDF5::Grouplike &group, size_t vertexCount);
+  RestartBodyIO(HDF5::Grouplike& group, const size_t vertexCount, const size_t id);
 
-  void write(ProgramState const &programState);
+  void write(const ProgramState& programState);
 
-  void read(size_t timeStep, ProgramState &programState);
+  void read(size_t timeStep, ProgramState& programState);
 
 private:
+  const size_t id_;
+
   HDF5::SequenceIO<2> displacementWriter_;
   HDF5::SequenceIO<2> velocityWriter_;
   HDF5::SequenceIO<2> accelerationWriter_;
   HDF5::SequenceIO<1> stateWriter_;
   HDF5::SequenceIO<1> weightedNormalStressWriter_;
+};
+
+template <class ProgramState> class RestartIO {
+
+public:
+  RestartIO(HDF5::Grouplike& group, const std::vector<size_t>& vertexCounts);
+  ~RestartIO();
+
+  void write(const ProgramState& programState);
+
+  void read(size_t timeStep, ProgramState& programState);
+
+private:
+  std::vector<RestartBodyIO<ProgramState>* > bodyIO_;
+
   HDF5::SequenceIO<0> relativeTimeWriter_;
   HDF5::SequenceIO<0> relativeTimeIncrementWriter_;
 };

src/hdf5/restart-io_tmpl.cc → dune/tectonic/io/hdf5/restart-io_tmpl.cc

-#include "../explicitvectors.hh"
+#include "../../explicitvectors.hh"
 
-#include "../program_state.hh"
+#include "../../data-structures/program_state.hh"
 
 using MyProgramState = ProgramState<Vector, ScalarVector>;
 
+template class RestartBodyIO<MyProgramState>;
 template class RestartIO<MyProgramState>;

src/hdf5/restrict.hh → dune/tectonic/io/hdf5/restrict.hh

@@ -5,7 +5,7 @@
 
 #include <dune/common/bitsetvector.hh>
 
-#include "../tobool.hh"
+#include "../../utils/tobool.hh"
 
 template <class Vector, class Patch>
 Vector restrictToSurface(Vector const &v1, Patch const &patch) {

src/hdf5/surface-writer.cc → dune/tectonic/io/hdf5/surface-writer.cc

@@ -7,8 +7,9 @@
 
 template <class ProgramState, class GridView>
 SurfaceWriter<ProgramState, GridView>::SurfaceWriter(
-    HDF5::Grouplike &file, Vector const &vertexCoordinates, Patch const &surface)
-    : group_(file, "surface"),
+    HDF5::Grouplike &file, Vector const &vertexCoordinates, Patch const &surface, size_t const id)
+    : id_(id),
+      group_(file, "surface"+std::to_string(id_)),
       surface_(surface),
       surfaceDisplacementWriter_(group_, "displacement", surface.numVertices(),
                                  Vector::block_type::dimension),
@@ -24,11 +25,12 @@ SurfaceWriter<ProgramState, GridView>::SurfaceWriter(
 template <class ProgramState, class GridView>
 void SurfaceWriter<ProgramState, GridView>::write(
     ProgramState const &programState) {
-  auto const surfaceDisplacements = restrictToSurface(programState.u, surface_);
+
+  auto const surfaceDisplacements = restrictToSurface(programState.u[id_], surface_);
   addEntry(surfaceDisplacementWriter_, programState.timeStep,
            surfaceDisplacements);
 
-  auto const surfaceVelocities = restrictToSurface(programState.v, surface_);
+  auto const surfaceVelocities = restrictToSurface(programState.v[id_], surface_);
   addEntry(surfaceVelocityWriter_, programState.timeStep, surfaceVelocities);
 }
 

src/hdf5/surface-writer.hh → dune/tectonic/io/hdf5/surface-writer.hh

@@ -12,11 +12,13 @@ template <class ProgramState, class GridView> class SurfaceWriter {
 
 public:
   SurfaceWriter(HDF5::Grouplike &file, Vector const &vertexCoordinates,
-                Patch const &surface);
+                Patch const &surface, size_t const id);
 
   void write(ProgramState const &programState);
 
 private:
+  size_t const id_;
+
   HDF5::Group group_;
 
   Patch const &surface_;

dune/tectonic/io/hdf5/surface-writer_tmpl.cc

+#include "../../explicitvectors.hh"
+#include "../../explicitgrid.hh"
+
+#include "../../data-structures/program_state.hh"
+
+using MyProgramState = ProgramState<Vector, ScalarVector>;
+
+template class SurfaceWriter<MyProgramState, DeformedGrid::LeafGridView>;

src/hdf5/time-writer_tmpl.cc → dune/tectonic/io/hdf5/time-writer_tmpl.cc

-#include "../explicitvectors.hh"
+#include "../../explicitvectors.hh"
 
-#include "../program_state.hh"
+#include "../../data-structures/program_state.hh"
 
 using MyProgramState = ProgramState<Vector, ScalarVector>;
 

src/uniform-grid-writer.cc → dune/tectonic/io/uniform-grid-writer.cc

@@ -11,10 +11,13 @@
 // #include <dune/common/parametertreeparser.hh>
 
 #include "assemblers.hh"
-#include "diameter.hh"
 #include "gridselector.hh"
+
+#include "io/vtk.hh"
+
 #include "one-body-problem-data/mygrid.hh"
-#include "vtk.hh"
+
+#include "utils/diameter.hh"
 
 size_t const dims = MY_DIM;
 size_t const refinements = 5; // FIXME?

src/vtk.cc → dune/tectonic/io/vtk.cc

@@ -8,51 +8,84 @@
 
 #include "vtk.hh"
 
+#include "../utils/debugutils.hh"
+
 template <class VertexBasis, class CellBasis>
-MyVTKWriter<VertexBasis, CellBasis>::MyVTKWriter(
-    CellBasis const &_cellBasis, VertexBasis const &_vertexBasis,
-    std::string _prefix)
-    : cellBasis(_cellBasis), vertexBasis(_vertexBasis), prefix(_prefix) {}
+BodyVTKWriter<VertexBasis, CellBasis>::BodyVTKWriter(
+    CellBasis const & cellBasis, VertexBasis const & vertexBasis,
+    std::string prefix)
+    : cellBasis_(cellBasis), vertexBasis_(vertexBasis), prefix_(prefix) {}
 
 template <class VertexBasis, class CellBasis>
 template <class Vector, class ScalarVector>
-void MyVTKWriter<VertexBasis, CellBasis>::write(
+void BodyVTKWriter<VertexBasis, CellBasis>::write(
     size_t record, Vector const &u, Vector const &v, ScalarVector const &alpha,
     ScalarVector const &stress) const {
   Dune::VTKWriter<typename VertexBasis::GridView> writer(
-      vertexBasis.getGridView());
+      vertexBasis_.getGridView());
 
   auto const displacementPointer =
       std::make_shared<VTKBasisGridFunction<VertexBasis, Vector> const>(
-          vertexBasis, u, "displacement");
+          vertexBasis_, u, "displacement");
   writer.addVertexData(displacementPointer);
 
   auto const velocityPointer =
       std::make_shared<VTKBasisGridFunction<VertexBasis, Vector> const>(
-          vertexBasis, v, "velocity");
+          vertexBasis_, v, "velocity");
   writer.addVertexData(velocityPointer);
 
   auto const AlphaPointer =
       std::make_shared<VTKBasisGridFunction<VertexBasis, ScalarVector> const>(
-          vertexBasis, alpha, "Alpha");
+          vertexBasis_, alpha, "Alpha");
   writer.addVertexData(AlphaPointer);
 
   auto const stressPointer =
       std::make_shared<VTKBasisGridFunction<CellBasis, ScalarVector> const>(
-          cellBasis, stress, "stress");
+          cellBasis_, stress, "stress");
   writer.addCellData(stressPointer);
 
-  std::string const filename = prefix + std::to_string(record);
+  std::string const filename = prefix_ + "_" + std::to_string(record);
   writer.write(filename.c_str(), Dune::VTK::appendedraw);
 }
 
 template <class VertexBasis, class CellBasis>
-void MyVTKWriter<VertexBasis, CellBasis>::writeGrid() const {
+void BodyVTKWriter<VertexBasis, CellBasis>::writeGrid() const {
   Dune::VTKWriter<typename VertexBasis::GridView> writer(
-      vertexBasis.getGridView());
+      vertexBasis_.getGridView());
 
-  std::string const filename = prefix + "_grid";
+  std::string const filename = prefix_ + "_grid";
   writer.write(filename.c_str(), Dune::VTK::appendedraw);
 }
 
+template <class VertexBasis, class CellBasis>
+MyVTKWriter<VertexBasis, CellBasis>::MyVTKWriter(const std::vector<const CellBasis* >& cellBases, const std::vector<const VertexBasis* >& vertexBases,
+  std::string prefix): bodyVTKWriters_(vertexBases.size()), prefix_(prefix) {
+
+  for (size_t i=0; i<bodyVTKWriters_.size(); i++) {
+    bodyVTKWriters_[i] = new BodyVTKWriter<VertexBasis, CellBasis>(*cellBases[i], *vertexBases[i], prefix_+std::to_string(i));
+  }
+}
+
+template <class VertexBasis, class CellBasis>
+MyVTKWriter<VertexBasis, CellBasis>::~MyVTKWriter() {
+  for (size_t i=0; i<bodyVTKWriters_.size(); i++) {
+    delete bodyVTKWriters_[i];
+  }
+}
+
+template <class VertexBasis, class CellBasis>
+template <class Vector, class ScalarVector>
+void MyVTKWriter<VertexBasis, CellBasis>::write(size_t record, const std::vector<Vector>& u, const std::vector<Vector>& v,
+                        const std::vector<ScalarVector>& alpha, const std::vector<ScalarVector>& stress) const {
+
+  for (size_t i=0; i<bodyVTKWriters_.size(); i++) {
+    bodyVTKWriters_[i]->write(record, u[i], v[i], alpha[i], stress[i]);
+  }
+}
+template <class VertexBasis, class CellBasis>
+void MyVTKWriter<VertexBasis, CellBasis>::writeGrids() const {
+  for (size_t i=0; i<bodyVTKWriters_.size(); i++) {
+    bodyVTKWriters_[i]->writeGrid();
+  }
+}
 #include "vtk_tmpl.cc"

src/vtk.hh → dune/tectonic/io/vtk.hh

@@ -3,13 +3,14 @@
 
 #include <string>
 
-template <class VertexBasis, class CellBasis> class MyVTKWriter {
-  CellBasis const &cellBasis;
-  VertexBasis const &vertexBasis;
-  std::string const prefix;
+template <class VertexBasis, class CellBasis> class BodyVTKWriter {
+private:
+  CellBasis const &cellBasis_;
+  VertexBasis const &vertexBasis_;
+  std::string const prefix_;
 
 public:
-  MyVTKWriter(CellBasis const &cellBasis, VertexBasis const &vertexBasis,
+  BodyVTKWriter(CellBasis const &cellBasis, VertexBasis const &vertexBasis,
               std::string prefix);
 
   template <class Vector, class ScalarVector>
@@ -18,4 +19,22 @@ template <class VertexBasis, class CellBasis> class MyVTKWriter {
 
   void writeGrid() const;
 };
+
+template <class VertexBasis, class CellBasis> class MyVTKWriter {
+private:
+  std::vector<BodyVTKWriter<VertexBasis, CellBasis>* > bodyVTKWriters_;
+  std::string const prefix_;
+
+public:
+  MyVTKWriter(const std::vector<const CellBasis* >& cellBases, const std::vector<const VertexBasis* >& vertexBases,
+              std::string prefix);
+
+  ~MyVTKWriter();
+
+  template <class Vector, class ScalarVector>
+  void write(size_t record, const std::vector<Vector>& u, const std::vector<Vector>& v,
+             const std::vector<ScalarVector>& alpha, const std::vector<ScalarVector>& stress) const;
+
+  void writeGrids() const;
+};
 #endif

src/vtk_tmpl.cc → dune/tectonic/io/vtk_tmpl.cc

@@ -2,17 +2,17 @@
 #error MY_DIM unset
 #endif
 
-#include "explicitgrid.hh"
-#include "explicitvectors.hh"
+#include "../explicitgrid.hh"
+#include "../explicitvectors.hh"
 
 #include <dune/fufem/functionspacebases/p0basis.hh>
 #include <dune/fufem/functionspacebases/p1nodalbasis.hh>
 
-using MyP0Basis = P0Basis<GridView, double>;
-using P1Basis = P1NodalBasis<GridView, double>;
+using MyP0Basis = P0Basis<DeformedGrid::LeafGridView, double>;
+using P1Basis = P1NodalBasis<DeformedGrid::LeafGridView, double>;
 
 template class MyVTKWriter<P1Basis, MyP0Basis>;
 
 template void MyVTKWriter<P1Basis, MyP0Basis>::write<Vector, ScalarVector>(
-    size_t record, Vector const &u, Vector const &v, ScalarVector const &alpha,
-    ScalarVector const &stress) const;
+    size_t record, const std::vector<Vector>& u, const std::vector<Vector>& v,
+        const std::vector<ScalarVector>& alpha, const std::vector<ScalarVector>& stress) const;

dune/tectonic/minimisation.hh

-#ifndef DUNE_TECTONIC_MINIMISATION_HH
-#define DUNE_TECTONIC_MINIMISATION_HH
-
-#include <dune/common/fmatrix.hh>
-#include <dune/common/fvector.hh>
-#include <dune/common/stdstreams.hh>
-
-#include <dune/fufem/arithmetic.hh>
-#include <dune/fufem/interval.hh>
-#include <dune/tnnmg/problem-classes/bisection.hh>
-
-#include <dune/tectonic/mydirectionalconvexfunction.hh>
-
-// Warning: this exploits the property v*x = 0
-template <class Functional>
-double lineSearch(Functional const &J,
-                  typename Functional::LocalVector const &x,
-                  typename Functional::LocalVector const &v,
-                  Bisection const &bisection) {
-  MyDirectionalConvexFunction<typename Functional::Nonlinearity> const JRest(
-      J.alpha * v.two_norm2(), J.b * v, J.phi, x, v);
-  int count;
-  return bisection.minimize(JRest, 0.0, 0.0, count);
-}
-
-/** Minimise a quadratic problem, for which both the quadratic and the
-    nonlinear term have gradients which point in the direction of
-    their arguments */
-template <class Functional>
-void minimise(Functional const &J, typename Functional::LocalVector &x,
-              Bisection const &bisection) {
-  auto v = J.b;
-  double const vnorm = v.two_norm();
-  if (vnorm <= 0.0)
-    return;
-  v /= vnorm;
-
-  double const alpha = lineSearch(J, x, v, bisection);
-  Arithmetic::addProduct(x, alpha, v);
-}
-#endif