diff --git a/src/tests/catch_test_alphaig.cxx b/src/tests/catch_test_alphaig.cxx index c19d24e749c521546bd06c39cb0b819c59b4521c..b541e8810b364faadfc66d637865f312ef4133f3 100644 --- a/src/tests/catch_test_alphaig.cxx +++ b/src/tests/catch_test_alphaig.cxx @@ -9,55 +9,55 @@ using Catch::Approx; using namespace teqp; TEST_CASE("Simplest case","[alphaig]") { - double a_1 = 1, a_2 = 2, T = 300; - nlohmann::json j0 = nlohmann::json::array(); - j0.push_back({ {"type", "Lead"}, { "a_1", 1 }, { "a_2", 2 } }); - nlohmann::json j = nlohmann::json::array(); - j.push_back(j0); - IdealHelmholtz ih(j); - std::valarray<double> molefrac{1.0}; - REQUIRE(ih.alphaig(T, 1, molefrac) == log(1) + a_1 + a_2 / T); + double a_1 = 1, a_2 = 2, T = 300; + nlohmann::json j0 = nlohmann::json::array(); + j0.push_back({ {"type", "Lead"}, { "a_1", 1 }, { "a_2", 2 } }); + nlohmann::json j = nlohmann::json::array(); + j.push_back(j0); + IdealHelmholtz ih(j); + std::valarray<double> molefrac{1.0}; + REQUIRE(ih.alphaig(T, 1, molefrac) == log(1) + a_1 + a_2 / T); } TEST_CASE("alphaig derivative", "[alphaig]") { - double a_1 = 1, a_2 = 2, T = 300, rho = 1; - nlohmann::json j0 = nlohmann::json::array(); - j0.push_back({ {"type", "Lead"}, { "a_1", 1 }, { "a_2", 2 } }); // For the first component - nlohmann::json j = nlohmann::json::array(); - j.push_back(j0); - IdealHelmholtz ih(j); - auto molefrac = (Eigen::ArrayXd(1) << 1.0).finished(); - auto wih = AlphaCallWrapper<1, decltype(ih)>(ih); - wih.alpha(T, rho, molefrac); - using tdx = TDXDerivatives<decltype(ih), double, Eigen::ArrayXd>; - SECTION("All cross derivatives should be zero") { - REQUIRE(tdx::get_Agenxy<1, 1, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); - REQUIRE(tdx::get_Agenxy<1, 2, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); - REQUIRE(tdx::get_Agenxy<1, 3, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); - REQUIRE(tdx::get_Agenxy<2, 1, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); - REQUIRE(tdx::get_Agenxy<2, 2, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); - REQUIRE(tdx::get_Agenxy<2, 3, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); - } + double a_1 = 1, a_2 = 2, T = 300, rho = 1; + nlohmann::json j0 = nlohmann::json::array(); + j0.push_back({ {"type", "Lead"}, { "a_1", 1 }, { "a_2", 2 } }); // For the first component + nlohmann::json j = nlohmann::json::array(); + j.push_back(j0); + IdealHelmholtz ih(j); + auto molefrac = (Eigen::ArrayXd(1) << 1.0).finished(); + auto wih = AlphaCallWrapper<1, decltype(ih)>(ih); + wih.alpha(T, rho, molefrac); + using tdx = TDXDerivatives<decltype(ih), double, Eigen::ArrayXd>; + SECTION("All cross derivatives should be zero") { + REQUIRE(tdx::get_Agenxy<1, 1, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); + REQUIRE(tdx::get_Agenxy<1, 2, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); + REQUIRE(tdx::get_Agenxy<1, 3, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); + REQUIRE(tdx::get_Agenxy<2, 1, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); + REQUIRE(tdx::get_Agenxy<2, 2, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); + REQUIRE(tdx::get_Agenxy<2, 3, ADBackends::autodiff>(wih, T, rho, molefrac) == 0); + } } TEST_CASE("Ammonia derivative", "[alphaig][NH3]") { - double T = 300, rho = 10; - double c0 = 4; - double a1 = -6.59406093943886, a2 = 5.60101151987913; - double Tcrit = 405.56, rhocrit = 13696.0; - std::valarray<double> n = { 2.224, 3.148, 0.9579 }, theta = { 1646, 3965, 7231 }; + double T = 300, rho = 10; + double c0 = 4; + double a1 = -6.59406093943886, a2 = 5.60101151987913; + double Tcrit = 405.56, rhocrit = 13696.0; + std::valarray<double> n = { 2.224, 3.148, 0.9579 }, theta = { 1646, 3965, 7231 }; - using o = nlohmann::json::object_t; - nlohmann::json j = { { - o{ {"type", "Lead"}, { "a_1", a1 - log(rhocrit) }, { "a_2", a2 * Tcrit } }, - o{ {"type", "LogT"}, { "a", -(c0 - 1) } }, - o{ {"type", "Constant"}, { "a", (c0 - 1) * log(Tcrit) } }, // Term from ln(tau) - o{ {"type", "PlanckEinstein"}, { "n", n}, {"theta", theta}} - } }; - IdealHelmholtz ih(j); - auto molefrac = (Eigen::ArrayXd(1) << 1.0).finished(); - auto wih = AlphaCallWrapper<1, decltype(ih)>(ih); - auto calc = wih.alpha(T, rho, molefrac); - auto expected = -5.3492909452728545; - REQUIRE(calc == Approx(expected)); + using o = nlohmann::json::object_t; + nlohmann::json j = { { + o{ {"type", "Lead"}, { "a_1", a1 - log(rhocrit) }, { "a_2", a2 * Tcrit } }, + o{ {"type", "LogT"}, { "a", -(c0 - 1) } }, + o{ {"type", "Constant"}, { "a", (c0 - 1) * log(Tcrit) } }, // Term from ln(tau) + o{ {"type", "PlanckEinstein"}, { "n", n}, {"theta", theta}} + } }; + IdealHelmholtz ih(j); + auto molefrac = (Eigen::ArrayXd(1) << 1.0).finished(); + auto wih = AlphaCallWrapper<1, decltype(ih)>(ih); + auto calc = wih.alpha(T, rho, molefrac); + auto expected = -5.3492909452728545; + REQUIRE(calc == Approx(expected)); } \ No newline at end of file diff --git a/src/tests/catch_test_ammonia_water.cxx b/src/tests/catch_test_ammonia_water.cxx index 0f17cbc0b4c97c336b2337b7914cdb59a66e5d4f..26a35a5155dd68215d014fc890863b468fbd91a1 100644 --- a/src/tests/catch_test_ammonia_water.cxx +++ b/src/tests/catch_test_ammonia_water.cxx @@ -10,58 +10,58 @@ using Catch::Approx; using namespace teqp; TEST_CASE("Trace critical curve w/ Tillner-Roth", "[NHfty3H2O]") { - auto model = AmmoniaWaterTillnerRoth(); - auto z = (Eigen::ArrayXd(2) << 0.7, 0.3).finished(); - auto Ar01 = teqp::TDXDerivatives<decltype(model)>::get_Ar01(model, 300, 300, z); + auto model = AmmoniaWaterTillnerRoth(); + auto z = (Eigen::ArrayXd(2) << 0.7, 0.3).finished(); + auto Ar01 = teqp::TDXDerivatives<decltype(model)>::get_Ar01(model, 300, 300, z); - double T0 = 405.40; - auto rhovec0 = (Eigen::ArrayXd(2) << 225/0.01703026, 0).finished(); - auto prc0 = IsochoricDerivatives<decltype(model)>::get_pr(model, T0, rhovec0); - auto pig0 = rhovec0.sum() * model.R(rhovec0/rhovec0.sum())*T0; - REQUIRE(prc0 + pig0 == Approx(11.33e6).margin(0.01e6)); + double T0 = 405.40; + auto rhovec0 = (Eigen::ArrayXd(2) << 225/0.01703026, 0).finished(); + auto prc0 = IsochoricDerivatives<decltype(model)>::get_pr(model, T0, rhovec0); + auto pig0 = rhovec0.sum() * model.R(rhovec0/rhovec0.sum())*T0; + REQUIRE(prc0 + pig0 == Approx(11.33e6).margin(0.01e6)); - TCABOptions opt; opt.polish = true; opt.integration_order = 1; opt.init_dt = 100; - CriticalTracing<decltype(model)>::trace_critical_arclength_binary(model, T0, rhovec0, "TillnerRoth_crit.csv", opt); + TCABOptions opt; opt.polish = true; opt.integration_order = 1; opt.init_dt = 100; + CriticalTracing<decltype(model)>::trace_critical_arclength_binary(model, T0, rhovec0, "TillnerRoth_crit.csv", opt); } TEST_CASE("Bell et al. REFPROP 10", "[NH3H2O]") { - auto model = build_multifluid_model({ "AMMONIA", "WATER" }, "../mycp", "", {{ "estimate","Lorentz-Berthelot" }}); + auto model = build_multifluid_model({ "AMMONIA", "WATER" }, "../mycp", "", {{ "estimate","Lorentz-Berthelot" }}); - std::string s = R"({ + std::string s = R"({ "0": { "1": { "BIP": { "betaT": 0.933585, - "gammaT": 1.015826, - "betaV": 1.044759, - "gammaV": 1.189754, + "gammaT": 1.015826, + "betaV": 1.044759, + "gammaV": 1.189754, "Fij": 1.0 }, "departure":{ - "type": "Gaussian+Exponential", - "n": [-2.00211,3.0813,-1.75352,2.9816,-3.82588,-1.7385,0.42008], - "t": [0.25,2.0,0.5,2.0,1.0,4.0,1.0], - "d": [1.0,1.0,1.0,1.0,1.0,1.0,3.0], - "l": [2.0,1.0,0.0,0.0,0.0,0.0,0.0], - "eta": [0.0,0.0,0.0,0.746,4.25,0.7,3.0], - "beta": [0.0,0.0,0.27,0.86,3.0,0.5,4.0], - "gamma": [0.0,0.0,2.8,1.8,1.5,0.8,1.3], - "epsilon": [0.0,0.0,0.0,2.0,-0.25,1.85,0.3], - "Npower": 2 + "type": "Gaussian+Exponential", + "n": [-2.00211,3.0813,-1.75352,2.9816,-3.82588,-1.7385,0.42008], + "t": [0.25,2.0,0.5,2.0,1.0,4.0,1.0], + "d": [1.0,1.0,1.0,1.0,1.0,1.0,3.0], + "l": [2.0,1.0,0.0,0.0,0.0,0.0,0.0], + "eta": [0.0,0.0,0.0,0.746,4.25,0.7,3.0], + "beta": [0.0,0.0,0.27,0.86,3.0,0.5,4.0], + "gamma": [0.0,0.0,2.8,1.8,1.5,0.8,1.3], + "epsilon": [0.0,0.0,0.0,2.0,-0.25,1.85,0.3], + "Npower": 2 } } } })"; - auto mutant = build_multifluid_mutant(model, nlohmann::json::parse(s)); - - double T0 = model.redfunc.Tc[0]; - auto rhovec0 = (Eigen::ArrayXd(2) << 1/model.redfunc.vc[0], 0).finished(); - auto der0 = CriticalTracing<decltype(mutant)>::get_drhovec_dT_crit(mutant, T0, rhovec0); - - auto prc0 = IsochoricDerivatives<decltype(mutant)>::get_pr(mutant, T0, rhovec0); - auto pig0 = rhovec0.sum() * mutant.R(rhovec0 / rhovec0.sum()) * T0; + auto mutant = build_multifluid_mutant(model, nlohmann::json::parse(s)); + + double T0 = model.redfunc.Tc[0]; + auto rhovec0 = (Eigen::ArrayXd(2) << 1/model.redfunc.vc[0], 0).finished(); + auto der0 = CriticalTracing<decltype(mutant)>::get_drhovec_dT_crit(mutant, T0, rhovec0); + + auto prc0 = IsochoricDerivatives<decltype(mutant)>::get_pr(mutant, T0, rhovec0); + auto pig0 = rhovec0.sum() * mutant.R(rhovec0 / rhovec0.sum()) * T0; - TCABOptions opt; opt.polish = true; opt.integration_order = 1; opt.init_dt = 100; opt.verbosity = 1000; + TCABOptions opt; opt.polish = true; opt.integration_order = 1; opt.init_dt = 100; opt.verbosity = 1000; - CriticalTracing<decltype(mutant)>::trace_critical_arclength_binary(mutant, T0, rhovec0, "BellREFPROP10_NH3.csv", opt); + CriticalTracing<decltype(mutant)>::trace_critical_arclength_binary(mutant, T0, rhovec0, "BellREFPROP10_NH3.csv", opt); } \ No newline at end of file diff --git a/src/tests/catch_test_mutant.cxx b/src/tests/catch_test_mutant.cxx index 6f4498eefb336aaf5d4ed83d98bc5901dd12daf4..bc8faf8c8eacf5ef1241ec7cf28b7be37ec11fdc 100644 --- a/src/tests/catch_test_mutant.cxx +++ b/src/tests/catch_test_mutant.cxx @@ -12,10 +12,10 @@ using namespace teqp; TEST_CASE("Test construction of mutant", "[mutant]") { - std::string coolprop_root = "../mycp"; - auto BIPcollection = coolprop_root + "/dev/mixtures/mixture_binary_pairs.json"; + std::string coolprop_root = "../mycp"; + auto BIPcollection = coolprop_root + "/dev/mixtures/mixture_binary_pairs.json"; - auto model = build_multifluid_model({ "Nitrogen", "Ethane" }, coolprop_root, BIPcollection); + auto model = build_multifluid_model({ "Nitrogen", "Ethane" }, coolprop_root, BIPcollection); std::string s0 = R"({"0": {} })"; nlohmann::json j0 = nlohmann::json::parse(s0); @@ -36,7 +36,7 @@ TEST_CASE("Test construction of mutant", "[mutant]") } } })"; - nlohmann::json j = nlohmann::json::parse(s); + nlohmann::json j = nlohmann::json::parse(s); auto mutant = build_multifluid_mutant(model, j); double T = 300, rho = 300;