|
version 0.6.0
|
|
|
version 0.6.0
|
This section describes the models and equations solved by Notus.
Fluid-flows models that can be solved are:
The International System of Units is used.
The incompressible Navier-Stokes equations take the following form:
\[ \rho \Big( \frac {\partial \mathbf{u}} {\partial t} + \mathbf{u} \cdot \nabla \mathbf{u} \Big) = - \nabla p + \nabla \cdot \Big( \mu \left( {\nabla \mathbf{u} + \nabla^\top \mathbf{u}} \right) \Big) + \mathbf{f} \\ \nabla \cdot \mathbf{u} = 0 \]
where \( \rho \) is the density, \( \mathbf{u} \) is the velocity, \( p \) is the pressure, \( \mu \) is the dynamic viscosity, and \( \mathbf{f} \) is a source term.
The compressible Navier-Stokes equations take the following form for the mass and momentum conservation equation, respectively,
\[ \frac{\partial \rho}{\partial t} + \boldsymbol{u} \cdot \boldsymbol{\nabla} \rho = -\rho \boldsymbol{\nabla} \cdot \boldsymbol{u} \\ \rho \Big( \frac {\partial \mathbf{u}} {\partial t} + \mathbf{u} \cdot \nabla \mathbf{u} \Big) = - \nabla p + \nabla \cdot \Big( \mu \left( {\nabla \mathbf{u} + \nabla^\top \mathbf{u}} \right) \Big) - \frac{2}{3} \boldsymbol{\nabla} \left(\mu\boldsymbol{\nabla}\cdot \boldsymbol{u} \right) + \mathbf{f} \]
These conservation equations are written in terms of primitive variables. This formulation considers the Stokes’ hypothesis as well as the Newtonian fluid hypothesis (as for incompressible flows).
Optionally, the following terms can be added under both incompressible and compressible hypothesis:
The gravity term
\[ \rho \mathbf{g} \]
where \( \mathbf{g} \) is the gravity field. The gravity field is assumed to be uniform. In this case, for incompressible flows, the Boussinesq approximation is used, i.e. \( \rho \) is taken as constant except in the gravity term. For compressible flows, an equation of state is used and \( \rho \) is kept variable.
The surface tension term
Following the Continuum Surface Force model [1], this term is written as:
\[ \sigma \kappa \delta_{\Gamma} \mathbf{n} = \sigma \kappa \nabla c \]
where \( \sigma \) is the surface tension between the two phases, \( \kappa \) the curvature of the interface \( \Gamma \), \( \mathbf{n} \) the normal to the interface, \( \delta_{\Gamma} \) the Dirac function that localizes the force on the interface.
The eddy-viscosity term
For turbulent flows, the dynamic viscosity \( \mu \) is replaced by \( (\mu + \mu_t) \) where \( \mu_t \) is the turbulent viscosity computed by a suited model.
A set of boundary conditions are available. If \( \mathbf{n} \) is the normal to a boundary and \( \mathbf{\tau} \) its tangent, they can be written as:
The energy equation takes the following form:
\[ \rho c_p \left( \frac {\partial T} {\partial t} + \mathbf{u} \cdot \nabla T \right) = \nabla \cdot \left( \lambda \nabla T \right) + \Phi_d \]
where \( T \) is the temperature, \( c_p \) is the heat capacity at constant pressure, \( \lambda \) is the thermal diffusion coefficient.
For anisothermal compressible flows, the conservation of energy, expressed in \( (c_p, T) \) formulation, is used
\[ \rho c_p \left( \frac{\partial T}{\partial t} + \boldsymbol{u} \cdot \boldsymbol{\nabla} T\right)= T \beta_p \left( \frac{\partial p}{\partial t} + \boldsymbol{u} \cdot \boldsymbol{\nabla} p\right) + \boldsymbol{\nabla}\cdot (\lambda\boldsymbol{\nabla}T) + \Phi_d \]
where \( T \) is the temperature field, \( \beta_p = -\rho^{-1} \left(\partial \rho\big/\partial T\right)_p \) is the isobaric thermal expansion coefficient, \( c_p \) is the heat capacity at constant pressure, \( \lambda \) is the thermal diffusion coefficient.
Optionally, the viscous dissipation rate of energy \( \Phi_d =\) can be added under both incompressible and compressible hypothesis:
\[ \Phi_d = - \frac{2 \mu}{3} \left( \boldsymbol{\nabla}\cdot \boldsymbol{u}\right)^2 + \frac{\mu}{2} \Big( \boldsymbol{\nabla} \boldsymbol{v} + \boldsymbol{\nabla} \boldsymbol{v}^T\Big) : \Big( \boldsymbol{\nabla} \boldsymbol{v} + \boldsymbol{\nabla} \boldsymbol{v}^T\Big) \, . \]
Available boundary conditions are:
The species transport equation takes the following form:
\[ \omega_i \left( \frac {\partial S_j} {\partial t} + \mathbf{u} \cdot \nabla S_j \right) = \nabla \cdot \left( D_j \nabla S_j \right) \]
where \(S_j\) is one of the species, and \(D_j\) is the corresponding diffusion coefficient. Both passive and active species are modeled this way.
Available boundary conditions are:
When solving multiphase flow of immiscible fluids, the advection of the \(i\)-th phase is done by:
\[ \frac {\partial C_i} {\partial t} + \mathbf{u} \cdot \nabla C_i = 0 \]
where \(C_i\) is the \(i\)-th fluid volume fraction.
Depending on the model, there may be a non-advected phase filling the remaining volume fraction \( (1 - \sum_{i=1}^{n} C_i) \).
Large Eddy Simulation approach is proposed to model the subgrid energy dissipation (to be validated). The only model available so far is the mixed scale one [7] which is a combination of the Smagorinsky and Turbulent Kinetic Energy ones.
Reynolds Average approach is proposed with two- and four-equations models: \(k-\omega~STD\), \(k-\omega~SST\), and \(v^2-f\) (to be validated).
The \(k-\omega \) model is a two-equation eddy-viscosity model designed to predict turbulent flows with improved accuracy near walls and in regions with adverse pressure gradients. It combines the advantages of the \(k-\omega \) model in the near-wall region with the \(k-\epsilon \) model in the free-stream region by using a blending function [6].
The equations [5,6] read for the turbulent kinetic energy \(k\)
\[ \rho \left( \frac{\partial k}{\partial t} + u_j \frac{\partial k}{\partial x_j} \right) = \tilde{P}_k - \beta^* \rho k \omega + \frac{\partial}{\partial x_j} \left[ \left( \mu + \sigma_k \mu_t \right) \frac{\partial k}{\partial x_j} \right] \]
with \( \tilde{P}_k = \min(P_k, 10 \beta^* \rho k \omega) \), \( P_k = \mu_t S_{ij} S_{ij} \), \( \nu_t = \frac{a_1 k}{\max(a_1 \omega, S F_2)} \), \( a_1 = 0.31 \), \( S = \sqrt{2S_{ij} S_{ij}} \), \( S_{ij} = \frac{1}{2} \left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) \), \( \sigma_{k,1} = 0.85 \), \( \sigma_{k,2} = 1 \), \( \beta^* = 0.09\),
and for \( \omega \) as
\[ \rho \left( \frac{\partial \omega}{\partial t} + u_j \frac{\partial \omega}{\partial x_j} \right) = \tilde{P}_\omega - \beta \rho \omega^2 + \frac{\partial}{\partial x_j} \left[ \left( \mu + \sigma_\omega \mu_t \right) \frac{\partial \omega}{\partial x_j} \right] + 2(1 - F_1) \rho \sigma_{\omega 2} \frac{1}{\omega} \frac{\partial k}{\partial x_j} \frac{\partial \omega}{\partial x_j} \]
with \( \tilde{P}_\omega = \gamma \frac{\tilde{P}_k}{\nu_t} \), \( \gamma_1 = 0.553166667\) and \( \gamma_2 = 0.440354667\), \( \sigma_{\omega,1} = 0.5 \), \( \sigma_{\omega,2} = 0.856 \), \( \sigma_{\omega2} = 0.856\) \( \beta_1 = 0.0750\) \( \beta_2 = 0.0828\)
The first blending function controlling the transition between near-wall and outer region is defined as:
\[ F_1 = \tanh \left( \left[ \min \left( \max \left( \frac{\sqrt{k}}{\beta^* \omega y}, \frac{500 \nu}{y^2 \omega} \right), \frac{4 \rho \sigma_{\omega 2} k}{CD_{k\omega} y^2} \right) \right]^4 \right) \]
with the cross-diffusion term \( CD_{k\omega} \) computed as
\[ CD_{k\omega} = \max \left( 2 \rho \sigma_{\omega 2} \frac{1}{\omega} \frac{\partial k}{\partial x_j} \frac{\partial \omega}{\partial x_j}, 10^{-10} \right) \, . \]
\(y\) is the wall distance and \( \nu = \mu / \rho \) the kinematic viscosity.
The second blending function is defined as:
\[ F_2 = \tanh \left( \left[ \max \left( \frac{2 \sqrt{k}}{\beta^* \omega y}, \frac{500 \nu}{y^2 \omega} \right) \right]^2 \right) \, . \]
The blending of a variable \( \phi \) is made by \( \phi = \phi_1 F_1 + \phi_2 (1 - F_1) \) with the two given parameters \( \phi_1\) and \( \phi_2\). The parameters \( \beta \sigma_k \sigma_\omega \) are computed with this blending.
All physical properties are considered as constant per phase, except density. Density may vary with temperature and species concentration:
\[ \rho = \rho_0 \Big( 1 - \beta_T ( T - T_{\text{ref}} ) + \sum_{j=1}^{m} \beta_{S_j} ( S_j - S_{j,\text{ref}} ) \Big) \]
where \( T_{\text{ref}}, S_{1,\text{ref}}, \ldots, S_{m,\text{ref}} \) are the reference temperature and references concentrations, respectively.
In multiphase flows, physical properties ( \( \rho_0, \mu, c_p, \lambda, \text{etc.}) \) depend on volume fractions \( C_i \) according to the following law:
\[ \phi = \phi_0 (1 - \sum_{i=1}^{n} C_i) + \sum_{i=1}^{n} \phi_i C_i \]
where \( \phi = \rho, \mu, c_p, \lambda, etc. \), the index 0 represent the non-advected phase.
To close the system of equations introduced above for compressible flows, in addition to specifying the initial and boundary conditions, an equation of state (EoS) for the density and other material properties \( x \), e.g. isothermal compressibility \( \chi_t \), thermal expansion \( \beta_p \), speed of sound \( c \) or heat capacity \( c_v, c_p \) must be chosen as
\[ x = \mathrm{EoS}(T,p) \, . \]
The choice of an EoS will induce the computation of all the thermophysical material properties of a fluid according to the chosen equation.
Constant EoS
The constant EoS considers constant material properties as defined by the user.
Perfect gas EoS
The perfect gas law is commonly used to model thermodynamic variations at atmospheric conditions. It links pressure, temperature, and density as
\[ \frac{P}{\rho} = R T \]
where \(P\) is the pressure, \(T\) is the absolute temperature, \(\rho\) is the density and \(R\) is the universal gas constant. This EoS may not accurately capture real-gas effects under high pressures, for example.
Material properties are updated according to this law as:
\[ \chi_T \!=\! 1/p \, , \beta_p \!=\! 1 / T \, , \, \text{and} \, , {c}^2 \!=\! \gamma p / \rho \, . \]
Peng-Robinson EoS
The Peng–Robinson (PR) EoS is a widely used cubic equation of state designed to accurately predict the thermodynamic properties of real gases, especially near the critical point.
\[ P = \frac{\rho RT}{1 - b \rho} - \frac{a(T) \rho^2}{1 + 2b \rho - b^2 \rho^2} \]
where \(P\) is the pressure, \(T\) is the absolute temperature, \(\rho\) is the density, \(R\) is the universal gas constant, \(a(T)\) is a temperature-dependent attractive parameter, \(b\) is the covolume parameter accounting for the finite size of molecules.
The temperature-dependent parameter \(a(T)\) is defined as:
\[ a(T) = a_c \cdot \alpha(T) \]
with
\[ a_c = 0.45724 \frac{R^2 T_c^2}{P_c} , \quad b = 0.07780 \frac{R T_c}{P_c} \, ,\\ \alpha(T) = \left[1 + \kappa \left(1 - \sqrt{\frac{T}{T_c}}\right)\right]^2 \,, \quad \kappa = 0.37464 + 1.54226\omega - 0.26992\omega^2 \, . \]
Here, \(T_c\) and \(P_c\) are the critical temperature and pressure of the fluid, respectively, and \(\omega\) is the acentric factor.
Material properties ( \(\chi_T \), \( \beta_p \), \( c_p, c_v\), \( c \)) are updated according to PR law. For more informations see [2]
Linear acoustic
In the context of linear acoustics (d'Alembert equation), a simplified equation of state is often used to relate pressure and density variations. It assumes small perturbations around a reference state \((T_0,P_0)\) and is given by
\[ \Delta p = c_0^2 \Delta \rho \, , \]
where \(\Delta p\) is the pressure perturbation, \(\Delta \rho\) is the density perturbation, and \(c_0=\sqrt{\gamma R T_0}\) is the constant speed of sound of the medium.
Refprop NIST library
For more accurate thermodynamic modeling, especially for real-fluids effects, the REFPROP (Reference Fluid Thermodynamic and Transport Properties) database [3] can be used. REFPROP provides high-precision equations of state and transport properties for a wide range of pure fluids and mixtures, based on experimental data and advanced thermodynamic models. REFPROP is a licensed software and is not freely available; it must be purchased from NIST for official use.
Enabling refprop in Notus, provided you have the operating license, replaces the calculation of all material properties made by an implemented EoS within Notus in favor of refprop. This choice is very accurate from a thermodynamic point of view, at the cost of a decrease in computational efficiency.
Coolprop library
CoolProp [4] is an open-source thermophysical property library that provides accurate thermodynamic and transport properties for a wide range of pure fluids and mixtures. It serves as a licence free alternative to commercial databases, like REFPROP for example.
[1] J. U. Brackbill, D. B. Korma, and C. Zemach, A Continuum Method for Modeling Surface Tension, Journal of Computational Physics, 100, 335-354 (1992).
[2] B. Poling, J. Prausnitz, and J. O'Connell, The properties of gases and liquids, McGraw-Hill, New York, (2001)
[3] E. W. Lemmon, I. H. Bell, M. L. Huber, and M. O. McLinden. NIST Standard Reference Database 23: Reference Fluid Thermodynamic and Transport Properties-REFPROP, Version 10., National Institute of Standards and Technology, (2018). url
[4] I.H. Bell, J. Wronski, S. Quoilin, V. Lemort, Pure and Pseudo-pure Fluid Thermophysical Property Evaluation and the Open-Source Thermophysical Property Library CoolProp, Industrial & Engineering Chemistry Research, 53-6, 2498-2508, (2014). url
[5] F. R. Menter. Two-equation eddy-viscosity turbulence models for engineering applications. AIAA Journal, 32(8):1598–1605, (1994).
[6] D. C. Wilcox. Turbulence Modeling for CFD. DCW Industries, La Canada, California, 3rd edition, (2008).
[7] P. Sagaut, Large Eddy Simulation for Incompressible Flows: An Introduction, Springer-Verlag Berlin, (2001), DOI : 10.1007/978-3-662-04416-2
[8] D. C. Wilcox, Formulation of the k-w Turbulence Model Revisited, AIAA Journal, Vol. 46, No. 11, (2008)
1.11.0