6#ifndef CT_PENGROBINSON_H
7#define CT_PENGROBINSON_H
31 explicit PengRobinson(
const string& infile=
"",
const string&
id=
"");
33 string type()
const override {
34 return "Peng-Robinson";
40 double cp_mole()
const override;
41 double cv_mole()
const override;
154 bool addSpecies(shared_ptr<Species> spec)
override;
173 void setBinaryCoeffs(
const string& species_i,
const string& species_j,
double a);
178 double sresid()
const override;
179 double hresid()
const override;
181 double liquidVolEst(
double T,
double& pres)
const override;
186 double dpdVCalc(
double T,
double molarVol,
double& presCalc)
const override;
231 void calculateAB(
double& aCalc,
double& bCalc,
double& aAlpha)
const;
233 void calcCriticalConditions(
double& pc,
double& tc,
double& vc)
const override;
236 int solveCubic(
double T,
double pres,
double a,
double b,
double aAlpha,
237 double Vroot[3])
const;
259 vector<double> m_b_coeffs;
260 vector<double> m_kappa;
262 mutable vector<double> m_dalphadT;
263 mutable vector<double> m_d2alphadT2;
264 vector<double> m_alpha;
276 double m_Vroot[3] = {0.0, 0.0, 0.0};
282 mutable vector<double> m_partialMolarVolumes;
305 enum class CoeffSource { EoS, CritProps, Database };
Header file for class Cantera::Array2D.
Header file for a derived class of ThermoPhase that handles non-ideal mixtures based on the fugacity ...
A map of string keys to values whose type can vary at runtime.
A class for 2D arrays stored in column-major (Fortran-compatible) form.
This is a filter class for ThermoPhase that implements some preparatory steps for efficiently handlin...
Implementation of a multi-species Peng-Robinson equation of state.
double thermalExpansionCoeff() const override
Return the volumetric thermal expansion coefficient. Units: 1/K.
double densSpinodalLiquid() const override
Return the value of the density at the liquid spinodal point (on the liquid side) for the current tem...
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
void setSpeciesCoeffs(const string &species, double a, double b, double w)
Set the pure fluid interaction parameters for a species.
double sresid() const override
Calculate the deviation terms for the total entropy of the mixture from the ideal gas mixture.
double soundSpeed() const override
Return the speed of sound. Units: m/s.
double pressure() const override
Return the thermodynamic pressure (Pa).
map< string, map< string, double > > m_binaryParameters
Explicitly-specified binary interaction parameters, to enable serialization.
void getSpeciesParameters(const string &name, AnyMap &speciesNode) const override
Get phase-specific parameters of a Species object such that an identical one could be reconstructed a...
static const double omega_b
Omega constant: b0 (= omega_b) used in Peng-Robinson equation of state.
vector< double > m_pp
Temporary storage - length = m_kk.
void setBinaryCoeffs(const string &species_i, const string &species_j, double a)
Set values for the interaction parameter between two species.
vector< double > m_dpdni
Vector of derivatives of pressure with respect to mole number.
string type() const override
String indicating the thermodynamic model implemented.
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
double densSpinodalGas() const override
Return the value of the density at the gas spinodal point (on the gas side) for the current temperatu...
double m_dpdT
The derivative of the pressure with respect to the temperature.
void getPartialMolarVolumes(double *vbar) const override
Return an array of partial molar volumes for the species in the mixture.
double cv_mole() const override
Molar heat capacity at constant volume. Units: J/kmol/K.
double liquidVolEst(double T, double &pres) const override
Estimate for the molar volume of the liquid.
int solveCubic(double T, double pres, double a, double b, double aAlpha, double Vroot[3]) const
Prepare variables and call the function to solve the cubic equation of state.
double daAlpha_dT() const
Calculate temperature derivative .
void calculatePressureDerivatives() const
Calculate and at the current conditions.
double isothermalCompressibility() const override
Returns the isothermal compressibility. Units: 1/Pa.
double hresid() const override
Calculate the deviation terms for the total enthalpy of the mixture from the ideal gas mixture.
double m_aAlpha_mix
Value of in the equation of state.
double m_a
Value of in the equation of state.
static const double omega_a
Omega constant: a0 (= omega_a) used in Peng-Robinson equation of state.
static const double omega_vc
Omega constant for the critical molar volume.
vector< CoeffSource > m_coeffSource
For each species, specifies the source of the a, b, and omega coefficients.
void calculateAB(double &aCalc, double &bCalc, double &aAlpha) const
Calculate the , , and parameters given the temperature.
void getPartialMolarIntEnergies(double *ubar) const override
Return an array of partial molar internal energies for the species in the mixture.
void updateMixingExpressions() override
Update the , , and parameters.
double cp_mole() const override
Molar heat capacity at constant pressure. Units: J/kmol/K.
double m_b
Value of in the equation of state.
void getPartialMolarCp(double *cpbar) const override
Calculate species-specific molar specific heats.
double standardConcentration(size_t k=0) const override
Returns the standard concentration , which is used to normalize the generalized concentration.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
double dpdVCalc(double T, double molarVol, double &presCalc) const override
Calculate the pressure and the pressure derivative given the temperature and the molar volume.
void getActivityCoefficients(double *ac) const override
Get the array of non-dimensional activity coefficients at the current solution temperature,...
double d2aAlpha_dT2() const
Calculate second derivative .
double densityCalc(double T, double pressure, int phase, double rhoguess) override
Calculates the density given the temperature and the pressure and a guess at the density.
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies of the species in the solution.
double m_dpdV
The derivative of the pressure with respect to the volume.
double speciesCritTemperature(double a, double b) const
Calculate species-specific critical temperature.
vector< double > m_acentric
acentric factor for each species, length m_kk
shared_ptr< Species > species(const string &name) const
Return the Species object for the named species.
string name() const
Return the name of the phase.
Namespace for the Cantera kernel.