MFC
Exascale flow solver
Loading...
Searching...
No Matches
m_viscous Module Reference

Computes viscous stress tensors and diffusive flux contributions for the Navier–Stokes equations. More...

Functions/Subroutines

impure subroutine, public s_initialize_viscous_module
 Allocates and populates the viscous Reynolds number arrays and transfers data to the GPU.
subroutine, public s_compute_viscous_stress_cylindrical_boundary (q_prim_vf, grad_x_vf, grad_y_vf, grad_z_vf, tau_re_vf, ix, iy, iz)
 The purpose of this subroutine is to compute the viscous.
subroutine, public s_get_viscous (ql_prim_rsx_vf, ql_prim_rsy_vf, ql_prim_rsz_vf, dql_prim_dx_n, dql_prim_dy_n, dql_prim_dz_n, ql_prim, qr_prim_rsx_vf, qr_prim_rsy_vf, qr_prim_rsz_vf, dqr_prim_dx_n, dqr_prim_dy_n, dqr_prim_dz_n, qr_prim, q_prim_qp, dq_prim_dx_qp, dq_prim_dy_qp, dq_prim_dz_qp, ix, iy, iz)
 Computes viscous terms.
subroutine s_reconstruct_cell_boundary_values_visc (v_vf, vl_x, vl_y, vl_z, vr_x, vr_y, vr_z, norm_dir, vl_prim_vf, vr_prim_vf, ix, iy, iz)
 Reconstructs left and right cell-boundary values of viscous primitive variables using WENO or MUSCL.
subroutine, public s_reconstruct_cell_boundary_values_visc_deriv (v_vf, vl_x, vl_y, vl_z, vr_x, vr_y, vr_z, norm_dir, vl_prim_vf, vr_prim_vf, ix, iy, iz)
 Reconstructs left and right cell-boundary values of viscous primitive variable derivatives using WENO or MUSCL.
subroutine s_apply_scalar_divergence_theorem (vl_vf, vr_vf, dv_ds_vf, norm_dir, ix, iy, iz, iv_in, dl, dim, buff_size_in)
 The purpose of this subroutine is to employ the inputted left and right cell-boundary integral-averaged variables to compute the relevant cell-average first-order spatial derivatives in the x-, y- or z-direction by means of the scalar divergence theorem.
subroutine s_compute_fd_gradient (var, grad_x, grad_y, grad_z)
 Computes the scalar gradient fields via finite differences.
subroutine, public s_compute_viscous_stress_tensor (viscous_stress_tensor, q_prim_vf, dynamic_viscosity, i, j, k)
 Computes the viscous stress tensor at a single grid cell using finite-difference velocity gradients.
impure subroutine, public s_finalize_viscous_module ()
 Deallocates the viscous Reynolds number arrays.

Variables

