Page MenuHomec4science

laws.cpp
No OneTemporary

File Metadata

Created
Thu, May 23, 17:30

laws.cpp

#include "catch.hpp"
#include "specmicp_common/physics/laws.hpp"
#include "specmicp_common/physics/unsaturated_laws.hpp"
#include <iostream>
using namespace specmicp;
TEST_CASE("Gas perfect") {
SECTION("Pressure") {
REQUIRE(laws::pressure_perfect_gas(1.0, 1.0, 1.0) == Approx(constants::gas_constant));
REQUIRE(laws::pressure_perfect_gas(2.0, 1.0, 1.0) == Approx(2.0*constants::gas_constant));
REQUIRE(laws::pressure_perfect_gas(1.0, 3.0, 1.0) == Approx(constants::gas_constant/3.0));
REQUIRE(laws::pressure_perfect_gas(1.0, 1.0, 4.0) == Approx(constants::gas_constant*4.0));
REQUIRE(laws::pressure_perfect_gas(2.0, 3.0, 4.0) == Approx(2.0*constants::gas_constant*4.0/3.0));
}
SECTION("Mole") {
REQUIRE(laws::mole_perfect_gas(1.0, 1.0, 1.0) == Approx(1.0/constants::gas_constant));
REQUIRE(laws::mole_perfect_gas(2.0, 1.0, 1.0) == Approx(2.0/constants::gas_constant));
REQUIRE(laws::mole_perfect_gas(1.0, 3.0, 1.0) == Approx(3.0/constants::gas_constant));
REQUIRE(laws::mole_perfect_gas(1.0, 1.0, 4.0) == Approx(1.0/(4.0*constants::gas_constant)));
REQUIRE(laws::mole_perfect_gas(2.0, 3.0, 4.0) == Approx(2.0*3.0/(4.0*constants::gas_constant)));
}
}
TEST_CASE("Debye-Huckel") {
SECTION("Debye-Huckel") {
REQUIRE(laws::debye_huckel(1.0, 1.0, 1.0) == Approx(-(constants::Adebye)/(1+constants::Bdebye)));
REQUIRE(laws::debye_huckel(1.0, 2.0, 1.0) == Approx(-(constants::Adebye*4)/(1+constants::Bdebye)));
REQUIRE(laws::debye_huckel(0.1, 1.0, 1.0) == Approx(-(constants::Adebye*0.1)/(1+constants::Bdebye*0.1)));
REQUIRE(laws::debye_huckel(0.1, 2.0, 1.0) == Approx(-(constants::Adebye*0.1*4)/(1+constants::Bdebye*0.1)));
REQUIRE(laws::debye_huckel(1.0, 1.0, 2.0) == Approx(-(constants::Adebye)/(1+constants::Bdebye*1.0*2.0)));
REQUIRE(laws::debye_huckel(0.1, 2.0, 3.0) == Approx(-(constants::Adebye*0.1*4)/(1+constants::Bdebye*0.1*3.0)));
}
SECTION("Extended Debye-Huckel") {
REQUIRE(laws::extended_debye_huckel(1.0, 1.0, 1.0, 1.0, 1.0) ==
Approx(-(constants::Adebye)/(1.0+constants::Bdebye)+1.0*1.0));
REQUIRE(laws::extended_debye_huckel(1.0, 1.0, 2.0, 1.0, 1.0) ==
Approx(-(constants::Adebye*4.0)/(1.0+constants::Bdebye)+1.0*1.0));
REQUIRE(laws::extended_debye_huckel(1.0, 1.0, 1.0, 1.0, 2.0) ==
Approx(-(constants::Adebye*1.0)/(1.0+constants::Bdebye)+1.0*2.0));
REQUIRE(laws::extended_debye_huckel(0.1, 0.01, 1.0, 1.0, 2.0) ==
Approx(-(constants::Adebye*1.0*0.01)/(1.0+constants::Bdebye*0.01)+2.0*0.1));
REQUIRE(laws::extended_debye_huckel(0.01, 0.1, -2.0, 3.0, 2.0) ==
Approx(-(4*constants::Adebye*1.0*0.1)/(1.0+constants::Bdebye*0.1*3.0)+2.0*0.01));
REQUIRE(laws::extended_debye_huckel(0.01, -2.0, 3.0, 2.0) ==
Approx(-(4*constants::Adebye*1.0*0.1)/(1.0+constants::Bdebye*0.1*3.0)+2.0*0.01));
}
}
TEST_CASE("Unsaturated laws", "[laws],[unsaturated],[physics]") {
SECTION("make_saturation_model") {
auto identity = [](scalar_t saturation) -> scalar_t {return saturation;};
auto model = laws::make_saturation_model(identity);
CHECK(model(1.0) == 1.0);
CHECK(model(0.5) == 0.5);
auto mult = [](scalar_t sat, scalar_t mult) -> scalar_t {
return mult*sat;
};
auto model2 = laws::make_saturation_model(mult, 2.0);
CHECK(model2(1.0) == Approx(2.0));
CHECK(model2(0.5) == Approx(1.0));
}
SECTION("Baroghel Bouny model") {
auto cap = laws::make_saturation_model(
&laws::capillary_pressure_BaroghelBouny,
1.0e8, 2.0);
CHECK(cap(1.0) == 0.0);
CHECK(cap(0.5) == Approx(1.732050808e8));
auto rh = laws::make_rh_model(cap);
CHECK(rh(1.0) == 1.0);
CHECK(rh(0.0) == 0.0);
CHECK(rh(0.5) == Approx(0.2832770066).epsilon(1e-3));
cap = laws::make_saturation_model(
&laws::capillary_pressure_BaroghelBouny,
2.0, 2.0);
CHECK(cap(0.5) == Approx(2.0*1.732050808));
}
SECTION("Mualem model") {
auto krl = laws::make_saturation_model(
&laws::relative_liquid_permeability_Mualem,
2.0
);
CHECK(krl(1.0) == Approx(1.0));
CHECK(krl(0.0) == Approx(0.0));
CHECK(krl(0.5) == Approx(0.5909902577));
auto krg = laws::make_saturation_model(
&laws::relative_gas_permeability_Mualem,
2.0
);
CHECK(krg(1.0) == Approx(0.0));
CHECK(krg(0.0) == Approx(1.0));
CHECK(krg(0.5) == Approx(0.005203820043));
}
}

Event Timeline