Have nonlinearities propagate their kinks

dune/tectonic/ellipticenergy.hh

@@ -21,6 +21,8 @@ template <int dim> class EllipticEnergy {
                  shared_ptr<NonlinearityType const> phi, int ignore = dim)
       : A(A), b(b), phi(phi), ignore(ignore) {}
 
+  std::vector<double> const &get_kinks() const { return phi->get_kinks(); }
+
   double operator()(SmallVector const &v) const {
     SmallVector y(0);
     A.usmv(0.5, v, y);        //  1/2 Av

dune/tectonic/localnonlinearity.hh

@@ -20,6 +20,8 @@ template <int dimension> class LocalNonlinearity {
 
   LocalNonlinearity(shared_ptr<NiceFunction const> func) : func(func) {}
 
+  std::vector<double> const &get_kinks() const { return func->get_kinks(); }
+
   double operator()(VectorType const &x) const {
     return func->evaluate(x.two_norm());
   }

dune/tectonic/minimisation.hh

@@ -147,8 +147,7 @@ void minimisationInitialiser(Functional const &J, Bisection const &bisection,
                              typename Functional::SmallVector &startingPoint) {
   typedef typename Functional::SmallVector SmallVector;
 
-  std::vector<double> const kinks = { 5, 10,
-                                      15 }; // FIXME: We happen to know these
+  std::vector<double> const &kinks = J.get_kinks();
 
   SmallVector x_old(0);
   double Jx_old = J(x_old);

dune/tectonic/nicefunction.hh

@@ -10,7 +10,17 @@
 
 namespace Dune {
 class NiceFunction {
+protected:
+  NiceFunction(std::vector<double> const &_kinks) : kinks(_kinks) {}
+
+private:
+  std::vector<double> const kinks;
+
 public:
+  std::vector<double> const &get_kinks() const { return kinks; }
+
+  NiceFunction() : kinks() {}
+
   virtual ~NiceFunction() {}
 
   double virtual leftDifferential(double s) const = 0;
@@ -32,7 +42,8 @@ class RuinaFunction : public NiceFunction {
 public:
   RuinaFunction(double coefficient, double a, double mu, double eta,
                 double normalStress, double b, double state, double L, double h)
-      : a(a),
+      : NiceFunction(),
+        a(a),
         eta(eta),
         h(h),
         coefficientProduct(coefficient * normalStress),
@@ -83,6 +94,8 @@ class RuinaFunction : public NiceFunction {
 
 class TrivialFunction : public NiceFunction {
 public:
+  TrivialFunction() : NiceFunction() {}
+
   double virtual evaluate(double) const { return 0; }
 
   double virtual leftDifferential(double) const { return 0; }

src/test-gradient-method-nicefunction.hh

@@ -7,6 +7,12 @@
 
 namespace Dune {
 class MyFunction : public NiceFunction {
+protected:
+  MyFunction(std::vector<double> const &_kinks) : NiceFunction(_kinks) {}
+
+public:
+  MyFunction() : NiceFunction() {}
+
   double virtual second_deriv(double) const {
     assert(false);
     return 0;
@@ -34,7 +40,7 @@ class Parabola : public MyFunction {
 
 class LinearFunction : public MyFunction {
 public:
-  LinearFunction(double a) : coefficient(a) {}
+  LinearFunction(double a) : MyFunction(), coefficient(a) {}
 
   double virtual evaluate(double x) const { return coefficient * x; }
 
@@ -52,6 +58,8 @@ class LinearFunction : public MyFunction {
 
 template <int slope> class SampleFunction : public MyFunction {
 public:
+  SampleFunction() : MyFunction({ 1 }) {}
+
   double virtual evaluate(double x) const {
     return (x < 1) ? x : (slope * (x - 1) + 1);
   }
@@ -76,6 +84,8 @@ class SteepFunction : public MyFunction {
 
 // slope in [n-1,n] is n
 class HorribleFunction : public MyFunction {
+  // TODO: Handle kinks
+
 public:
   double virtual evaluate(double x) const {
     double const fl = floor(x);
@@ -103,6 +113,8 @@ class HorribleFunction : public MyFunction {
 // slope in [n-1,n] is log(n+1)
 class HorribleFunctionLogarithmic : public MyFunction {
 public:
+  // TODO: Handle kinks
+
   double virtual evaluate(double x) const {
     double y = 0;
     size_t const fl = floor(x);
@@ -131,14 +143,14 @@ class HorribleFunctionLogarithmic : public MyFunction {
 };
 
 class ThreeKinkFunction : public MyFunction {
-private:
-  std::vector<double> kinks = { 5, 10, 15 };
-
 public:
+  ThreeKinkFunction() : MyFunction({ 5, 10, 15 }) {}
+
   double virtual evaluate(double x) const {
     double acc = 0;
     double last_kink = 0;
     int i;
+    auto const &kinks = this->get_kinks();
     for (i = 0; i < kinks.size(); ++i) {
       if (x <= kinks[i])
         break;