type(int_bounds_infoiv
type(int_bounds_infois1_viscous
type(int_bounds_infois2_viscous
type(int_bounds_infois3_viscous
real(wp), dimension(:, :), allocatable res_viscous

Detailed Description

Computes viscous stress tensors and diffusive flux contributions for the Navier–Stokes equations.

Function/Subroutine Documentation

◆ s_apply_scalar_divergence_theorem()

subroutine m_viscous::s_apply_scalar_divergence_theorem ( type(scalar_field), dimension(iv%beg:iv%end), intent(in) vl_vf,
type(scalar_field), dimension(iv%beg:iv%end), intent(in) vr_vf,
type(scalar_field), dimension(iv%beg:iv%end), intent(inout) dv_ds_vf,
integer, intent(in) norm_dir,
type(int_bounds_info), intent(in) ix,
type(int_bounds_info), intent(in) iy,
type(int_bounds_info), intent(in) iz,
type(int_bounds_info), intent(in) iv_in,
real(wp), dimension(-buff_size_in:dim + buff_size_in), intent(in) dl,
integer, intent(in) dim,
integer, intent(in) buff_size_in )

The purpose of this subroutine is to employ the inputted left and right cell-boundary integral-averaged variables to compute the relevant cell-average first-order spatial derivatives in the x-, y- or z-direction by means of the scalar divergence theorem.

Parameters
vL_vfLeft cell-boundary integral averages
vR_vfRight cell-boundary integral averages
dv_ds_vfCell-average first-order spatial derivatives
norm_dirSplitting coordinate direction
ixIndex bounds in the x-direction
iyIndex bounds in the y-direction
izIndex bounds in the z-direction
iv_inVariable index bounds
dLCell width array
dimDimension size
buff_size_inBuffer layer size

Definition at line 3372 of file m_viscous.fpp.f90.

Here is the caller graph for this function:

◆ s_compute_fd_gradient()

subroutine m_viscous::s_compute_fd_gradient ( type(scalar_field), intent(in) var,
type(scalar_field), intent(inout) grad_x,
type(scalar_field), intent(inout) grad_y,
type(scalar_field), intent(inout) grad_z )

Computes the scalar gradient fields via finite differences.

Parameters
varVariable to compute derivative of
grad_xFirst coordinate direction component of the derivative
grad_ySecond coordinate direction component of the derivative
grad_zThird coordinate direction component of the derivative

Definition at line 3650 of file m_viscous.fpp.f90.

Here is the caller graph for this function:

◆ s_compute_viscous_stress_cylindrical_boundary()

subroutine, public m_viscous::s_compute_viscous_stress_cylindrical_boundary ( type(scalar_field), dimension(sys_size), intent(in) q_prim_vf,
type(scalar_field), dimension(num_dims), intent(in) grad_x_vf,
type(scalar_field), dimension(num_dims), intent(in) grad_y_vf,
type(scalar_field), dimension(num_dims), intent(in) grad_z_vf,
type(scalar_field), dimension(1:sys_size), intent(inout) tau_re_vf,
type(int_bounds_info), intent(in) ix,
type(int_bounds_info), intent(in) iy,
type(int_bounds_info), intent(in) iz )

The purpose of this subroutine is to compute the viscous.

Definition at line 433 of file m_viscous.fpp.f90.

◆ s_compute_viscous_stress_tensor()

subroutine, public m_viscous::s_compute_viscous_stress_tensor ( real(wp), dimension(1:3, 1:3), intent(inout) viscous_stress_tensor,
type(scalar_field), dimension(1:sys_size), intent(in) q_prim_vf,
real(wp), intent(in) dynamic_viscosity,
integer, intent(in) i,
integer, intent(in) j,
integer, intent(in) k )

Computes the viscous stress tensor at a single grid cell using finite-difference velocity gradients.

Definition at line 4294 of file m_viscous.fpp.f90.

◆ s_finalize_viscous_module()

impure subroutine, public m_viscous::s_finalize_viscous_module

Deallocates the viscous Reynolds number arrays.

Definition at line 4370 of file m_viscous.fpp.f90.

Here is the caller graph for this function:

◆ s_get_viscous()

subroutine, public m_viscous::s_get_viscous ( real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) ql_prim_rsx_vf,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) ql_prim_rsy_vf,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) ql_prim_rsz_vf,
type(vector_field), dimension(1:num_dims), intent(inout) dql_prim_dx_n,
type(vector_field), dimension(1:num_dims), intent(inout) dql_prim_dy_n,
type(vector_field), dimension(1:num_dims), intent(inout) dql_prim_dz_n,
type(vector_field), dimension(num_dims), intent(inout) ql_prim,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) qr_prim_rsx_vf,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) qr_prim_rsy_vf,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) qr_prim_rsz_vf,
type(vector_field), dimension(1:num_dims), intent(inout) dqr_prim_dx_n,
type(vector_field), dimension(1:num_dims), intent(inout) dqr_prim_dy_n,
type(vector_field), dimension(1:num_dims), intent(inout) dqr_prim_dz_n,
type(vector_field), dimension(num_dims), intent(inout) qr_prim,
type(vector_field), intent(in) q_prim_qp,
type(vector_field), dimension(1), intent(inout) dq_prim_dx_qp,
type(vector_field), dimension(1), intent(inout) dq_prim_dy_qp,
type(vector_field), dimension(1), intent(inout) dq_prim_dz_qp,
type(int_bounds_info), intent(in) ix,
type(int_bounds_info), intent(in) iy,
type(int_bounds_info), intent(in) iz )

Computes viscous terms.

Parameters
qL_prim_rsx_vfLeft reconstructed primitive variables in x
qL_prim_rsy_vfLeft reconstructed primitive variables in y
qL_prim_rsz_vfLeft reconstructed primitive variables in z
dqL_prim_dx_nLeft primitive x-derivative
dqL_prim_dy_nLeft primitive y-derivative
dqL_prim_dz_nLeft primitive z-derivative
qL_primLeft cell-boundary primitive variables
qR_prim_rsx_vfRight reconstructed primitive variables in x
qR_prim_rsy_vfRight reconstructed primitive variables in y
qR_prim_rsz_vfRight reconstructed primitive variables in z
dqR_prim_dx_nRight primitive x-derivative
dqR_prim_dy_nRight primitive y-derivative
dqR_prim_dz_nRight primitive z-derivative
qR_primRight cell-boundary primitive variables
q_prim_qpCell-averaged primitive variables
dq_prim_dx_qpCell-averaged primitive x-derivative
dq_prim_dy_qpCell-averaged primitive y-derivative
dq_prim_dz_qpCell-averaged primitive z-derivative
ixIndex bounds in the x-direction
iyIndex bounds in the y-direction
izIndex bounds in the z-direction

Definition at line 1430 of file m_viscous.fpp.f90.

Here is the call graph for this function:

◆ s_initialize_viscous_module()

impure subroutine, public m_viscous::s_initialize_viscous_module

Allocates and populates the viscous Reynolds number arrays and transfers data to the GPU.

Definition at line 356 of file m_viscous.fpp.f90.

Here is the caller graph for this function:

◆ s_reconstruct_cell_boundary_values_visc()

subroutine m_viscous::s_reconstruct_cell_boundary_values_visc ( type(scalar_field), dimension(iv%beg:iv%end), intent(in) v_vf,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) vl_x,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) vl_y,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) vl_z,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) vr_x,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) vr_y,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, 1:), intent(inout) vr_z,
integer, intent(in) norm_dir,
type(scalar_field), dimension(iv%beg:iv%end), intent(inout) vl_prim_vf,
type(scalar_field), dimension(iv%beg:iv%end), intent(inout) vr_prim_vf,
type(int_bounds_info), intent(in) ix,
type(int_bounds_info), intent(in) iy,
type(int_bounds_info), intent(in) iz )

Reconstructs left and right cell-boundary values of viscous primitive variables using WENO or MUSCL.

Definition at line 2797 of file m_viscous.fpp.f90.

Here is the caller graph for this function:

◆ s_reconstruct_cell_boundary_values_visc_deriv()

subroutine, public m_viscous::s_reconstruct_cell_boundary_values_visc_deriv ( type(scalar_field), dimension(iv%beg:iv%end), intent(in) v_vf,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, iv%beg:), intent(inout) vl_x,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, iv%beg:), intent(inout) vl_y,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, iv%beg:), intent(inout) vl_z,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, iv%beg:), intent(inout) vr_x,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, iv%beg:), intent(inout) vr_y,
real(wp), dimension(idwbuff(1)%beg:, idwbuff(2)%beg:, idwbuff(3)%beg:, iv%beg:), intent(inout) vr_z,
integer, intent(in) norm_dir,
type(scalar_field), dimension(iv%beg:iv%end), intent(inout) vl_prim_vf,
type(scalar_field), dimension(iv%beg:iv%end), intent(inout) vr_prim_vf,
type(int_bounds_info), intent(in) ix,
type(int_bounds_info), intent(in) iy,
type(int_bounds_info), intent(in) iz )

Reconstructs left and right cell-boundary values of viscous primitive variable derivatives using WENO or MUSCL.

Definition at line 3077 of file m_viscous.fpp.f90.

Variable Documentation

◆ is1_viscous

type(int_bounds_info) m_viscous::is1_viscous

Definition at line 327 of file m_viscous.fpp.f90.

◆ is2_viscous

type(int_bounds_info) m_viscous::is2_viscous

Definition at line 327 of file m_viscous.fpp.f90.

◆ is3_viscous

type(int_bounds_info) m_viscous::is3_viscous

Definition at line 327 of file m_viscous.fpp.f90.

◆ iv

type(int_bounds_info) m_viscous::iv

Definition at line 326 of file m_viscous.fpp.f90.

◆ res_viscous

real(wp), dimension(:, :), allocatable m_viscous::res_viscous

Definition at line 340 of file m_viscous.fpp.f90.