[Cleanup] Move weights into normal stress

dune/tectonic/frictionpotential.hh

@@ -29,8 +29,9 @@ class FrictionPotential {
 
 class TruncatedRateState : public FrictionPotential {
 public:
-  TruncatedRateState(double coefficient, double _normalStress, FrictionData _fd)
-      : fd(_fd), weight(coefficient), normalStress(_normalStress) {}
+  TruncatedRateState(double _weight, double _weightedNormalStress,
+                     FrictionData _fd)
+      : fd(_fd), weight(_weight), weightedNormalStress(_weightedNormalStress) {}
 
   double coefficientOfFriction(double V) const {
     if (V <= Vmin)
@@ -40,14 +41,14 @@ class TruncatedRateState : public FrictionPotential {
   }
 
   double differential(double V) const {
-    return weight * (fd.C - normalStress * coefficientOfFriction(V));
+    return weight * fd.C - weightedNormalStress * coefficientOfFriction(V);
   }
 
   double second_deriv(double V) const {
     if (V <= Vmin)
       return 0;
 
-    return weight * (-normalStress) * (fd.a / V);
+    return -weightedNormalStress * (fd.a / V);
   }
 
   double regularity(double V) const {
@@ -65,26 +66,26 @@ class TruncatedRateState : public FrictionPotential {
 private:
   FrictionData const fd;
   double const weight;
-  double const normalStress;
+  double const weightedNormalStress;
   double Vmin;
 };
 
 class RegularisedRateState : public FrictionPotential {
 public:
-  RegularisedRateState(double coefficient, double _normalStress,
+  RegularisedRateState(double _weight, double _weightedNormalStress,
                        FrictionData _fd)
-      : fd(_fd), weight(coefficient), normalStress(_normalStress) {}
+      : fd(_fd), weight(_weight), weightedNormalStress(_weightedNormalStress) {}
 
   double coefficientOfFriction(double V) const {
     return fd.a * std::asinh(0.5 * V / Vmin);
   }
 
   double differential(double V) const {
-    return weight * (fd.C - normalStress * coefficientOfFriction(V));
+    return weight * fd.C - weightedNormalStress * coefficientOfFriction(V);
   }
 
   double second_deriv(double V) const {
-    return weight * (-normalStress) * fd.a / std::hypot(2.0 * Vmin, V);
+    return -weightedNormalStress * fd.a / std::hypot(2.0 * Vmin, V);
   }
 
   double regularity(double V) const { return std::abs(second_deriv(V)); }
@@ -97,7 +98,7 @@ class RegularisedRateState : public FrictionPotential {
 private:
   FrictionData const fd;
   double const weight;
-  double const normalStress;
+  double const weightedNormalStress;
   double Vmin;
 };
 

dune/tectonic/globalratestatefriction.hh

@@ -26,10 +26,9 @@ class GlobalRateStateFriction : public GlobalFriction<Matrix, Vector> {
   GlobalRateStateFriction(BoundaryPatch<GridView> const &frictionalBoundary,
                           GridView const &gridView,
                           GlobalFrictionData<block_size> const &frictionInfo,
-                          // Note: passing the following two makes no sense
                           ScalarVector const &weights,
-                          ScalarVector const &normalStress)
-      : restrictions(normalStress.size()) {
+                          ScalarVector const &weightedNormalStress)
+      : restrictions(weightedNormalStress.size()) {
     auto trivialNonlinearity =
         std::make_shared<Friction>(std::make_shared<TrivialFunction>());
 
@@ -44,7 +43,7 @@ class GlobalRateStateFriction : public GlobalFriction<Matrix, Vector> {
         continue;
       }
       auto const fp = std::make_shared<ScalarFriction>(
-          weights[i], normalStress[i], frictionInfo(coordinate));
+          weights[i], weightedNormalStress[i], frictionInfo(coordinate));
       restrictions[i] = std::make_shared<Friction>(fp);
     }
   }

src/assemblers.cc

@@ -142,15 +142,21 @@ auto MyAssembler<GridView, dimension>::assembleFrictionNonlinearity(
     vertexAssembler.assembleBoundaryFunctional(frictionalBoundaryAssembler,
                                                weights, frictionalBoundary);
   }
+  ScalarVector weightedNormalStress(weights.size());
+  for (size_t i = 0; i < weights.size(); ++i)
+    weightedNormalStress[i] = weights[i] * normalStress[i];
+
   switch (frictionModel) {
     case Config::Truncated:
       return std::make_shared<GlobalRateStateFriction<
           Matrix, Vector, TruncatedRateState, GridView>>(
-          frictionalBoundary, gridView, frictionInfo, weights, normalStress);
+          frictionalBoundary, gridView, frictionInfo, weights,
+          weightedNormalStress);
     case Config::Regularised:
       return std::make_shared<GlobalRateStateFriction<
           Matrix, Vector, RegularisedRateState, GridView>>(
-          frictionalBoundary, gridView, frictionInfo, weights, normalStress);
+          frictionalBoundary, gridView, frictionInfo, weights,
+          weightedNormalStress);
     default:
       assert(false);
   }