Handle choice Ruina/Laursen in the parset file

src/one-body-sample.cc

@@ -216,34 +216,44 @@ int main(int argc, char *argv[]) {
       stiffnessMatrix.umv(u3, b3);
 
       // {{{ Assemble terms for the nonlinearity
-      std::vector<double> a;
-      a.resize(grid.size(grid.maxLevel(), dim));
-      std::fill(a.begin(), a.end(), parset.get<double>("boundary.a"));
 
       std::vector<double> coefficientOfFriction;
       coefficientOfFriction.resize(grid.size(grid.maxLevel(), dim));
       std::fill(coefficientOfFriction.begin(), coefficientOfFriction.end(),
                 parset.get<double>("boundary.mu"));
 
-      std::vector<double> eta;
-      eta.resize(grid.size(grid.maxLevel(), dim));
-      std::fill(eta.begin(), eta.end(), parset.get<double>("boundary.eta"));
-
       std::vector<double> normalStress;
       normalStress.resize(grid.size(grid.maxLevel(), dim));
       std::fill(normalStress.begin(), normalStress.end(),
                 parset.get<double>("boundary.normalstress"));
 
-      Dune::GlobalLaursenNonlinearity<dim, Dune::LinearFunction>
-      myGlobalNonlinearity(coefficientOfFriction, normalStress, nodalIntegrals);
-      // Dune::GlobalRuinaNonlinearity<dim>
-      //   myGlobalNonlinearity(nodalIntegrals, a, coefficientOfFriction, eta,
-      // normalStress);
+      Dune::GlobalNonlinearity<dim> *myGlobalNonlinearity;
+
+      std::string const friction_model =
+          parset.get<std::string>("boundary.friction.model");
+      if (friction_model == std::string("Ruina")) {
+        std::vector<double> a;
+        a.resize(grid.size(grid.maxLevel(), dim));
+        std::fill(a.begin(), a.end(), parset.get<double>("boundary.a"));
+
+        std::vector<double> eta;
+        eta.resize(grid.size(grid.maxLevel(), dim));
+        std::fill(eta.begin(), eta.end(), parset.get<double>("boundary.eta"));
+
+        myGlobalNonlinearity = new Dune::GlobalRuinaNonlinearity<dim>(
+            nodalIntegrals, a, coefficientOfFriction, eta, normalStress);
+      } else if (friction_model == std::string("Laursen")) {
+        myGlobalNonlinearity =
+            new Dune::GlobalLaursenNonlinearity<dim, Dune::LinearFunction>(
+                coefficientOfFriction, normalStress, nodalIntegrals);
+      } else {
+        assert(false);
+      }
       // }}}
 
       {
         MyConvexProblemType myConvexProblem(stiffnessMatrix,
-                                            myGlobalNonlinearity, b1, u1);
+                                            *myGlobalNonlinearity, b1, u1);
         MyBlockProblemType myBlockProblem(myConvexProblem);
         nonlinearGSStep.setProblem(u1, myBlockProblem);
 
@@ -251,6 +261,7 @@ int main(int argc, char *argv[]) {
                                       solver_tolerance, &energyNorm,
                                       Solver::QUIET);
         solver.solve();
+        free(myGlobalNonlinearity);
       }
 
       auto *displacement =

src/one-body-sample.parset

@@ -15,4 +15,7 @@ tolerance = 1e-6
 a = 0.0015
 normalstress = 0.1
 mu = 0.75
-eta = 1
\ No newline at end of file
+eta = 1
+
+[boundary.friction]
+model = Laursen