From 5b25edc7755e4814590788ee68bacffb155b3945 Mon Sep 17 00:00:00 2001
From: Ian Bell <ian.bell@nist.gov>
Date: Thu, 8 Apr 2021 09:01:31 -0400
Subject: [PATCH] Allow PCSAFT models to be specified by constants directly

---
 include/teqp/models/pcsaft.hpp | 27 ++++++++++++++++++++++++---
 interface/pybind11_wrapper.cpp | 16 ++++++++++++++++
 2 files changed, 40 insertions(+), 3 deletions(-)

diff --git a/include/teqp/models/pcsaft.hpp b/include/teqp/models/pcsaft.hpp
index b55178e..9b39958 100644
--- a/include/teqp/models/pcsaft.hpp
+++ b/include/teqp/models/pcsaft.hpp
@@ -176,7 +176,7 @@ private:
     std::vector<std::string> names;
     double k_ij; ///< binary interaction parameter
 public:
-    PCSAFTMixture(const std::vector<std::string> names) : names(names)
+    PCSAFTMixture(const std::vector<std::string> &names) : names(names)
     {
         m.resize(names.size());
         mminus1.resize(names.size());
@@ -191,9 +191,30 @@ public:
             epsilon_over_k[i] = coeff.epsilon_over_k;
             i++;
         }
-
         k_ij = 0;
     };
+    PCSAFTMixture(const std::vector<SAFTCoeffs> &coeffs) 
+    {
+        m.resize(coeffs.size());
+        mminus1.resize(coeffs.size());
+        sigma_Angstrom.resize(coeffs.size());
+        epsilon_over_k.resize(coeffs.size());
+        names.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;
+            i++;
+        }
+        k_ij = 0;
+    };
+    auto get_m(){ return m; }
+    auto get_sigma_Angstrom() { return sigma_Angstrom; }
+    auto get_epsilon_over_k_K() { return epsilon_over_k; }
+
     void print_info() {
         std::cout << "i m sigma / A e/kB / K \n  ++++++++++++++" << std::endl;
         for (auto i = 0; i < m.size(); ++i) {
@@ -207,7 +228,7 @@ public:
         for (auto i = 0; i < N; ++i) {
             d[i] = sigma_Angstrom[i] * (1.0 - 0.12 * exp(-3.0 * epsilon_over_k[i] / T));
         }
-        return 6 * 0.74 / EIGEN_PI / (mole_fractions*m*powvec(d, 3)).sum();
+        return 6 * 0.74 / EIGEN_PI / (mole_fractions*m*powvec(d, 3)).sum()*1e30; // particles/m^3
     }
     const double R = get_R_gas<double>();
 
diff --git a/interface/pybind11_wrapper.cpp b/interface/pybind11_wrapper.cpp
index 2f7da8c..c52b413 100644
--- a/interface/pybind11_wrapper.cpp
+++ b/interface/pybind11_wrapper.cpp
@@ -63,8 +63,24 @@ void init_teqp(py::module& m) {
         ;
     add_derivatives<vdWEOS1>(m);
 
+    py::class_<SAFTCoeffs>(m, "SAFTCoeffs")
+        .def(py::init<>())
+        .def_readwrite("name", &SAFTCoeffs::name)
+        .def_readwrite("m", &SAFTCoeffs::m)
+        .def_readwrite("sigma_Angstrom", &SAFTCoeffs::sigma_Angstrom)
+        .def_readwrite("epsilon_over_k", &SAFTCoeffs::epsilon_over_k)
+        .def_readwrite("BibTeXKey", &SAFTCoeffs::BibTeXKey)
+        ;
+
     py::class_<PCSAFTMixture>(m, "PCSAFTEOS")
         .def(py::init<const std::vector<std::string> &>(), py::arg("names"))
+        .def(py::init<const std::vector<SAFTCoeffs> &>(), py::arg("coeffs"))
+        .def("print_info", &PCSAFTMixture::print_info)
+        .def("max_rhoN", &PCSAFTMixture::max_rhoN<Eigen::ArrayXd>)
+        .def("get_m", &PCSAFTMixture::get_m)
+        .def("get_sigma_Angstrom", &PCSAFTMixture::get_sigma_Angstrom)
+        .def("get_epsilon_over_k_K", &PCSAFTMixture::get_epsilon_over_k_K)
+
         ;
     add_derivatives<PCSAFTMixture>(m);
 
-- 
GitLab