diff --git a/dune/tectonic/globallaursennonlinearity.hh b/dune/tectonic/globallaursennonlinearity.hh
index 4b24e83c1a7c6e96ca0c955c76507086eedb16f1..6d9d1ddcc023d59b561940acc1fa39f2a1c477fc 100644
--- a/dune/tectonic/globallaursennonlinearity.hh
+++ b/dune/tectonic/globallaursennonlinearity.hh
@@ -7,6 +7,7 @@
#include <dune/common/fvector.hh>
+#include "localnonlinearity.hh"
#include "globalnonlinearity.hh"
#include "nicefunction.hh"
@@ -33,11 +34,15 @@ class GlobalLaursenNonlinearity : public Dune::GlobalNonlinearity<dim> {
sigma_n [id + mu id] = sigma_n (1 + mu) id
*/
- virtual Dune::NiceFunction *restriction(int i) const {
+ virtual shared_ptr<LocalNonlinearity<dim> const> restriction(int i) const {
double coefficient = nodalIntegrals[i][0];
coefficient *= normalStress[i];
coefficient *= 1 + mu[i];
- return new OuterFunctionType(coefficient);
+
+ shared_ptr<NiceFunction const> const func(
+ new OuterFunctionType(coefficient));
+ return shared_ptr<LocalNonlinearity<dim> const>(
+ new LocalNonlinearity<dim>(func));
}
private:
diff --git a/dune/tectonic/globalnonlinearity.hh b/dune/tectonic/globalnonlinearity.hh
index 7915f9ae636d09c4ff3de0bb973842b2e02e4f62..15d15da8a19bed63f2e835e519e8b85353feafef 100644
--- a/dune/tectonic/globalnonlinearity.hh
+++ b/dune/tectonic/globalnonlinearity.hh
@@ -6,14 +6,15 @@
#include <dune/common/fvector.hh>
#include "nicefunction.hh"
+#include "localnonlinearity.hh"
namespace Dune {
template <int dim> class GlobalNonlinearity {
public:
/*
- Return a restriction of the outer function to the i'th node. If
+ Return a restriction of the outer function to the i'th node.
*/
- virtual NiceFunction* restriction(int i) const = 0;
+ virtual shared_ptr<LocalNonlinearity<dim> const> restriction(int i) const = 0;
};
}
#endif
diff --git a/dune/tectonic/globalruinanonlinearity.hh b/dune/tectonic/globalruinanonlinearity.hh
index a2b2cb68dbb7c19a746a46917bc6991cb27df5df..249be621e63271a939615082aff8cd55f7ea23a9 100644
--- a/dune/tectonic/globalruinanonlinearity.hh
+++ b/dune/tectonic/globalruinanonlinearity.hh
@@ -22,20 +22,25 @@ class GlobalRuinaNonlinearity : public Dune::GlobalNonlinearity<dim> {
a(a),
mu(mu),
eta(eta),
- normalStress(normalStress) {}
+ normalStress(normalStress),
+ trivialNonlinearity(new LocalNonlinearity<dim>(
+ shared_ptr<NiceFunction const>(new TrivialFunction))) {}
/*
Return a restriction of the outer function to the i'th node.
*/
- virtual Dune::NiceFunction *restriction(int i) const {
+ virtual shared_ptr<LocalNonlinearity<dim> const> restriction(int i) const {
if (nodalIntegrals[i][0] == 0)
- return new Dune::TrivialFunction();
- else
- return new Dune::RuinaFunction(nodalIntegrals[i][0], a[i], mu[i], eta[i],
- normalStress[i]);
+ return trivialNonlinearity;
+
+ shared_ptr<NiceFunction const> const func(new Dune::RuinaFunction(
+ nodalIntegrals[i][0], a[i], mu[i], eta[i], normalStress[i]));
+ return shared_ptr<LocalNonlinearity<dim>>(new LocalNonlinearity<dim>(func));
}
private:
+ shared_ptr<LocalNonlinearity<dim> const> const trivialNonlinearity;
+
std::vector<Dune::FieldVector<double, 1>> nodalIntegrals;
std::vector<double> a;
std::vector<double> mu;
diff --git a/dune/tectonic/localnonlinearity.hh b/dune/tectonic/localnonlinearity.hh
index fbdd63257359a3df461d7396b48a6c8841be9f1e..0d5b5b78c5d45745f81f4b33254b2b13b6924ba6 100644
--- a/dune/tectonic/localnonlinearity.hh
+++ b/dune/tectonic/localnonlinearity.hh
@@ -18,7 +18,7 @@ template <int dimension> class LocalNonlinearity {
typedef FieldVector<double, dimension> VectorType;
typedef FieldMatrix<double, dimension, dimension> MatrixType;
- LocalNonlinearity(Dune::shared_ptr<NiceFunction> &func) : func_(func) {}
+ LocalNonlinearity(Dune::shared_ptr<NiceFunction const> func) : func_(func) {}
double operator()(VectorType const x) const {
double ret;
@@ -63,7 +63,7 @@ template <int dimension> class LocalNonlinearity {
}
private:
- Dune::shared_ptr<NiceFunction> &func_;
+ Dune::shared_ptr<NiceFunction const> const func_;
};
}
#endif
diff --git a/dune/tectonic/myblockproblem.hh b/dune/tectonic/myblockproblem.hh
index dd085f9a4825519a51217b1361fc2843dacb65dd..cb96d7d0a7f6aa98fc5a8e6b6e097de31ca12d0b 100644
--- a/dune/tectonic/myblockproblem.hh
+++ b/dune/tectonic/myblockproblem.hh
@@ -9,7 +9,6 @@
#include <dune/tnnmg/problem-classes/onedconvexfunction.hh>
#include "localnonlinearity.hh"
-#include "nicefunction.hh"
#include "samplefunctional.hh"
/** \brief Base class for problems where each block can be solved with a
@@ -118,8 +117,7 @@ class MyBlockProblem<MyConvexProblemTypeTEMPLATE>::IterateObject {
}
assert(localA != NULL);
- Dune::shared_ptr<Dune::NiceFunction> f(problem.phi.restriction(m));
- Dune::LocalNonlinearity<block_size> const phi(f);
+ auto phi = problem.phi.restriction(m);
Dune::SampleFunctional<block_size> localJ(*localA, localb, phi,
ignore_component);
diff --git a/dune/tectonic/samplefunctional.hh b/dune/tectonic/samplefunctional.hh
index f5de28a3311d9d2b60932987df095cce11f15f48..d042d9e7edaffa673707932d52b68689946671f6 100644
--- a/dune/tectonic/samplefunctional.hh
+++ b/dune/tectonic/samplefunctional.hh
@@ -24,15 +24,15 @@ template <int dim> class SampleFunctional {
typedef LocalNonlinearity<dim> NonlinearityType;
SampleFunctional(SmallMatrix const &A, SmallVector const &b,
- NonlinearityType const &phi, int ignore = dim)
+ shared_ptr<NonlinearityType const> phi, int ignore = dim)
: A(A), b(b), phi(phi), ignore(ignore) {}
double operator()(SmallVector const &v) const {
SmallVector y;
- A.mv(v, y); // Av
- y /= 2; // 1/2 Av
- y -= b; // 1/2 Av - b
- return y * v + phi(v); // <1/2 Av - b,v> + H(|v|)
+ A.mv(v, y); // Av
+ y /= 2; // 1/2 Av
+ y -= b; // 1/2 Av - b
+ return y * v + *phi(v); // <1/2 Av - b,v> + H(|v|)
}
bool descentDirection(SmallVector const x, SmallVector &ret) const {
@@ -44,7 +44,7 @@ template <int dim> class SampleFunctional {
smoothGradient(x, d);
Interval<double> D;
- phi.directionalSubDiff(x, d, D);
+ phi->directionalSubDiff(x, d, D);
double const nonlinearDecline = D[1];
double const smoothDecline = -(d * d);
double const combinedDecline = smoothDecline + nonlinearDecline;
@@ -94,7 +94,7 @@ template <int dim> class SampleFunctional {
SmallMatrix const &A;
SmallVector const &b;
- NonlinearityType const φ
+ shared_ptr<NonlinearityType const> const phi;
int const ignore;
private:
@@ -107,7 +107,7 @@ template <int dim> class SampleFunctional {
}
void upperGradient(SmallVector const &x, SmallVector &y) const {
- phi.upperGradient(x, y);
+ phi->upperGradient(x, y);
SmallVector z;
smoothGradient(x, z);
y += z;
@@ -116,7 +116,7 @@ template <int dim> class SampleFunctional {
}
void lowerGradient(SmallVector const &x, SmallVector &y) const {
- phi.lowerGradient(x, y);
+ phi->lowerGradient(x, y);
SmallVector z;
smoothGradient(x, z);
y += z;
@@ -156,7 +156,7 @@ void minimise(Functional const J, typename Functional::SmallVector &x,
return;
typedef typename Functional::NonlinearityType LocalNonlinearityType;
- LocalNonlinearityType phi = J.phi;
+ LocalNonlinearityType phi = *J.phi; // copy, see below
if (linesearchp) {
// {{{ Construct a restriction of J to the line x + t * descDir