|
MFC: Simulation
High-fidelity multiphase flow simulation
|
The module contains the subroutines used to compute viscous terms. More...
Functions/Subroutines | |
| impure subroutine, public | s_initialize_viscous_module |
| subroutine, public | s_compute_viscous_stress_tensor (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) |
| 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) |
| 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. | |
| impure subroutine, public | s_finalize_viscous_module () |
Variables | |
| type(int_bounds_info) | iv |
| type(int_bounds_info) | is1_viscous |
| type(int_bounds_info) | is2_viscous |
| type(int_bounds_info) | is3_viscous |
| real(wp), dimension(:, :), allocatable | res_viscous |
The module contains the subroutines used to compute viscous terms.
| 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.
| vL_vf | Left cell-boundary integral averages |
| vR_vf | Right cell-boundary integral averages |
| dv_ds_vf | Cell-average first-order spatial derivatives |
| norm_dir | Splitting coordinate direction |
| 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.
| var | Variable to compute derivative of |
| grad_x | First coordinate direction component of the derivative |
| grad_y | Second coordinate direction component of the derivative |
| grad_z | Third coordinate direction component of the derivative |
| norm | Norm of the gradient vector |
| subroutine, public m_viscous::s_compute_viscous_stress_tensor | ( | 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.
| impure subroutine, public m_viscous::s_finalize_viscous_module |
| 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.
| q_cons_vf | Cell-averaged conservative variables |
| q_prim_vf | Cell-averaged primitive variables |
| rhs_vf | Cell-averaged RHS variables |
| impure subroutine, public m_viscous::s_initialize_viscous_module |
| 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 ) |
| 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 ) |
| type(int_bounds_info) m_viscous::is1_viscous |
| type(int_bounds_info) m_viscous::is2_viscous |
| type(int_bounds_info) m_viscous::is3_viscous |
| type(int_bounds_info) m_viscous::iv |
| real(wp), dimension(:, :), allocatable m_viscous::res_viscous |