Feature/stress plot

src/film-on-substrate.cc

@@ -230,6 +230,12 @@ int main (int argc, char *argv[]) try
     power<dim>(
         lagrange<displacementOrder>()
   ));
+  auto orientationPowerBasis = makeBasis(
+    gridView,
+    power<dim>(
+        lagrange<rotationOrder>()
+  ));
+
   typedef Dune::Functions::LagrangeBasis<GridView,displacementOrder> DeformationFEBasis;
   typedef Dune::Functions::LagrangeBasis<GridView,rotationOrder> OrientationFEBasis;
   DeformationFEBasis deformationFEBasis(gridView);
@@ -628,7 +634,6 @@ int main (int argc, char *argv[]) try
     /////////////////////////////////
 
     // Compute the displacement
-    BlockVector<FieldVector<double,dim> > displacement(compositeBasis.size({0}));
     for (int i = 0; i < compositeBasis.size({0}); i++) {
        for (int j = 0; j  < dim; j++) {
         displacement[i][j] = x[_0][i][j];
@@ -642,6 +647,51 @@ int main (int argc, char *argv[]) try
     vtkWriter.addVertexData(displacementFunction, VTK::FieldInfo("displacement", VTK::FieldInfo::Type::scalar, dim));
     vtkWriter.write(resultPath + "finite-strain_homotopy_" + parameterSet.get<std::string>("energy") + "_" + std::to_string(neumannValues[0]) + "_" + std::to_string(i+1));
   }
+  std::string ending = grid->leafGridView().comm().size() > 1 ? std::to_string(mpiHelper.rank()) : "";
+  std::ofstream file;
+  file.open("deformation" + ending);
+  for (int i = 0; i < identity.size(); i++){
+    file << identity[i] << ":" << displacement[i] << "\n";
+  }
+
+  file.close();
+  
+  BlockVector<FieldVector<double,dim> > identityRotation(orientationFEBasis.size());
+  Dune::Functions::interpolate(orientationPowerBasis, identityRotation, [](FieldVector<double,dim> x){ return x; });
+
+  file.open("rotation" + ending);
+  for (int i = 0; i < identityRotation.size(); i++){
+    file << identityRotation[i] << ":" << x[_1][i] << "\n";
+  }
+
+  file.close();
+
+  MPI_Barrier(grid->leafGridView().comm());
+
+  if (grid->leafGridView().comm().size() > 1 && mpiHelper.rank() == 0) {
+    file.open("deformation");
+    for (int i = 0; i < grid->leafGridView().comm().size(); i++) {
+      std::ifstream deformationInput("deformation" + std::to_string(i));
+      if (deformationInput.is_open()) {
+        file << deformationInput.rdbuf();
+      }
+      deformationInput.close();
+      std::remove(("deformation" + std::to_string(i)).c_str());
+    }
+    file.close();
+
+    file.open("rotation");
+    for (int i = 0; i < grid->leafGridView().comm().size(); i++) {
+      std::ifstream rotationInput("rotation" + std::to_string(i));
+      if (rotationInput.is_open()) {
+        file << rotationInput.rdbuf();
+      }
+      rotationInput.close();
+      std::remove(("rotation" + std::to_string(i)).c_str());
+    }
+    file.close();
+  }
+
 } catch (Exception& e) {
     std::cout << e.what() << std::endl;
 }