MFC:Post_process  v1.0
m_derived_variables Module Reference

This module features subroutines that allow for the derivation of numerous flow variables from the conservative and primitive ones. Currently, the available derived variables include the unadvected volume fraction, specific heat ratio, liquid stiffness, speed of sound, vorticity and the numerical Schlieren function. More...

Functions/Subroutines

subroutine, public s_initialize_derived_variables_module ()
 Computation of parameters, allocation procedures, and/or any other tasks needed to properly setup the module. More...
 
The purpose of this subroutine is to compute the finite-

difference coefficients for the centered schemes utilized in computations of first order spatial derivatives in the s-coordinate direction. The s-coordinate direction refers to the x-, y- or z-coordinate direction, depending on the subroutine's inputs. Note that coefficients of up to 4th order accuracy are available.

Parameters
qNumber of cells in the s-coordinate direction
offset_sSize of the ghost zone layer in the s-coordinate direction
s_ccLocations of the cell-centers in the s-coordinate direction
fd_coeff_sFinite-diff. coefficients in the s-coordinate direction
subroutine, public s_compute_finite_difference_coefficients (q, offset_s, s_cc, fd_coeff_s)
 
subroutine, public s_derive_unadvected_volume_fraction (q_cons_vf, q_sf)
 This subroutine is used together with the volume fraction model and when called upon, it computes the values of the unadvected volume fraction from the inputted conservative variables, q_cons_vf. The calculated values are stored in the derived flow quantity storage variable, q_sf. More...
 
subroutine, public s_derive_specific_heat_ratio (gamma_sf, q_sf)
 This subroutine receives as input the specific heat ratio function, gamma_sf, and derives from it the specific heat ratio. The latter is stored in the derived flow quantity storage variable, q_sf. More...
 
subroutine, public s_derive_liquid_stiffness (gamma_sf, pi_inf_sf, q_sf)
 This subroutine admits as inputs the specific heat ratio function and the liquid stiffness function, gamma_sf and pi_inf_sf, respectively. These are used to calculate the values of the liquid stiffness, which are stored in the derived flow quantity storage variable, q_sf. More...
 
subroutine, public s_derive_sound_speed (q_prim_vf, rho_sf, gamma_sf, pi_inf_sf, q_sf)
 This subroutine admits as inputs the primitive variables, the density, the specific heat ratio function and liquid stiffness function. It then computes from those variables the values of the speed of sound, which are stored in the derived flow quantity storage variable, q_sf. More...
 
subroutine, public s_derive_flux_limiter (i, q_prim_vf, q_sf)
 This subroutine derives the flux_limiter at cell boundary i+1/2. This is an approximation because the velocity used to determine the upwind direction is the velocity at the cell center i instead of the contact velocity at the cell boundary from the Riemann solver. More...
 
subroutine, public s_derive_curvature (i, q_prim_vf, q_sf)
 Computes the local curvatures. More...
 
subroutine s_solve_linear_system (A, b, sol, ndim)
 Computes the solution to the linear system Ax=b w/ sol = x. More...
 
subroutine, public s_derive_vorticity_component (i, q_prim_vf, q_sf)
 This subroutine receives as inputs the indicator of the component of the vorticity that should be outputted and the primitive variables. From those inputs, it proceeds to calculate values of the desired vorticity component, which are subsequently stored in derived flow quantity storage variable, q_sf. More...
 
subroutine, public s_derive_numerical_schlieren_function (q_cons_vf, rho_sf, q_sf)
 This subroutine gets as inputs the conservative variables and density. From those inputs, it proceeds to calculate the values of the numerical Schlieren function, which are subsequently stored in the derived flow quantity storage variable, q_sf. More...
 
subroutine, public s_finalize_derived_variables_module ()
 Deallocation procedures for the module. More...
 

Variables

real(kind(0d0)), dimension(:,:,:), allocatable gm_rho_sf
 Gradient magnitude (gm) of the density for each cell of the computational sub-domain. This variable is employed in the calculation of the numerical Schlieren function. More...
 
integer, private flg
 Flagging (flg) variable used to annotate the dimensionality of the dataset that is undergoing the post-process. A flag value of 1 indicates that the dataset is 3D, while a flag value of 0 indicates that it is not. This flg variable is necessary to avoid cycling through the third dimension of the flow variable(s) when the simulation is not 3D and the size of the buffer is non-zero. Note that a similar procedure does not have to be applied to the second dimension since in 1D, the buffer size is always zero. More...
 
