From 7c1bf5328dc0b39f83102184fee7f1f91c378bd9 Mon Sep 17 00:00:00 2001 From: Ian Bell Date: Sat, 23 Sep 2023 13:15:57 -0400 Subject: [PATCH] Implement get_names and get_BibTeXKeys for PC-SAFT and SAFT-VR-Mie Closes #51 --- include/teqp/models/pcsaft.hpp | 5 +++- include/teqp/models/saftvrmie.hpp | 41 ++++++++++++++++++++---------- interface/pybind11_wrapper.cpp | 3 +++ src/tests/catch_test_PCSAFT.cxx | 3 ++- src/tests/catch_test_SAFTVRMie.cxx | 8 ++++++ 5 files changed, 45 insertions(+), 15 deletions(-) diff --git a/include/teqp/models/pcsaft.hpp b/include/teqp/models/pcsaft.hpp index 0fb5d468..0072e998 100644 --- a/include/teqp/models/pcsaft.hpp +++ b/include/teqp/models/pcsaft.hpp @@ -286,7 +286,7 @@ class PCSAFTMixture { mminus1, ///< m-1 sigma_Angstrom, ///< epsilon_over_k; ///< depth of pair potential divided by Boltzman constant - std::vector names; + std::vector names, bibtex; Eigen::ArrayXXd kmat; ///< binary interaction parameter matrix PCSAFTHardChainContribution hardchain; @@ -316,6 +316,7 @@ class PCSAFTMixture { sigma_Angstrom.resize(coeffs.size()); epsilon_over_k.resize(coeffs.size()); names.resize(coeffs.size()); + bibtex.resize(coeffs.size()); auto i = 0; for (const auto &coeff : coeffs) { m[i] = coeff.m; @@ -323,6 +324,7 @@ class PCSAFTMixture { sigma_Angstrom[i] = coeff.sigma_Angstrom; epsilon_over_k[i] = coeff.epsilon_over_k; names[i] = coeff.name; + bibtex[i] = coeff.BibTeXKey; i++; } return PCSAFTHardChainContribution(m, mminus1, sigma_Angstrom, epsilon_over_k, kmat); @@ -374,6 +376,7 @@ class PCSAFTMixture { auto get_epsilon_over_k_K() const { return epsilon_over_k; } auto get_kmat() const { return kmat; } auto get_names() const { return names;} + auto get_BibTeXKeys() const { return bibtex;} auto print_info() { std::string s = std::string("i m sigma / A e/kB / K \n ++++++++++++++") + "\n"; diff --git a/include/teqp/models/saftvrmie.hpp b/include/teqp/models/saftvrmie.hpp index 29e52d80..6c742199 100644 --- a/include/teqp/models/saftvrmie.hpp +++ b/include/teqp/models/saftvrmie.hpp @@ -741,7 +741,7 @@ struct SAFTVRMieChainContributionTerms{ class SAFTVRMieMixture { private: - std::vector names; + std::vector names, bibtex; const SAFTVRMieChainContributionTerms terms; const std::optional polar; // Can be present or not @@ -760,6 +760,20 @@ class SAFTVRMieMixture { SAFTVRMieLibrary library; return library.get_coeffs(names); } + auto get_names(const std::vector &coeffs){ + std::vector names_; + for (auto c : coeffs){ + names_.push_back(c.name); + } + return names_; + } + auto get_bibtex(const std::vector &coeffs){ + std::vector keys_; + for (auto c : coeffs){ + keys_.push_back(c.BibTeXKey); + } + return keys_; + } public: static auto build_chain(const std::vector &coeffs, const std::optional& kmat, const std::optional& flags = std::nullopt){ if (kmat){ @@ -807,8 +821,8 @@ class SAFTVRMieMixture { public: SAFTVRMieMixture(const std::vector &names, const std::optional& kmat = std::nullopt, const std::optional&flags = std::nullopt) : SAFTVRMieMixture(get_coeffs_from_names(names), kmat, flags){}; - SAFTVRMieMixture(const std::vector &coeffs, const std::optional &kmat = std::nullopt, const std::optional&flags = std::nullopt) : terms(build_chain(coeffs, kmat, flags)), polar(build_polar(coeffs)) {}; - SAFTVRMieMixture(SAFTVRMieChainContributionTerms&& terms, std::optional &&polar = std::nullopt) : terms(std::move(terms)), polar(std::move(polar)) {}; + SAFTVRMieMixture(const std::vector &coeffs, const std::optional &kmat = std::nullopt, const std::optional&flags = std::nullopt) : names(get_names(coeffs)), bibtex(get_bibtex(coeffs)), terms(build_chain(coeffs, kmat, flags)), polar(build_polar(coeffs)) {}; + SAFTVRMieMixture(SAFTVRMieChainContributionTerms&& terms, const std::vector &coeffs, std::optional &&polar = std::nullopt) : names(get_names(coeffs)), bibtex(get_bibtex(coeffs)), terms(std::move(terms)), polar(std::move(polar)) {}; // PCSAFTMixture( const PCSAFTMixture& ) = delete; // non construction-copyable @@ -855,7 +869,8 @@ class SAFTVRMieMixture { {"alphar_chain", val.alphar_chain} }; } - + auto get_names() const { return names; } + auto get_BibTeXKeys() const { return bibtex; } auto get_m() const { return terms.m; } auto get_sigma_Angstrom() const { return (terms.sigma_A).eval(); } auto get_sigma_m() const { return terms.sigma_A/1e10; } @@ -965,7 +980,7 @@ inline auto SAFTVRMiefactory(const nlohmann::json & spec){ if (!something_polar){ // Nonpolar, just m, epsilon, sigma and possibly a kmat matrix with kij coefficients - return SAFTVRMieMixture(SAFTVRMieMixture::build_chain(coeffs, kmat)); + return SAFTVRMieMixture(SAFTVRMieMixture::build_chain(coeffs, kmat), coeffs); } else{ // Polar term is also provided, along with the chain terms @@ -1046,34 +1061,34 @@ inline auto SAFTVRMiefactory(const nlohmann::json & spec){ auto mustar2 = (mustar2factor*mu_Cm.pow(2)/(ms*epsks*sigma_ms.pow(3))).eval(); auto Qstar2 = (Qstar2factor*Q_Cm2.pow(2)/(ms*epsks*sigma_ms.pow(5))).eval(); auto polar = MultipolarContributionGrossVrabec(ms, sigma_ms*1e10, epsks, mustar2, nmu, Qstar2, nQ); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } if (polar_model == "GubbinsTwu+Luckas"){ using MCGTL = MultipolarContributionGubbinsTwu; auto mubar2 = (mustar2factor*mu_Cm.pow(2)/(epsks*sigma_ms.pow(3))).eval(); auto Qbar2 = (Qstar2factor*Q_Cm2.pow(2)/(epsks*sigma_ms.pow(5))).eval(); auto polar = MCGTL(sigma_ms, epsks, mubar2, Qbar2, multipolar_rhostar_approach::use_packing_fraction); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } if (polar_model == "GubbinsTwu+GubbinsTwu"){ using MCGG = MultipolarContributionGubbinsTwu; auto mubar2 = (mustar2factor*mu_Cm.pow(2)/(epsks*sigma_ms.pow(3))).eval(); auto Qbar2 = (Qstar2factor*Q_Cm2.pow(2)/(epsks*sigma_ms.pow(5))).eval(); auto polar = MCGG(sigma_ms, epsks, mubar2, Qbar2, multipolar_rhostar_approach::use_packing_fraction); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } if (polar_model == "GubbinsTwu+Gottschalk"){ using MCGG = MultipolarContributionGubbinsTwu; auto mubar2 = (mustar2factor*mu_Cm.pow(2)/(epsks*sigma_ms.pow(3))).eval(); auto Qbar2 = (Qstar2factor*Q_Cm2.pow(2)/(epsks*sigma_ms.pow(5))).eval(); auto polar = MCGG(sigma_ms, epsks, mubar2, Qbar2, multipolar_rhostar_approach::use_packing_fraction); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } if (polar_model == "GrayGubbins+GubbinsTwu"){ using MCGG = MultipolarContributionGrayGubbins; auto polar = MCGG(sigma_ms, epsks, SIGMAIJ, EPSKIJ, mu_Cm, Q_Cm2, polar_flags); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } // if (polar_model == "GrayGubbins+Gottschalk"){ // using MCGG = MultipolarContributionGrayGubbins; @@ -1083,7 +1098,7 @@ inline auto SAFTVRMiefactory(const nlohmann::json & spec){ if (polar_model == "GrayGubbins+Luckas"){ using MCGG = MultipolarContributionGrayGubbins; auto polar = MCGG(sigma_ms, epsks, SIGMAIJ, EPSKIJ, mu_Cm, Q_Cm2, polar_flags); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } @@ -1092,14 +1107,14 @@ inline auto SAFTVRMiefactory(const nlohmann::json & spec){ auto mubar2 = (mustar2factor*mu_Cm.pow(2)/(epsks*sigma_ms.pow(3))).eval(); auto Qbar2 = (Qstar2factor*Q_Cm2.pow(2)/(epsks*sigma_ms.pow(5))).eval(); auto polar = MCGTL(sigma_ms, epsks, mubar2, Qbar2, multipolar_rhostar_approach::calculate_Gubbins_rhostar); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } if (polar_model == "GubbinsTwu+GubbinsTwu+GubbinsTwuRhostar"){ using MCGG = MultipolarContributionGubbinsTwu; auto mubar2 = (mustar2factor*mu_Cm.pow(2)/(epsks*sigma_ms.pow(3))).eval(); auto Qbar2 = (Qstar2factor*Q_Cm2.pow(2)/(epsks*sigma_ms.pow(5))).eval(); auto polar = MCGG(sigma_ms, epsks, mubar2, Qbar2, multipolar_rhostar_approach::calculate_Gubbins_rhostar); - return SAFTVRMieMixture(std::move(chain), std::move(polar)); + return SAFTVRMieMixture(std::move(chain), coeffs, std::move(polar)); } throw teqp::InvalidArgument("didn't understand this polar_model:"+polar_model); } diff --git a/interface/pybind11_wrapper.cpp b/interface/pybind11_wrapper.cpp index dd473563..43fbdfa7 100644 --- a/interface/pybind11_wrapper.cpp +++ b/interface/pybind11_wrapper.cpp @@ -138,10 +138,13 @@ void attach_model_specific_methods(py::object& obj){ setattr("get_sigma_Angstrom", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_sigma_Angstrom(); }), obj)); setattr("get_epsilon_over_k_K", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_epsilon_over_k_K(); }), obj)); setattr("get_names", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_names(); }), obj)); + setattr("get_BibTeXKeys", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_BibTeXKeys(); }), obj)); setattr("max_rhoN", MethodType(py::cpp_function([](py::object& o, double T, REArrayd& molefrac){ return get_typed(o).max_rhoN(T, molefrac); }, "self"_a, "T"_a, "molefrac"_a), obj)); setattr("get_kmat", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_kmat(); }), obj)); } else if (index == SAFTVRMie_i){ + setattr("get_names", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_names(); }), obj)); + setattr("get_BibTeXKeys", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_BibTeXKeys(); }), obj)); setattr("get_m", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_m(); }), obj)); setattr("get_sigma_Angstrom", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_sigma_Angstrom(); }), obj)); setattr("get_sigma_m", MethodType(py::cpp_function([](py::object& o){ return get_typed(o).get_sigma_m(); }), obj)); diff --git a/src/tests/catch_test_PCSAFT.cxx b/src/tests/catch_test_PCSAFT.cxx index e5d9f16a..b5afc50c 100644 --- a/src/tests/catch_test_PCSAFT.cxx +++ b/src/tests/catch_test_PCSAFT.cxx @@ -26,11 +26,12 @@ TEST_CASE("Single alphar check value", "[PCSAFT]") CHECK(tdx::get_Ar00(model, T, Dmolar, z) == Approx(-0.032400020930842724)); } -TEST_CASE("Check get_names", "[PCSAFT]") +TEST_CASE("Check get_names and get_BibTeXKeys", "[PCSAFT]") { std::vector names = { "Methane" }; auto model = PCSAFTMixture(names); CHECK(model.get_names()[0] == "Methane"); + CHECK(model.get_BibTeXKeys()[0] == "Gross-IECR-2001"); } TEST_CASE("Check 0n derivatives", "[PCSAFT]") diff --git a/src/tests/catch_test_SAFTVRMie.cxx b/src/tests/catch_test_SAFTVRMie.cxx index e9bfba7a..fd01a8c4 100644 --- a/src/tests/catch_test_SAFTVRMie.cxx +++ b/src/tests/catch_test_SAFTVRMie.cxx @@ -22,6 +22,14 @@ using namespace teqp; #include "boost/multiprecision/cpp_complex.hpp" #include "teqp/math/quadrature.hpp" +TEST_CASE("Check get_names and get_BibTeXKeys", "[SAFTVRMIE]") +{ + std::vector names = { "Methane", "Ethane" }; + auto model = SAFTVRMieMixture(names); + CHECK(model.get_names()[0] == "Methane"); + CHECK(model.get_BibTeXKeys()[0] == "Lafitte-JCP-2001"); +} + TEST_CASE("Check integration", "[SAFTVRMIE]"){ std::function f = [](const double&x){ return x*sin(x); }; auto exact = -2*cos(1) + 2*sin(1);