Skip to content
Snippets Groups Projects
Commit 7485f0a2 authored by Timo Koch's avatar Timo Koch
Browse files

[fluidmatrixinteraction] Use std::clamp for improved readability

parent 77d04e3e
No related branches found
No related tags found
1 merge request!1992[fluidmatrixinteraction] Use std::clamp for improved readability
...@@ -72,10 +72,9 @@ public: ...@@ -72,10 +72,9 @@ public:
static Scalar pc(const Params &params, Scalar swe) static Scalar pc(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
return params.pe()*pow(swe, -1.0/params.lambda()); return params.pe()*pow(swe, -1.0/params.lambda());
} }
...@@ -136,10 +135,9 @@ public: ...@@ -136,10 +135,9 @@ public:
static Scalar dpc_dswe(const Params &params, Scalar swe) static Scalar dpc_dswe(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
return - params.pe()/params.lambda() * pow(swe, -1/params.lambda() - 1); return - params.pe()/params.lambda() * pow(swe, -1/params.lambda() - 1);
} }
...@@ -184,10 +182,9 @@ public: ...@@ -184,10 +182,9 @@ public:
static Scalar krw(const Params &params, Scalar swe) static Scalar krw(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
return pow(swe, 2.0/params.lambda() + 3); return pow(swe, 2.0/params.lambda() + 3);
} }
...@@ -210,10 +207,9 @@ public: ...@@ -210,10 +207,9 @@ public:
static Scalar dkrw_dswe(const Params &params, Scalar swe) static Scalar dkrw_dswe(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
return (2.0/params.lambda() + 3)*pow(swe, 2.0/params.lambda() + 2); return (2.0/params.lambda() + 3)*pow(swe, 2.0/params.lambda() + 2);
} }
...@@ -235,10 +231,9 @@ public: ...@@ -235,10 +231,9 @@ public:
static Scalar krn(const Params &params, Scalar swe) static Scalar krn(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const Scalar exponent = 2.0/params.lambda() + 1; const Scalar exponent = 2.0/params.lambda() + 1;
const Scalar tmp = 1.0 - swe; const Scalar tmp = 1.0 - swe;
...@@ -264,10 +259,9 @@ public: ...@@ -264,10 +259,9 @@ public:
static Scalar dkrn_dswe(const Params &params, Scalar swe) static Scalar dkrn_dswe(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const auto lambdaInv = 1.0/params.lambda(); const auto lambdaInv = 1.0/params.lambda();
const auto swePow = pow(swe, 2*lambdaInv); const auto swePow = pow(swe, 2*lambdaInv);
......
...@@ -144,9 +144,8 @@ public: ...@@ -144,9 +144,8 @@ public:
*/ */
static Scalar krw(const Params &params, Scalar swe) static Scalar krw(const Params &params, Scalar swe)
{ {
using std::max; using std::clamp;
using std::min; return clamp(swe, 0.0, 1.0);
return max(min(swe,1.0),0.0);
} }
/*! /*!
...@@ -160,10 +159,8 @@ public: ...@@ -160,10 +159,8 @@ public:
*/ */
static Scalar krn(const Params &params, Scalar swe) static Scalar krn(const Params &params, Scalar swe)
{ {
Scalar sne = 1 - swe; using std::clamp;
using std::max; return clamp(1.0-swe, 0.0, 1.0);
using std::min;
return max(min(sne,1.0),0.0);
} }
}; };
} // end namespace Dumux } // end namespace Dumux
......
...@@ -72,10 +72,9 @@ public: ...@@ -72,10 +72,9 @@ public:
static Scalar pc(const Params &params, Scalar swe) static Scalar pc(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const Scalar pc = pow(pow(swe, -1.0/params.vgm()) - 1, 1.0/params.vgn())/params.vgAlpha(); const Scalar pc = pow(pow(swe, -1.0/params.vgm()) - 1, 1.0/params.vgn())/params.vgAlpha();
return pc; return pc;
...@@ -142,10 +141,9 @@ public: ...@@ -142,10 +141,9 @@ public:
static Scalar dpc_dswe(const Params &params, Scalar swe) static Scalar dpc_dswe(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const Scalar powSwe = pow(swe, -1/params.vgm()); const Scalar powSwe = pow(swe, -1/params.vgm());
return - 1.0/params.vgAlpha() * pow(powSwe - 1, 1.0/params.vgn() - 1)/params.vgn() return - 1.0/params.vgAlpha() * pow(powSwe - 1, 1.0/params.vgn() - 1)/params.vgn()
...@@ -193,10 +191,9 @@ public: ...@@ -193,10 +191,9 @@ public:
static Scalar krw(const Params &params, Scalar swe) static Scalar krw(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const Scalar r = 1.0 - pow(1.0 - pow(swe, 1.0/params.vgm()), params.vgm()); const Scalar r = 1.0 - pow(1.0 - pow(swe, 1.0/params.vgm()), params.vgm());
return pow(swe, params.vgl())*r*r; return pow(swe, params.vgl())*r*r;
...@@ -219,10 +216,9 @@ public: ...@@ -219,10 +216,9 @@ public:
static Scalar dkrw_dswe(const Params &params, Scalar swe) static Scalar dkrw_dswe(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const Scalar x = 1.0 - pow(swe, 1.0/params.vgm()); const Scalar x = 1.0 - pow(swe, 1.0/params.vgm());
const Scalar xToM = pow(x, params.vgm()); const Scalar xToM = pow(x, params.vgm());
...@@ -242,15 +238,14 @@ public: ...@@ -242,15 +238,14 @@ public:
* *
* \note Instead of undefined behaviour if pc is not in the valid range, we return a valid number, * \note Instead of undefined behaviour if pc is not in the valid range, we return a valid number,
* by clamping the input. * by clamping the input.
* \note See e.g. Dury, Fischer, Schulin (1999) for application of Mualem model to non-wetting rel. perm. * \note See e.g. Dury, Fischer, Schulin (1999) for application of Mualem model to non-wetting rel. perm.
*/ */
static Scalar krn(const Params &params, Scalar swe) static Scalar krn(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
return pow(1 - swe, params.vgl()) * pow(1 - pow(swe, 1.0/params.vgm()), 2*params.vgm()); return pow(1 - swe, params.vgl()) * pow(1 - pow(swe, 1.0/params.vgm()), 2*params.vgm());
} }
...@@ -273,10 +268,9 @@ public: ...@@ -273,10 +268,9 @@ public:
static Scalar dkrn_dswe(const Params &params, Scalar swe) static Scalar dkrn_dswe(const Params &params, Scalar swe)
{ {
using std::pow; using std::pow;
using std::min; using std::clamp;
using std::max;
swe = min(max(swe, 0.0), 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0 swe = clamp(swe, 0.0, 1.0); // the equation below is only defined for 0.0 <= sw <= 1.0
const auto sne = 1.0 - swe; const auto sne = 1.0 - swe;
const auto x = 1.0 - pow(swe, 1.0/params.vgm()); const auto x = 1.0 - pow(swe, 1.0/params.vgm());
......
...@@ -235,13 +235,12 @@ public: ...@@ -235,13 +235,12 @@ public:
krn -= pow(1 - pow(ste, 1/params.vgm()), params.vgm()); krn -= pow(1 - pow(ste, 1/params.vgm()), params.vgm());
krn *= krn; krn *= krn;
using std::max; using std::clamp;
using std::min;
using std::sqrt; using std::sqrt;
if (params.krRegardsSnr()) if (params.krRegardsSnr())
{ {
// regard Snr in the permeability of the n-phase, see Helmig1997 // regard Snr in the permeability of the n-phase, see Helmig1997
Scalar resIncluded = max(min((sn - params.snr()/ (1-params.swr())), 1.0), 0.0); const Scalar resIncluded = clamp(sn - params.snr()/ (1-params.swr()), 0.0, 1.0);
krn *= sqrt(resIncluded ); krn *= sqrt(resIncluded );
} }
else else
......
...@@ -25,6 +25,8 @@ ...@@ -25,6 +25,8 @@
#ifndef REGULARIZED_PARKERVANGEN_3P_HH #ifndef REGULARIZED_PARKERVANGEN_3P_HH
#define REGULARIZED_PARKERVANGEN_3P_HH #define REGULARIZED_PARKERVANGEN_3P_HH
#include <algorithm>
#include "parkervangen3p.hh" #include "parkervangen3p.hh"
#include "regularizedparkervangen3pparams.hh" #include "regularizedparkervangen3pparams.hh"
...@@ -93,13 +95,9 @@ public: ...@@ -93,13 +95,9 @@ public:
static Scalar pcgw(const Params &params, Scalar swe) static Scalar pcgw(const Params &params, Scalar swe)
{ {
//if specified, a constant value is used for regularization //if specified, a constant value is used for regularization
if(params.constRegularization()) using std::clamp;
{ if (params.constRegularization())
if(swe < 0.0) swe = clamp(swe, 0.0, 1.0);
swe = 0.0;
if(swe > 1.0)
swe = 1.0;
}
// retrieve the low and the high threshold saturations for the // retrieve the low and the high threshold saturations for the
// unregularized capillary pressure curve from the parameters // unregularized capillary pressure curve from the parameters
...@@ -333,10 +331,9 @@ public: ...@@ -333,10 +331,9 @@ public:
*/ */
static Scalar krn(const Params &params, Scalar swe, Scalar sn, Scalar ste) static Scalar krn(const Params &params, Scalar swe, Scalar sn, Scalar ste)
{ {
using std::max; using std::clamp;
using std::min; swe = clamp(swe, 0.0, 1.0);
swe = max(min(swe, 1.0), 0.0); ste = clamp(ste, 0.0, 1.0);
ste = max(min(ste, 1.0), 0.0);
if(ste - swe <= 0.0) if(ste - swe <= 0.0)
return 0.0; return 0.0;
......
...@@ -90,10 +90,9 @@ public: ...@@ -90,10 +90,9 @@ public:
const FluidState &state, const FluidState &state,
int wPhaseIdx = 0) int wPhaseIdx = 0)
{ {
using std::max; using std::clamp;
using std::min;
for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
values[phaseIdx] = max(min(state.saturation(phaseIdx),1.0),0.0); values[phaseIdx] = clamp(state.saturation(phaseIdx), 0.0, 1.0);
} }
}; };
} // end namespace Dumux } // end namespace Dumux
......
...@@ -76,9 +76,8 @@ public: ...@@ -76,9 +76,8 @@ public:
for (int dir = 0; dir < FVGridGeom::GridView::dimension; ++dir) for (int dir = 0; dir < FVGridGeom::GridView::dimension; ++dir)
divU += gradU[dir][dir]; divU += gradU[dir][dir];
using std::min; using std::clamp;
using std::max; return clamp((refPoro+divU)/(1.0+divU), minPoro, maxPoro);
return min( maxPoro, max(minPoro, (refPoro+divU)/(1.0+divU)) );
} }
/*! /*!
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment