diff --git a/dune/tectonic/nicefunction.hh b/dune/tectonic/nicefunction.hh index 626c59e60a027448dd7ae033bc1007cd36d5b3fa..93234e1d8f668f58e4c01d565d57180e6bacddf7 100644 --- a/dune/tectonic/nicefunction.hh +++ b/dune/tectonic/nicefunction.hh @@ -56,11 +56,13 @@ class RuinaFunction : public NiceFunction { rho = exp(-mu/a) */ void virtual evaluate(double const &x, double &y) const { - double const arg = std::max(eta * x / h, rho); - double const r = arg * (a * (std::log(arg) - 1) + mu + compound_state); - y = 1 / eta * normalStress * (r + a * rho + mu + compound_state); - y *= coefficient; - y *= h; + double const arg = eta * x / h; + if (arg <= rho) { + y = -log(rho); // TODO: We can write this out explicitly + return; + } + y = arg * (log(arg) - 1) - arg * log(rho) + rho - log(rho); + y *= coefficient * h * normalStress / eta; } /*