Use the VirtualFunction class

src/finitescalarincreasingconvexfunction.hh

@@ -3,17 +3,21 @@
 #ifndef FINITE_SCALAR_INCREASING_CONVEX_FUNCTION_HH
 #define FINITE_SCALAR_INCREASING_CONVEX_FUNCTION_HH
 
+#include <dune/common/function.hh>
+
 namespace Dune {
-class FiniteScalarIncreasingConvexFunction {
+class FiniteScalarIncreasingConvexFunction
+    : public VirtualFunction<double, double> {
 public:
-  virtual double operator()(const double s) const = 0;
   virtual double leftDifferential(const double s) const = 0;
   virtual double rightDifferential(const double s) const = 0;
 };
 
 class SampleFunction : public FiniteScalarIncreasingConvexFunction {
 public:
-  double operator()(const double s) const { return (s < 1) ? s : (2 * s - 1); }
+  void evaluate(const double& x, double& y) const {
+    y = (x < 1) ? x : (2 * x - 1);
+  }
 
   double leftDifferential(const double s) const { return (s <= 1) ? 1 : 2; }
 
@@ -22,7 +26,7 @@ class SampleFunction : public FiniteScalarIncreasingConvexFunction {
 
 class TrivialFunction : public FiniteScalarIncreasingConvexFunction {
 public:
-  double operator()(const double) const { return 0; }
+  void evaluate(const double& x, double& y) const { y = 0; }
 
   double leftDifferential(const double) const { return 0; }
 

src/mynonlinearity.hh

@@ -22,7 +22,11 @@ class MyNonlinearity {
   typedef SmallVector VectorType;
   typedef SmallMatrix MatrixType;
 
-  double operator()(VectorType const x) const { return func_(x.two_norm()); }
+  double operator()(VectorType const x) const {
+    double ret;
+    func_.evaluate(x.two_norm(), ret);
+    return ret;
+  }
 
   // directional subdifferential: at u on the line u + t*v
   // u and v are assumed to be non-zero