CubicEOS.hpp

// Reaktoro is a unified framework for modeling chemically reactive systems.
//
// Copyright © 2014-2024 Allan Leal
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public License
// along with this library. If not, see <http://www.gnu.org/licenses/>.
 
#pragma once
 
// Reaktoro includes
#include <Reaktoro/Common/Constants.hpp>
#include <Reaktoro/Common/Matrix.hpp>
#include <Reaktoro/Common/Types.hpp>
#include <Reaktoro/Core/Model.hpp>
#include <Reaktoro/Core/StateOfMatter.hpp>
 
namespace Reaktoro {
namespace CubicEOS {
 
struct Substance
{
    String formula;    
    real Tcr = NaN;   
    real Pcr = NaN;   
    real omega = NaN; 
};
 
struct Props
{
    real V;            
    real VT;           
    real VP;           
    real Gres;         
    real Hres;         
    real Cpres;        
    real Cvres;        
    ArrayXr Vi;        
    ArrayXr ln_phi;    
    StateOfMatter som; 
};
 
struct Bip
{
    MatrixXr k;   
    MatrixXr kT;  
    MatrixXr kTT; 
};
 
// Note: The list of arguments below for the computation of @eq{k_{ij}} and its
// first and second order temperature derivatives is inspired by the binary
// interaction parameter model of Jaubert et al. (2004):
// \eqc{k_{ij}=\dfrac{-\dfrac{1}{2}\sum_{k=1}^{N_{g}}\sum_{l=1}^{N_{g}}(\alpha_{ik}-\alpha_{jk})(\alpha_{il}-\alpha_{jl})A_{kl}\cdot(298.15T^{-1})^{\left(\frac{B_{kl}}{A_{kl}}-1\right)}-\left(\dfrac{\sqrt{a_{i}}}{b_{i}}-\dfrac{\sqrt{a_{j}}}{b_{j}}\right)}{2\dfrac{\sqrt{a_{i}a_{j}}}{b_{i}b_{j}}}.}
 
struct BipModelArgs
{
    Strings const& substances; 
    real const& T;             
    ArrayXrConstRef Tcr;       
    ArrayXrConstRef Pcr;       
    ArrayXrConstRef omega;     
    ArrayXrConstRef a;         
    ArrayXrConstRef aT;        
    ArrayXrConstRef aTT;       
    ArrayXrConstRef alpha;     
    ArrayXrConstRef alphaT;    
    ArrayXrConstRef alphaTT;   
    ArrayXrConstRef b;         
};
 
using BipModel = Model<Bip(BipModelArgs const&)>;
 
struct Alpha
{
    real alpha;   
    real alphaT;  
    real alphaTT; 
};
 
struct AlphaModelArgs
{
    real const& Tr;    
    real const& TrT;   
    real const& omega; 
};
 
using AlphaModel = Model<Alpha(AlphaModelArgs const&)>;
 
struct EquationModel
{
    real epsilon;      
    real sigma;        
    real Omega;        
    real Psi;          
    AlphaModel alphafn; 
};
 
auto EquationModelVanDerWaals() -> EquationModel;
 
auto EquationModelRedlichKwong() -> EquationModel;
 
auto EquationModelSoaveRedlichKwong() -> EquationModel;
 
auto EquationModelPengRobinson76() -> EquationModel;
 
auto EquationModelPengRobinson78() -> EquationModel;
 
auto EquationModelPengRobinson() -> EquationModel;
 
struct EquationSpecs
{
    Vec<Substance> substances; 
    EquationModel eqmodel;     
    BipModel bipmodel;         
};
 
class Equation
{
public:
    explicit Equation(EquationSpecs const& specs);
 
    Equation(Equation const& other);
 
    virtual ~Equation();
 
    auto operator=(Equation other) -> Equation&;
 
    auto equationSpecs() const -> EquationSpecs const&;
 
    auto compute(Props& props, real const& T, real const& P, ArrayXrConstRef const& x) -> void;
 
private:
    struct Impl;
 
    Ptr<Impl> pimpl;
};
 
struct BipModelParamsPhreeqc
{
    real kH2O_CO2  = 0.19; 
    real kH2O_H2S  = 0.19; 
    real kH2O_CH4  = 0.49; 
    real kH2O_N2   = 0.49; 
    real kH2O_C2H6 = 0.49; 
    real kH2O_C3H8 = 0.55; 
};
 
struct BipModelParamsSoreideWhitson
{
    real kH2O_CO2 = 0.1896;
 
    real kH2O_N2 = 0.4778;
 
    real kH2O_CH4 = 0.4850;
 
    real kH2O_C2H6 = 0.4920;
 
    real kH2O_C3H8 = 0.5525;
 
    real kH2O_nC4H10 = 0.5525;
 
    real kH2O_H2S_a1 = 0.19031;
 
    real kH2O_H2S_a2 = -0.05965;
};
 
auto BipModelPhreeqc(Strings const& substances, BipModelParamsPhreeqc const& params = {}) -> BipModel;
 
auto BipModelSoreideWhitson(Strings const& substances, BipModelParamsSoreideWhitson const& params = {}) -> BipModel;
 
// DEPRECATED METHODS TO BE REMOVED IN THE NEAR FUTURE
 
struct [[deprecated("BipModelParamsPHREEQC has been deprecated and renamed to BipModelParamsPhreeqc.")]] BipModelParamsPHREEQC
{};
 
[[deprecated("BipModelPHREEQC has been deprecated and renamed to BipModelPhreeqc.")]]
auto BipModelPHREEQC(Strings const& substances, BipModelParamsPHREEQC const& params = {}) -> BipModel;
 
} // namespace CubicEOS
} // namespace Reaktoro