Skip to content
GitLab
Commit d3d6f0c9 authored by Martin Utz's avatar Martin Utz
Browse files

[swe][frictionlaws] Restructure source calculation 

The source is now entirely calculated within the friction law.
parent f563b6ae
Hide whitespace changes
Inline Side-by-side
dumux/material/fluidmatrixinteractions/frictionlaws/frictionlaw.hh
View file @ d3d6f0c9
......@@ -35,17 +35,22 @@ namespace Dumux {
* The LET mobility model is used to limit the friction for small water depths.
*/
template <typename Scalar>
template <typename Scalar, typename NumEqVector>
class FrictionLaw
{
public:
/*!
* \brief Compute the friction ustar_h.
* \brief Compute the friction source term.
*
* \return ustar_h friction used for the source term in shallow water models.
* The friction source term contains the losses due to friction.
*
* \return Friction source term.
*/
virtual Scalar computeUstarH(const Scalar waterDept, const Scalar frictionValue) const = 0;
virtual NumEqVector computeSource(const Scalar waterDept,
const Scalar frictionValue,
const Scalar u,
const Scalar v) const = 0;
/*!
* \brief Limit the friction for small water depth.
......@@ -81,12 +86,12 @@ public:
using std::min;
using std::max;
Scalar mobility_max = 1.0; //!< maximal mobility
Scalar mobilityMax = 1.0; //!< maximal mobility
Scalar minUpperH = roughnessHeight * 2.0;
Scalar sw = min(waterDepth * (1.0/minUpperH),1.0);
sw = max(0.0,sw);
auto mobility = mobility_max /(1 + (1.0-sw)*(1.0-sw));
auto mobility = mobilityMax /(1 + (1.0-sw)*(1.0-sw));
return roughnessHeight * (1.0 - mobility);
}
};
......
dumux/material/fluidmatrixinteractions/frictionlaws/manning.hh
View file @ d3d6f0c9
......@@ -36,32 +36,43 @@ namespace Dumux {
* The LET mobility model is used to limit the friction for small water depths.
*/
template <typename Scalar>
class FrictionLawManning : public FrictionLaw<Scalar>
template <typename Scalar, typename NumEqVector>
class FrictionLawManning : public FrictionLaw<Scalar, NumEqVector>
{
public:
FrictionLawManning(const Scalar gravity)
: gravity_(gravity) {}
/*!
* \brief Compute the friction ustar_h.
* \brief Compute the friction source term.
*
* \param waterDepth water depth.
* \param frictionValue The Manning friction value.
* \param frictionValue Manning friction value.
* \param u velocity in x-direction.
* \param v velocity in y-direction.
*
* \return ustar_h friction used for the source term in shallow water models.
* \return Friction source term.
*/
Scalar computeUstarH(const Scalar waterDepth, const Scalar frictionValue) const final
NumEqVector computeSource(const Scalar waterDepth,
const Scalar frictionValue,
const Scalar u,
const Scalar v) const final
{
using std::pow;
using std::hypot;
NumEqVector source(0.0);
Scalar ustar_h = 0.0;
Scalar rough_h = pow(25.68/(1.0/frictionValue),6.0);
Scalar roughnessHeight = pow(25.68/(1.0/frictionValue),6.0);
roughnessHeight = this->limitRoughH(roughnessHeight, waterDepth);
const Scalar c = pow((waterDepth + roughnessHeight),1.0/6.0) * 1.0/(frictionValue);
const Scalar uv = hypot(u,v);
rough_h = this->limitRoughH(rough_h, waterDepth);
source[0] = 0.0;
source[1] = -gravity_/(c*c) * u * uv;
source[2] = -gravity_/(c*c) * v * uv;
auto cfric = pow((waterDepth + rough_h),1.0/6.0) * 1.0/(frictionValue);
ustar_h = gravity_ / pow(cfric,2.0);
return ustar_h;
return source;
}
private:
Scalar gravity_;
......
dumux/material/fluidmatrixinteractions/frictionlaws/nikuradse.hh
View file @ d3d6f0c9
......@@ -36,35 +36,43 @@ namespace Dumux {
* The LET mobility model is used to limit the friction for small water depths.
*/
template <typename Scalar>
class FrictionLawNikuradse : public FrictionLaw<Scalar>
template <typename Scalar, typename NumEqVector>
class FrictionLawNikuradse : public FrictionLaw<Scalar, NumEqVector>
{
public:
FrictionLawNikuradse(const Scalar gravity)
: gravity_(gravity) {}
/*!
* \brief Compute the friction ustar_h.
* \brief Compute the friction source term.
*
* \param waterDepth water depth.
* \param frictionValue The equivalent sand roughness.
* \param u velocity in x-direction.
* \param v velocity in y-direction.
*
* \return ustar_h friction used for the source term in shallow water models.
* \return Friction source term.
*/
Scalar computeUstarH(const Scalar waterDepth, const Scalar frictionValue) const final
NumEqVector computeSource(const Scalar waterDepth,
const Scalar frictionValue,
const Scalar u,
const Scalar v) const final
{
using std::pow;
using std::log;
using std::hypot;
Scalar ustar_h = 0.0;
Scalar rough_h = frictionValue;
NumEqVector source(0.0);
rough_h = this->limitRoughH(rough_h, waterDepth);
ustar_h = pow(0.41,2.0)/pow(log((12*(waterDepth + rough_h))/frictionValue),2.0);
return ustar_h;
Scalar roughnessHeight = frictionValue;
roughnessHeight = this->limitRoughH(roughnessHeight, waterDepth);
const Scalar ustarH = pow(0.41,2.0)/pow(log((12*(waterDepth + roughnessHeight))/frictionValue),2.0);
const Scalar uv = hypot(u,v);
source[0] = 0.0;
source[1] = -ustarH * u * uv;
source[2] = -ustarH * v * uv;
return source;
}
private:
Scalar gravity_;
};
} // end namespace Dumux
......
test/freeflow/shallowwater/roughchannel/problem.hh
View file @ d3d6f0c9
......@@ -147,11 +147,11 @@ public:
const auto gravity = this->spatialParams().gravity();
if (getParam<std::string>("Problem.FrictionLaw") == ("Manning"))
{
frictionLaw_ = std::make_shared<FrictionLawManning<Scalar>>(gravity);
frictionLaw_ = std::make_shared<FrictionLawManning<Scalar,NumEqVector>>(gravity);
}
else if (getParam<std::string>("Problem.FrictionLaw") == ("Nikuradse"))
{
frictionLaw_ = std::make_shared<FrictionLawNikuradse<Scalar>>(gravity);
frictionLaw_ = std::make_shared<FrictionLawNikuradse<Scalar,NumEqVector>>();
std::cout<<"\nWARNING: This test is meant to be run for the friction law after Manning.\n";
std::cout<<" You are running it with Nikuradse, although the friction values in\n";
std::cout<<" the input file and the analytic solution don't fit to Nikuradse!\n\n.";
......@@ -240,8 +240,6 @@ public:
{
using std::hypot;
NumEqVector source(0.0);
const auto& globalPos = scv.center();
const auto& volVars = elemVolVars[scv];
const Scalar manningN = this->spatialParams().frictionValue(globalPos);
......@@ -249,12 +247,7 @@ public:
const Scalar h = volVars.waterDepth();
const Scalar u = volVars.velocity(0);
const Scalar v = volVars.velocity(1);
const Scalar ustarH = frictionLaw_->computeUstarH(h, manningN);
const Scalar uv = hypot(u,v);
source[0] = 0.0;
source[1] = -ustarH * u * uv;
source[2] = -ustarH * v * uv;
NumEqVector source = frictionLaw_->computeSource(h, manningN, u, v);
return source;
}
......@@ -387,7 +380,7 @@ private:
Scalar constManningN_; // analytic solution is only available for const friction.
Scalar bedSlope_;
Scalar discharge_; // discharge at the inflow boundary
std::shared_ptr<FrictionLaw<Scalar>> frictionLaw_;
std::shared_ptr<FrictionLaw<Scalar,NumEqVector>> frictionLaw_;
static constexpr Scalar eps_ = 1.0e-6;
std::string name_;
};
......
test/freeflow/shallowwater/roughchannel/spatialparams.hh
View file @ d3d6f0c9
......@@ -76,6 +76,11 @@ public:
/*! \brief Define the friction value.
*
* The unit of the friction value depends on the used friction law.
* The Mannng friction value is just a coefficient and has no unit $[-]$
* In the Nikuradse law the friction value is the equivalent sand roughness
* with the unit $[m]$.
*
* \return friction value
*/
Scalar frictionValue(const GlobalPosition& globalPos) const
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment