Nektar++
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | Friends | List of all members
Nektar::VanDerWaalsEoS Class Reference

van der Waals equation of state: p = RT/(1/rho - b) - a * rho^2 with a = 27/64 * (R*Tc)^2 / Pc b = 1/8 * (R*Tc) / Pc More...

#include <VanDerWaalsEoS.h>

Inheritance diagram for Nektar::VanDerWaalsEoS:
[legend]

Static Public Member Functions

static EquationOfStateSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

NekDouble GetTemperature (const NekDouble &rho, const NekDouble &e) final
 Calculate the temperature. More...
 
vec_t GetTemperature (const vec_t &rho, const vec_t &e) final
 
NekDouble GetPressure (const NekDouble &rho, const NekDouble &e) final
 Calculate the pressure. More...
 
vec_t GetPressure (const vec_t &rho, const vec_t &e) final
 
NekDouble v_GetEntropy (const NekDouble &rho, const NekDouble &e) final
 
NekDouble v_GetDPDrho_e (const NekDouble &rho, const NekDouble &e) final
 
NekDouble v_GetDPDe_rho (const NekDouble &rho, const NekDouble &e) final
 
NekDouble v_GetEFromRhoP (const NekDouble &rho, const NekDouble &p) final
 
NekDouble v_GetRhoFromPT (const NekDouble &rho, const NekDouble &p) final
 
- Protected Member Functions inherited from Nektar::EquationOfState
 EquationOfState (const LibUtilities::SessionReaderSharedPtr &pSession)
 Constructor. More...
 
 EquationOfState (const NekDouble &gamma, const NekDouble &gasConstant)
 Programmatic Constructor. More...
 
virtual NekDouble v_GetSoundSpeed (const NekDouble &rho, const NekDouble &e)
 

Protected Attributes

NekDouble m_a
 
NekDouble m_b
 
- Protected Attributes inherited from Nektar::EquationOfState
NekDouble m_gamma
 
NekDouble m_gasConstant
 
NekDouble m_gammaMone
 
NekDouble m_gammaMoneOgasConst
 

Private Member Functions

 VanDerWaalsEoS (const LibUtilities::SessionReaderSharedPtr &pSession)
 
 ~VanDerWaalsEoS (void)
 
template<class T , typename = typename std::enable_if< std::is_floating_point<T>::value || tinysimd::is_vector_floating_point<T>::value>::type>
GetTemperatureKernel (const T &rho, const T &e)
 
template<class T , typename = typename std::enable_if< std::is_floating_point<T>::value || tinysimd::is_vector_floating_point<T>::value>::type>
GetPressureKernel (const T &rho, const T &e)
 

Friends

class MemoryManager< VanDerWaalsEoS >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::EquationOfState
virtual ~EquationOfState ()
 
NekDouble GetSoundSpeed (const NekDouble &rho, const NekDouble &e)
 Calculate the sound speed. More...
 
NekDouble GetEntropy (const NekDouble &rho, const NekDouble &e)
 Calculate the entropy. More...
 
NekDouble GetDPDrho_e (const NekDouble &rho, const NekDouble &e)
 Calculate the partial derivative of P(rho,e) with respect to rho. More...
 
NekDouble GetDPDe_rho (const NekDouble &rho, const NekDouble &e)
 Calculate the partial derivative of P(rho,e) with respect to e. More...
 
NekDouble GetEFromRhoP (const NekDouble &rho, const NekDouble &p)
 Obtain the internal energy from rho and P. More...
 
NekDouble GetRhoFromPT (const NekDouble &p, const NekDouble &T)
 Obtain the density from P and T. More...
 

Detailed Description

van der Waals equation of state: p = RT/(1/rho - b) - a * rho^2 with a = 27/64 * (R*Tc)^2 / Pc b = 1/8 * (R*Tc) / Pc

Definition at line 49 of file VanDerWaalsEoS.h.

Constructor & Destructor Documentation

◆ VanDerWaalsEoS()

Nektar::VanDerWaalsEoS::VanDerWaalsEoS ( const LibUtilities::SessionReaderSharedPtr pSession)
private

Definition at line 49 of file VanDerWaalsEoS.cpp.

51  : EquationOfState(pSession)
52 {
53  NekDouble Tcrit, Pcrit;
54  pSession->LoadParameter("Tcrit", Tcrit);
55  pSession->LoadParameter("Pcrit", Pcrit);
56 
57  m_a = 27.0 / 64.0 * m_gasConstant * m_gasConstant * Tcrit * Tcrit / Pcrit;
58  m_b = 1.0 / 8.0 * m_gasConstant * Tcrit / Pcrit;
59 }
EquationOfState(const LibUtilities::SessionReaderSharedPtr &pSession)
Constructor.
double NekDouble

References m_a, m_b, and Nektar::EquationOfState::m_gasConstant.

◆ ~VanDerWaalsEoS()

Nektar::VanDerWaalsEoS::~VanDerWaalsEoS ( void  )
inlineprivate

Definition at line 91 of file VanDerWaalsEoS.h.

91 {};

Member Function Documentation

◆ create()

static EquationOfStateSharedPtr Nektar::VanDerWaalsEoS::create ( const LibUtilities::SessionReaderSharedPtr pSession)
inlinestatic

Creates an instance of this class.

Definition at line 55 of file VanDerWaalsEoS.h.

57  {
60  return p;
61  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< EquationOfState > EquationOfStateSharedPtr
A shared pointer to an equation of state object.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

◆ GetPressure() [1/2]

NekDouble Nektar::VanDerWaalsEoS::GetPressure ( const NekDouble rho,
const NekDouble e 
)
finalprotectedvirtual

Calculate the pressure.

Implements Nektar::EquationOfState.

Definition at line 72 of file VanDerWaalsEoS.cpp.

73 {
74  return GetPressureKernel(rho, e);
75 }
T GetPressureKernel(const T &rho, const T &e)

References GetPressureKernel().

◆ GetPressure() [2/2]

vec_t Nektar::VanDerWaalsEoS::GetPressure ( const vec_t rho,
const vec_t e 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 77 of file VanDerWaalsEoS.cpp.

78 {
79  return GetPressureKernel(rho, e);
80 }

References GetPressureKernel().

◆ GetPressureKernel()

template<class T , typename = typename std::enable_if< std::is_floating_point<T>::value || tinysimd::is_vector_floating_point<T>::value>::type>
T Nektar::VanDerWaalsEoS::GetPressureKernel ( const T &  rho,
const T &  e 
)
inlineprivate

Definition at line 104 of file VanDerWaalsEoS.h.

105  {
106  return (e + m_a * rho) * (m_gamma - 1) / (1.0 / rho - m_b) -
107  m_a * rho * rho;
108  }

References m_a, m_b, and Nektar::EquationOfState::m_gamma.

Referenced by GetPressure().

◆ GetTemperature() [1/2]

NekDouble Nektar::VanDerWaalsEoS::GetTemperature ( const NekDouble rho,
const NekDouble e 
)
finalprotectedvirtual

Calculate the temperature.

Implements Nektar::EquationOfState.

Definition at line 61 of file VanDerWaalsEoS.cpp.

63 {
64  return GetTemperatureKernel(rho, e);
65 }
T GetTemperatureKernel(const T &rho, const T &e)

References GetTemperatureKernel().

Referenced by v_GetEntropy().

◆ GetTemperature() [2/2]

vec_t Nektar::VanDerWaalsEoS::GetTemperature ( const vec_t rho,
const vec_t e 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 67 of file VanDerWaalsEoS.cpp.

68 {
69  return GetTemperatureKernel(rho, e);
70 }

References GetTemperatureKernel().

◆ GetTemperatureKernel()

template<class T , typename = typename std::enable_if< std::is_floating_point<T>::value || tinysimd::is_vector_floating_point<T>::value>::type>
T Nektar::VanDerWaalsEoS::GetTemperatureKernel ( const T &  rho,
const T &  e 
)
inlineprivate

Definition at line 96 of file VanDerWaalsEoS.h.

97  {
98  return (e + m_a * rho) * (m_gamma - 1) / m_gasConstant;
99  }

References m_a, Nektar::EquationOfState::m_gamma, and Nektar::EquationOfState::m_gasConstant.

Referenced by GetTemperature().

◆ v_GetDPDe_rho()

NekDouble Nektar::VanDerWaalsEoS::v_GetDPDe_rho ( const NekDouble rho,
const NekDouble e 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 103 of file VanDerWaalsEoS.cpp.

105 {
106  boost::ignore_unused(e);
107  return (m_gamma - 1) / (1.0 / rho - m_b);
108 }

References m_b, and Nektar::EquationOfState::m_gamma.

◆ v_GetDPDrho_e()

NekDouble Nektar::VanDerWaalsEoS::v_GetDPDrho_e ( const NekDouble rho,
const NekDouble e 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 91 of file VanDerWaalsEoS.cpp.

93 {
94  NekDouble result;
95 
96  result = (m_gamma - 1) * (e + 2 * m_a * rho - m_a * m_b * rho * rho);
97  result = result / ((1 - m_b * rho) * (1 - m_b * rho));
98  result = result - 2 * m_a * rho;
99 
100  return result;
101 }

References m_a, m_b, and Nektar::EquationOfState::m_gamma.

◆ v_GetEFromRhoP()

NekDouble Nektar::VanDerWaalsEoS::v_GetEFromRhoP ( const NekDouble rho,
const NekDouble p 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 110 of file VanDerWaalsEoS.cpp.

112 {
113  return (p + m_a * rho * rho) * (1.0 / rho - m_b) / (m_gamma - 1) -
114  m_a * rho;
115 }

References m_a, m_b, Nektar::EquationOfState::m_gamma, and CellMLToNektar.cellml_metadata::p.

◆ v_GetEntropy()

NekDouble Nektar::VanDerWaalsEoS::v_GetEntropy ( const NekDouble rho,
const NekDouble e 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 82 of file VanDerWaalsEoS.cpp.

83 {
84  NekDouble T = GetTemperature(rho, e);
85  NekDouble sIg =
86  m_gasConstant / (m_gamma - 1) * log(T) - m_gasConstant * log(rho);
87 
88  return sIg + m_gasConstant * log(1 - m_b * rho);
89 }
NekDouble GetTemperature(const NekDouble &rho, const NekDouble &e) final
Calculate the temperature.
scalarT< T > log(scalarT< T > in)
Definition: scalar.hpp:300

References GetTemperature(), tinysimd::log(), m_b, Nektar::EquationOfState::m_gamma, and Nektar::EquationOfState::m_gasConstant.

◆ v_GetRhoFromPT()

NekDouble Nektar::VanDerWaalsEoS::v_GetRhoFromPT ( const NekDouble rho,
const NekDouble p 
)
finalprotectedvirtual

Implements Nektar::EquationOfState.

Definition at line 117 of file VanDerWaalsEoS.cpp.

118 {
119  // First solve for the compressibility factor Z using the cubic equation
120  // Z^3 + k1 * Z^2 + k2 * Z + k3 = 0
121  // for van der Waals:
122  // k1 = -(B+1), k2 = A, k3 = -AB
123  // where A = aP/(RT)^2, B = bP/(RT)
124  NekDouble A = m_a * p / (m_gasConstant * m_gasConstant * T * T);
125  NekDouble B = m_b * p / (m_gasConstant * T);
126 
127  NekDouble k1 = -(B + 1);
128  NekDouble k2 = A;
129  NekDouble k3 = -A * B;
130 
131  // Use ideal gas (Z=1) as starting guess for iteration
132  NekDouble Z = 1.0;
133  // Newton-Raphson iteration to find Z
134  NekDouble tol = 1e-6;
135  NekDouble maxIter = 100;
136  NekDouble residual = 1;
137  NekDouble f, df;
138  unsigned int cnt = 0;
139  while (abs(residual) > tol && cnt < maxIter)
140  {
141  f = Z * Z * Z + k1 * Z * Z + k2 * Z + k3;
142  df = 3 * Z * Z + 2 * k1 * Z + k2;
143  residual = f / df;
144  Z -= residual;
145  ++cnt;
146  }
147  if (cnt == maxIter)
148  {
149  cout << "Newton-Raphson in VanDerWaalsEoS::v_GetRhoFromPT did not "
150  "converge in "
151  << maxIter << " iterations (residual = " << residual << ")"
152  << endl;
153  }
154 
155  // Now calculate rho = p/(ZRT)
156  return p / (Z * m_gasConstant * T);
157 }
scalarT< T > abs(scalarT< T > in)
Definition: scalar.hpp:295

References tinysimd::abs(), m_a, m_b, Nektar::EquationOfState::m_gasConstant, and CellMLToNektar.cellml_metadata::p.

Friends And Related Function Documentation

◆ MemoryManager< VanDerWaalsEoS >

friend class MemoryManager< VanDerWaalsEoS >
friend

Definition at line 1 of file VanDerWaalsEoS.h.

Member Data Documentation

◆ className

std::string Nektar::VanDerWaalsEoS::className
static
Initial value:
=
"VanDerWaals", VanDerWaalsEoS::create,
"Van der Waals equation of state.")
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:198
static EquationOfStateSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession)
Creates an instance of this class.
EquationOfStateFactory & GetEquationOfStateFactory()
Declaration of the equation of state factory singleton.

Name of the class.

Definition at line 64 of file VanDerWaalsEoS.h.

◆ m_a

NekDouble Nektar::VanDerWaalsEoS::m_a
protected

◆ m_b

NekDouble Nektar::VanDerWaalsEoS::m_b
protected