[Transport] Use muParser in transport parameter selection

parent ec89d0ea
......@@ -2,15 +2,7 @@ volumes:
sand:
index: 0
transport:
type: superposition
eff_diff:
type: milligton_quirk_1
eff_hydromechanic_disp:
type: isotropic
long_disp:
type: const
trans_disp:
type: const
type: Deff_MQ1 + Dhm_iso
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......@@ -21,15 +13,7 @@ volumes:
silt:
index: 1
transport:
type: superposition
eff_diff:
type: milligton_quirk_1
eff_hydromechanic_disp:
type: isotropic
long_disp:
type: const
trans_disp:
type: const
type: Deff_MQ1 + Dhm_iso
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......
......@@ -19,6 +19,8 @@
#include <dune/dorie/model/transport/parameterization/hydrodynamic_dispersion/eff_diffusion/millington_quirk_2.hh>
#include <dune/dorie/model/transport/parameterization/hydrodynamic_dispersion/superposition.hh>
#include <muParser.h>
namespace Dune {
namespace Dorie {
namespace Parameterization {
......@@ -46,83 +48,111 @@ struct TransportFactory
std::shared_ptr<Transport<Traits>> hd_dips;
const auto type = type_node["type"].as<std::string>();
uint type;
try {
mu::Parser parser;
parser.DefineConst(ConstHydrodynamicDispersion<Traits>::type, 1);
parser.DefineConst(PowerLawDispersion<Traits>::type, 2);
parser.DefineConst(ConstEffectiveDiffusion<Traits>::type, 10);
parser.DefineConst(MillingtonQuirk1<Traits>::type, 20);
parser.DefineConst(MillingtonQuirk2<Traits>::type, 30);
parser.DefineConst(ConstEffectiveHydromechanicDispersion<Traits>::type, 100);
parser.DefineConst(IsotropicEffectiveHydromechanicDispersion<Traits>::type, 200);
if (type == ConstHydrodynamicDispersion<Traits>::type) {
hd_dips = std::make_shared<ConstHydrodynamicDispersion<Traits>>(name);
} else if (type == PowerLawDispersion<Traits>::type) {
hd_dips = std::make_shared<PowerLawDispersion<Traits>>(name);
} else if (type == HydrodynamicDispersionSuperposition<Traits>::type) {
parser.SetExpr(type_node["type"].as<std::string>());
type = parser.Eval();
// Select type for effective diffusion
const YAML::Node& eff_diff_node = type_node["eff_diff"];
auto eff_diff_type = eff_diff_node["type"].as<std::string>();
} catch (mu::Parser::exception_type& e) {
log->error("Evaluating transport parameterization type failed:");
log->error(" Parsed expression: {}",e.GetExpr());
log->error(" Token: {}",e.GetToken());
log->error(" Error position: {}",e.GetPos());
log->error(" Error code: {}",e.GetCode());
log->error(" Error message: {}",e.GetMsg());
DUNE_THROW(IOError, "Error evaluating transport parameterization type");
}
// get digit x from right to left of a uint.
// uint ....... xxxxxx
// ^^ ^^
// digit ...... 54..10
auto get_digit = [](uint number, int x)
{return (uint)(number*std::pow(10,-1*x))%10;};
if (type == 1) {
hd_dips = std::make_shared<ConstHydrodynamicDispersion<Traits>>(name);
} else if (type == 2) {
hd_dips = std::make_shared<PowerLawDispersion<Traits>>(name);
} else if (get_digit(type,0) == 0) {
std::shared_ptr<EffectiveDiffusion<Traits>> eff_diff;
// Select effective diffusion according to digit 1 of 'type'
uint type_1 = get_digit(type,1);
// Build parametrization for effective diffusion
if (eff_diff_type == ConstEffectiveDiffusion<Traits>::type) {
if (type_1 == 1) {
eff_diff = std::make_shared<ConstEffectiveDiffusion<Traits>>(name);
} else if (eff_diff_type == MillingtonQuirk1<Traits>::type) {
} else if (type_1 == 2) {
eff_diff = std::make_shared<MillingtonQuirk1<Traits>>(name);
} else if (eff_diff_type == MillingtonQuirk2<Traits>::type) {
} else if (type_1 == 3) {
eff_diff = std::make_shared<MillingtonQuirk2<Traits>>(name);
} else {
log->error("Effective Diffusion parameterization '{}' "
"has unknown type '{}'",
name, eff_diff_type);
name, type_node["type"].as<std::string>());
DUNE_THROW(IOError, "Unknown Transport parameterization type");
}
// Select type for effective hydromechanic dispersion
const YAML::Node& eff_hm_disp_node = type_node["eff_hydromechanic_disp"];
auto eff_hm_disp_type = eff_hm_disp_node["type"].as<std::string>();
std::shared_ptr<EffectiveHydromechanicDispersion<Traits>> eff_hm_disp;
// Select effective hydromechanic dispersion according to digit 2 of 'type'
uint type_2 = get_digit(type,2);
// Build parametrization for effective hydromechanic dispersion
if (eff_hm_disp_type == ConstEffectiveHydromechanicDispersion<Traits>::type) {
if (type_2 == 1) {
eff_hm_disp = std::make_shared<ConstEffectiveHydromechanicDispersion<Traits>>(name);
} else if (eff_hm_disp_type == IsotropicEffectiveHydromechanicDispersion<Traits>::type) {
} else if (type_2 == 2) {
// Select type for longitudinal and transverse dispersivity
const YAML::Node& long_disp_node = eff_hm_disp_node["long_disp"];
const YAML::Node& trans_disp_node = eff_hm_disp_node["trans_disp"];
auto long_disp_type = long_disp_node["type"].as<std::string>();
auto trans_disp_type = trans_disp_node["type"].as<std::string>();
std::shared_ptr<ConstLongitudinalDispersivity<Traits>> lambda_l;
std::shared_ptr<ConstTransverseDispersivity<Traits>> lambda_t;
// Build parametrization for longitudinal dispersivity
if (long_disp_type == ConstLongitudinalDispersivity<Traits>::type) {
// if (...) {
lambda_l = std::make_shared<ConstLongitudinalDispersivity<Traits>>(name);
} else {
log->error("Longitudinal Dispersivity parameterization '{}' "
"has unknown type '{}'",
name, long_disp_type);
DUNE_THROW(IOError, "Unknown Transport parameterization type");
}
// } else {
// log->error("Longitudinal Dispersivity parameterization '{}' "
// "has unknown type '{}'",
// name, long_disp_type);
// DUNE_THROW(IOError, "Unknown Transport parameterization type");
// }
// Build parametrization for transverse dispersivity
if (trans_disp_type == ConstTransverseDispersivity<Traits>::type) {
// if (...) {
lambda_t = std::make_shared<ConstTransverseDispersivity<Traits>>(name);
} else {
log->error("Transverse Dispersivity parameterization '{}' "
"has unknown type '{}'",
name, trans_disp_type);
DUNE_THROW(IOError, "Unknown Transport parameterization type");
}
// } else {
// log->error("Transverse Dispersivity parameterization '{}' "
// "has unknown type '{}'",
// name, trans_disp_type);
// DUNE_THROW(IOError, "Unknown Transport parameterization type");
// }
eff_hm_disp = std::make_shared<IsotropicEffectiveHydromechanicDispersion<Traits>>(name,lambda_l,lambda_t);
} else {
log->error("Effective Hydromechanic Dispersion parameterization '{}' "
"has unknown type '{}'",
name, eff_hm_disp_type);
name, type_node["type"].as<std::string>());
DUNE_THROW(IOError, "Unknown Transport parameterization type");
}
hd_dips = std::make_shared<HydrodynamicDispersionSuperposition<Traits>>(name,eff_diff,eff_hm_disp);
} else {
log->error("Transport parameterization '{}' has unknown type '{}'",
name, type);
name, type_node["type"].as<std::string>());
DUNE_THROW(IOError, "Unknown Transport parameterization type");
}
return hd_dips;
......
......@@ -41,7 +41,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "const";
static inline std::string type = "Dhd_const";
ConstHydrodynamicDispersionType _const_hydrodynamic_dispersion;
......
......@@ -42,7 +42,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "const";
static inline std::string type = "Deff_const";
ConstEffectiveDiffusionType _const_eff_diff;
......
......@@ -53,7 +53,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "milligton_quirk_1";
static inline std::string type = "Deff_MQ1";
ModelcularDiffusionType _mol_diff;
PorosityType _porosity;
......
......@@ -53,7 +53,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "milligton_quirk_2";
static inline std::string type = "Deff_MQ2";
ModelcularDiffusionType _mol_diff;
PorosityType _porosity;
......
......@@ -41,7 +41,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "const";
static inline std::string type = "Dhm_const";
ConstEffectiveHydromechanicDispersionType _const_eff_hydromechanic_dispersion;
......
......@@ -58,7 +58,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "isotropic";
static inline std::string type = "Dhm_iso";
std::shared_ptr<LambdaL> _long_disp;
std::shared_ptr<LambdaT> _trans_disp;
......
......@@ -55,7 +55,7 @@ public:
};
/// The name of this parameterization type (This is not the instance name)
static inline std::string type = "power_law";
static inline std::string type = "Dhd_pl";
GammaType _gamma;
AlphaType _alpha;
......
......@@ -2,7 +2,7 @@ volumes:
material_0:
index: 0
transport:
type: const
type: Dhd_const
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......@@ -12,7 +12,7 @@ volumes:
material_1:
index: 1
transport:
type: power_law
type: Dhd_pl
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......@@ -22,11 +22,7 @@ volumes:
material_2:
index: 2
transport:
type: superposition
eff_diff:
type: const
eff_hydromechanic_disp:
type: const
type: Deff_const + Dhm_const
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......@@ -37,15 +33,7 @@ volumes:
material_3:
index: 3
transport:
type: superposition
eff_diff:
type: milligton_quirk_1
eff_hydromechanic_disp:
type: isotropic
long_disp:
type: const
trans_disp:
type: const
type: Deff_MQ1 + Dhm_iso
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......@@ -56,11 +44,7 @@ volumes:
material_4:
index: 4
transport:
type: superposition
eff_diff:
type: milligton_quirk_2
eff_hydromechanic_disp:
type: const
type: Deff_MQ2 + Dhm_const
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......
......@@ -2,15 +2,7 @@ volumes:
sand:
index: 0
transport:
type: superposition
eff_diff:
type: milligton_quirk_1
eff_hydromechanic_disp:
type: isotropic
long_disp:
type: const
trans_disp:
type: const
type: Deff_MQ1 + Dhm_iso
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......@@ -21,15 +13,7 @@ volumes:
silt:
index: 1
transport:
type: superposition
eff_diff:
type: milligton_quirk_1
eff_hydromechanic_disp:
type: isotropic
long_disp:
type: const
trans_disp:
type: const
type: Deff_MQ1 + Dhm_iso
parameters:
char_length: 1.5E-11
mol_diff: 2.E-9
......
Markdown is supported
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