Finite-difference (fd) coefficients in x-, y- and z-coordinate directions.

Note that because sufficient boundary information is available for all the active coordinate directions, the centered family of the finite-difference schemes is used.

real(kind(0d0)), dimension(:,:), allocatable, public fd_coeff_x
 
real(kind(0d0)), dimension(:,:), allocatable, public fd_coeff_y
 
real(kind(0d0)), dimension(:,:), allocatable, public fd_coeff_z
 

Detailed Description

This module features subroutines that allow for the derivation of numerous flow variables from the conservative and primitive ones. Currently, the available derived variables include the unadvected volume fraction, specific heat ratio, liquid stiffness, speed of sound, vorticity and the numerical Schlieren function.

Function/Subroutine Documentation

◆ s_compute_finite_difference_coefficients()

subroutine, public m_derived_variables::s_compute_finite_difference_coefficients ( integer, intent(in)  q,
type(bounds_info), intent(in)  offset_s,
real(kind(0d0)), dimension(-buff_size:q+buff_size), intent(in)  s_cc,
real(kind(0d0)), dimension(-fd_number:fd_number, -offset_s%beg:q+offset_s%end), intent(inout)  fd_coeff_s 
)

Definition at line 166 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_curvature()

subroutine, public m_derived_variables::s_derive_curvature ( integer, intent(in)  i,
type(scalar_field), dimension(sys_size), intent(in)  q_prim_vf,
real(kind(0d0)), dimension (-offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end), intent(inout)  q_sf 
)

Computes the local curvatures.

Parameters
iFluid indicator
q_prim_vfPrimitive variables
q_sfCurvature

Definition at line 513 of file m_derived_variables.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ s_derive_flux_limiter()

subroutine, public m_derived_variables::s_derive_flux_limiter ( integer, intent(in)  i,
type(scalar_field), dimension(sys_size), intent(in)  q_prim_vf,
real(kind(0d0)), dimension ( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end), intent(inout)  q_sf 
)

This subroutine derives the flux_limiter at cell boundary i+1/2. This is an approximation because the velocity used to determine the upwind direction is the velocity at the cell center i instead of the contact velocity at the cell boundary from the Riemann solver.

Parameters
iComponent indicator
q_prim_vfPrimitive variables
q_sfFlux limiter

Definition at line 417 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_liquid_stiffness()

subroutine, public m_derived_variables::s_derive_liquid_stiffness ( real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  gamma_sf,
real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  pi_inf_sf,
real(kind(0d0)), dimension( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end ), intent(inout)  q_sf 
)

This subroutine admits as inputs the specific heat ratio function and the liquid stiffness function, gamma_sf and pi_inf_sf, respectively. These are used to calculate the values of the liquid stiffness, which are stored in the derived flow quantity storage variable, q_sf.

Parameters
gamma_sfSpecific heat ratio
pi_inf_sfLiquid stiffness function
q_sfLiquid stiffness

Definition at line 305 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_numerical_schlieren_function()

subroutine, public m_derived_variables::s_derive_numerical_schlieren_function ( type(scalar_field), dimension(sys_size), intent(in)  q_cons_vf,
real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  rho_sf,
real(kind(0d0)), dimension( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end ), intent(inout)  q_sf 
)

This subroutine gets as inputs the conservative variables and density. From those inputs, it proceeds to calculate the values of the numerical Schlieren function, which are subsequently stored in the derived flow quantity storage variable, q_sf.

Parameters
q_cons_vfConservative variables
rho_sfDensity
q_sfNumerical Schlieren function

Definition at line 840 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_sound_speed()

subroutine, public m_derived_variables::s_derive_sound_speed ( type(scalar_field), dimension(sys_size), intent(in)  q_prim_vf,
real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  rho_sf,
real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  gamma_sf,
real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  pi_inf_sf,
real(kind(0d0)), dimension( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end ), intent(inout)  q_sf 
)

This subroutine admits as inputs the primitive variables, the density, the specific heat ratio function and liquid stiffness function. It then computes from those variables the values of the speed of sound, which are stored in the derived flow quantity storage variable, q_sf.

Parameters
q_prim_vfPrimitive variables
rho_sfDensity
gamma_sfSpecific heat ratio function
pi_inf_sfLiquid stiffness function
q_sfSpeed of sound

Definition at line 349 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_specific_heat_ratio()

