Skip to content

Commit

Permalink
Implement get_names and get_BibTeXKeys for PC-SAFT and SAFT-VR-Mie
Browse files Browse the repository at this point in the history 
Closes #51
  • Loading branch information
@ianhbell
ianhbell committed Sep 23, 2023
1 parent 4565948 commit 7c1bf53
Show file tree
Hide file tree
Showing 5 changed files with 45 additions and 15 deletions.
5 changes: 4 additions & 1 deletion include/teqp/models/pcsaft.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -286,7 +286,7 @@ class PCSAFTMixture {
mminus1, ///< m-1
sigma_Angstrom, ///<
epsilon_over_k; ///< depth of pair potential divided by Boltzman constant
std::vector<std::string> names;
std::vector<std::string> names, bibtex;
Eigen::ArrayXXd kmat; ///< binary interaction parameter matrix

PCSAFTHardChainContribution hardchain;
Expand Down Expand Up @@ -316,13 +316,15 @@ 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;
mminus1[i] = m[i] - 1;
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);
Expand Down Expand Up @@ -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";
Expand Down
41 changes: 28 additions & 13 deletions include/teqp/models/saftvrmie.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -741,7 +741,7 @@ struct SAFTVRMieChainContributionTerms{
class SAFTVRMieMixture {
private:

std::vector<std::string> names;
std::vector<std::string> names, bibtex;
const SAFTVRMieChainContributionTerms terms;
const std::optional<SAFTpolar::multipolar_contributions_variant> polar; // Can be present or not

Expand All @@ -760,6 +760,20 @@ class SAFTVRMieMixture {
SAFTVRMieLibrary library;
return library.get_coeffs(names);
}
auto get_names(const std::vector<SAFTVRMieCoeffs> &coeffs){
std::vector<std::string> names_;
for (auto c : coeffs){
names_.push_back(c.name);
}
return names_;
}
auto get_bibtex(const std::vector<SAFTVRMieCoeffs> &coeffs){
std::vector<std::string> keys_;
for (auto c : coeffs){
keys_.push_back(c.BibTeXKey);
}
return keys_;
}
public:
static auto build_chain(const std::vector<SAFTVRMieCoeffs> &coeffs, const std::optional<Eigen::ArrayXXd>& kmat, const std::optional<nlohmann::json>& flags = std::nullopt){
if (kmat){
Expand Down Expand Up @@ -807,8 +821,8 @@ class SAFTVRMieMixture {

public:
SAFTVRMieMixture(const std::vector<std::string> &names, const std::optional<Eigen::ArrayXXd>& kmat = std::nullopt, const std::optional<nlohmann::json>&flags = std::nullopt) : SAFTVRMieMixture(get_coeffs_from_names(names), kmat, flags){};
SAFTVRMieMixture(const std::vector<SAFTVRMieCoeffs> &coeffs, const std::optional<Eigen::ArrayXXd> &kmat = std::nullopt, const std::optional<nlohmann::json>&flags = std::nullopt) : terms(build_chain(coeffs, kmat, flags)), polar(build_polar(coeffs)) {};
SAFTVRMieMixture(SAFTVRMieChainContributionTerms&& terms, std::optional<SAFTpolar::multipolar_contributions_variant> &&polar = std::nullopt) : terms(std::move(terms)), polar(std::move(polar)) {};
SAFTVRMieMixture(const std::vector<SAFTVRMieCoeffs> &coeffs, const std::optional<Eigen::ArrayXXd> &kmat = std::nullopt, const std::optional<nlohmann::json>&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<SAFTVRMieCoeffs> &coeffs, std::optional<SAFTpolar::multipolar_contributions_variant> &&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
Expand Down Expand Up @@ -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; }
Expand Down Expand Up @@ -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
Expand Down Expand Up @@ -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<LuckasJIntegral, LuckasKIntegral>;
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<GubbinsTwuJIntegral, GubbinsTwuKIntegral>;
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<GottschalkJIntegral, GottschalkKIntegral>;
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<GubbinsTwuJIntegral, GubbinsTwuKIntegral>;
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<GottschalkJIntegral, GottschalkKIntegral>;
Expand All @@ -1083,7 +1098,7 @@ inline auto SAFTVRMiefactory(const nlohmann::json & spec){
if (polar_model == "GrayGubbins+Luckas"){
using MCGG = MultipolarContributionGrayGubbins<LuckasJIntegral, LuckasKIntegral>;
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));
}


Expand All @@ -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<GubbinsTwuJIntegral, GubbinsTwuKIntegral>;
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);
}
Expand Down
3 changes: 3 additions & 0 deletions interface/pybind11_wrapper.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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<PCSAFT_t>(o).get_sigma_Angstrom(); }), obj));
setattr("get_epsilon_over_k_K", MethodType(py::cpp_function([](py::object& o){ return get_typed<PCSAFT_t>(o).get_epsilon_over_k_K(); }), obj));
setattr("get_names", MethodType(py::cpp_function([](py::object& o){ return get_typed<PCSAFT_t>(o).get_names(); }), obj));
setattr("get_BibTeXKeys", MethodType(py::cpp_function([](py::object& o){ return get_typed<PCSAFT_t>(o).get_BibTeXKeys(); }), obj));
setattr("max_rhoN", MethodType(py::cpp_function([](py::object& o, double T, REArrayd& molefrac){ return get_typed<PCSAFT_t>(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<PCSAFT_t>(o).get_kmat(); }), obj));
}
else if (index == SAFTVRMie_i){
setattr("get_names", MethodType(py::cpp_function([](py::object& o){ return get_typed<SAFTVRMie_t>(o).get_names(); }), obj));
setattr("get_BibTeXKeys", MethodType(py::cpp_function([](py::object& o){ return get_typed<SAFTVRMie_t>(o).get_BibTeXKeys(); }), obj));
setattr("get_m", MethodType(py::cpp_function([](py::object& o){ return get_typed<SAFTVRMie_t>(o).get_m(); }), obj));
setattr("get_sigma_Angstrom", MethodType(py::cpp_function([](py::object& o){ return get_typed<SAFTVRMie_t>(o).get_sigma_Angstrom(); }), obj));
setattr("get_sigma_m", MethodType(py::cpp_function([](py::object& o){ return get_typed<SAFTVRMie_t>(o).get_sigma_m(); }), obj));
Expand Down
3 changes: 2 additions & 1 deletion src/tests/catch_test_PCSAFT.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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<std::string> 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]")
Expand Down
8 changes: 8 additions & 0 deletions src/tests/catch_test_SAFTVRMie.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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<std::string> 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<double(double)> f = [](const double&x){ return x*sin(x); };
auto exact = -2*cos(1) + 2*sin(1);
Expand Down

0 comments on commit 7c1bf53

Please sign in to comment.