programmer's documentation
Functions/Subroutines | Variables
cavitation Module Reference

Functions/Subroutines

subroutine init_cavitation
 Default initialization of the module variables. More...
 
subroutine cavitation_compute_source_term (pressure, voidf)
 Compute the vaporization source term $ \Gamma_V \left(\alpha, p\right) = m^+ + m^- $ using the Merkle model:

\[ m^+ = -\dfrac{C_{prod} \rho_l \min \left( p-p_V,0 \right)\alpha(1-\alpha)} {0.5\rho_lu_\infty^2t_\infty}, \]

\[ m^- = -\dfrac{C_{dest} \rho_v \max \left( p-p_V,0 \right)\alpha(1-\alpha)} {0.5\rho_lu_\infty^2t_\infty}, \]

with $ C_{prod}, C_{dest} $ empirical constants, $ t_\infty=l_\infty/u_\infty $ a reference time scale and $p_V$ the reference saturation pressure. $ l_\infty $, $ u_\infty $ and $p_V$ may be provided by the user (user function). Note that the r.h.s. of the void fraction transport equation is $ \Gamma_V/\rho_v $. More...

 
subroutine cavitation_correct_visc_turb (crom, voidf, visct)
 Modify eddy viscosity using the Reboud correction:

\[ \mu_t'= \dfrac{\rho_v + (1-\alpha)^{mcav}(\rho_l-\rho_v)}{\rho}\mu_t. \]

. More...

 

Variables

double precision, save presat
 reference saturation pressure in kg/(m s2) More...
 
double precision, save linf
 reference length scale of the flow in meters More...
 
double precision, save uinf
 reference velocity of the flow in m/s More...
 
double precision, save cprod
 constant Cprod of the vaporization source term (Merkle model) More...
 
double precision, save cdest
 constant Cdest of the condensation source term (Merkle model) More...
 
integer, save icvevm
 activation of the eddy-viscosity correction (Reboud correction) More...
 
double precision, save mcav
 constant mcav of the eddy-viscosity correction (Reboud correction) More...
 
integer, save itscvi
 implicitation in pressure of the vaporization/condensation model More...
 

Function/Subroutine Documentation

◆ cavitation_compute_source_term()

subroutine cavitation::cavitation_compute_source_term ( double precision, dimension(ncelet)  pressure,
double precision, dimension(ncelet)  voidf 
)

Compute the vaporization source term $ \Gamma_V \left(\alpha, p\right) = m^+ + m^- $ using the Merkle model:

\[ m^+ = -\dfrac{C_{prod} \rho_l \min \left( p-p_V,0 \right)\alpha(1-\alpha)} {0.5\rho_lu_\infty^2t_\infty}, \]

\[ m^- = -\dfrac{C_{dest} \rho_v \max \left( p-p_V,0 \right)\alpha(1-\alpha)} {0.5\rho_lu_\infty^2t_\infty}, \]

with $ C_{prod}, C_{dest} $ empirical constants, $ t_\infty=l_\infty/u_\infty $ a reference time scale and $p_V$ the reference saturation pressure. $ l_\infty $, $ u_\infty $ and $p_V$ may be provided by the user (user function). Note that the r.h.s. of the void fraction transport equation is $ \Gamma_V/\rho_v $.

Parameters
[in]pressurePressure array
[in]voidfVoid fraction array

◆ cavitation_correct_visc_turb()

subroutine cavitation::cavitation_correct_visc_turb ( double precision, dimension(ncelet)  crom,
double precision, dimension(ncelet)  voidf,
double precision, dimension(ncelet)  visct 
)

Modify eddy viscosity using the Reboud correction:

\[ \mu_t'= \dfrac{\rho_v + (1-\alpha)^{mcav}(\rho_l-\rho_v)}{\rho}\mu_t. \]

.

Parameters
[in]cromDensity array
[in]voidfVoid fraction array
[in,out]voidfturbulent viscosity

◆ init_cavitation()

subroutine cavitation::init_cavitation ( )

Default initialization of the module variables.