subroutine, public m_derived_variables::s_derive_specific_heat_ratio ( real(kind(0d0)), dimension( -buff_size : m+buff_size , -buff_size : n+buff_size , -buff_size*flg : (p+buff_size)*flg ), intent(in)  gamma_sf,
real(kind(0d0)), dimension( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end ), intent(inout)  q_sf 
)

This subroutine receives as input the specific heat ratio function, gamma_sf, and derives from it the specific heat ratio. The latter is stored in the derived flow quantity storage variable, q_sf.

Parameters
gamma_sfSpecific heat ratio function
q_sfSpecific heat ratio

Definition at line 265 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_unadvected_volume_fraction()

subroutine, public m_derived_variables::s_derive_unadvected_volume_fraction ( type(scalar_field), dimension(sys_size), intent(in)  q_cons_vf,
real(kind(0d0)), dimension( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end ), intent(inout)  q_sf 
)

This subroutine is used together with the volume fraction model and when called upon, it computes the values of the unadvected volume fraction from the inputted conservative variables, q_cons_vf. The calculated values are stored in the derived flow quantity storage variable, q_sf.

Parameters
q_cons_vfConservative variables
q_sfUnadvected volume fraction

Definition at line 223 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_derive_vorticity_component()

subroutine, public m_derived_variables::s_derive_vorticity_component ( integer, intent(in)  i,
type(scalar_field), dimension(sys_size), intent(in)  q_prim_vf,
real(kind(0d0)), dimension( -offset_x%beg : m+offset_x%end , -offset_y%beg : n+offset_y%end , -offset_z%beg : p+offset_z%end ), intent(inout)  q_sf 
)

This subroutine receives as inputs the indicator of the component of the vorticity that should be outputted and the primitive variables. From those inputs, it proceeds to calculate values of the desired vorticity component, which are subsequently stored in derived flow quantity storage variable, q_sf.

Parameters
iVorticity component indicator
q_prim_vfPrimitive variables
q_sfVorticity component

Definition at line 728 of file m_derived_variables.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ s_finalize_derived_variables_module()

subroutine, public m_derived_variables::s_finalize_derived_variables_module ( )

Deallocation procedures for the module.

Definition at line 1003 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_initialize_derived_variables_module()

subroutine, public m_derived_variables::s_initialize_derived_variables_module ( )

Computation of parameters, allocation procedures, and/or any other tasks needed to properly setup the module.

Definition at line 98 of file m_derived_variables.f90.

Here is the caller graph for this function:

◆ s_solve_linear_system()

subroutine m_derived_variables::s_solve_linear_system ( real(kind(0d0)), dimension(ndim,ndim), intent(inout)  A,
real(kind(0d0)), dimension(ndim), intent(inout)  b,
real(kind(0d0)), dimension(ndim), intent(out)  sol,
integer, intent(in)  ndim 
)

Computes the solution to the linear system Ax=b w/ sol = x.

Parameters
AInput matrix
bright-hand-side
solSolution
ndimProblem size

Definition at line 658 of file m_derived_variables.f90.

Here is the caller graph for this function:

Variable Documentation

◆ fd_coeff_x

real(kind(0d0)), dimension(:,:), allocatable, public m_derived_variables::fd_coeff_x

Definition at line 75 of file m_derived_variables.f90.

◆ fd_coeff_y

real(kind(0d0)), dimension(:,:), allocatable, public m_derived_variables::fd_coeff_y

Definition at line 76 of file m_derived_variables.f90.

◆ fd_coeff_z

real(kind(0d0)), dimension(:,:), allocatable, public m_derived_variables::fd_coeff_z

Definition at line 77 of file m_derived_variables.f90.

◆ flg

integer, private m_derived_variables::flg
private

Flagging (flg) variable used to annotate the dimensionality of the dataset that is undergoing the post-process. A flag value of 1 indicates that the dataset is 3D, while a flag value of 0 indicates that it is not. This flg variable is necessary to avoid cycling through the third dimension of the flow variable(s) when the simulation is not 3D and the size of the buffer is non-zero. Note that a similar procedure does not have to be applied to the second dimension since in 1D, the buffer size is always zero.

Definition at line 80 of file m_derived_variables.f90.

◆ gm_rho_sf

real(kind(0d0)), dimension(:,:,:), allocatable m_derived_variables::gm_rho_sf

Gradient magnitude (gm) of the density for each cell of the computational sub-domain. This variable is employed in the calculation of the numerical Schlieren function.

Definition at line 65 of file m_derived_variables.f90.