|
MFC
Exascale flow solver
|
Shared Riemann-solver module state and the per-sweep setup, state-buffer population, viscous source flux, and finalization helpers. More...
Functions/Subroutines | |
| subroutine | s_compute_viscous_source_flux (vell_vf, dvell_dx_vf, dvell_dy_vf, dvell_dz_vf, velr_vf, dvelr_dx_vf, dvelr_dy_vf, dvelr_dz_vf, flux_src_vf, q_prim_vf, norm_dir, ix, iy, iz) |
| Dispatch to the subroutines that are utilized to compute the viscous source fluxes for either Cartesian or cylindrical geometries. For more information please refer to: 1) s_compute_cartesian_viscous_source_flux 2) s_compute_cylindrical_viscous_source_flux. | |
| subroutine | s_populate_riemann_states_variables_buffers (ql_prim_rsx_vf, dql_prim_dx_vf, dql_prim_dy_vf, dql_prim_dz_vf, qr_prim_rsx_vf, dqr_prim_dx_vf, dqr_prim_dy_vf, dqr_prim_dz_vf, norm_dir, ix, iy, iz) |
| Populate the left and right Riemann state variable buffers based on boundary conditions. | |
| subroutine | s_initialize_riemann_solver (flux_src_vf, norm_dir) |
| Set up the chosen Riemann solver algorithm for the current direction. | |
| subroutine | s_compute_cylindrical_viscous_source_flux (vell_vf, dvell_dx_vf, dvell_dy_vf, dvell_dz_vf, velr_vf, dvelr_dx_vf, dvelr_dy_vf, dvelr_dz_vf, flux_src_vf, q_prim_vf, norm_dir, ix, iy, iz) |
| Compute cylindrical viscous source flux contributions for momentum and energy. | |
| subroutine | s_compute_cartesian_viscous_source_flux (dvell_dx_vf, dvell_dy_vf, dvell_dz_vf, dvelr_dx_vf, dvelr_dy_vf, dvelr_dz_vf, flux_src_vf, q_prim_vf, norm_dir) |
| Compute Cartesian viscous source flux contributions for momentum and energy. | |
| subroutine | s_calculate_shear_stress_tensor (vel_grad_avg, re_shear, divergence_v, tau_shear_out) |
| Compute shear stress tensor components. | |
| subroutine | s_calculate_bulk_stress_tensor (re_bulk, divergence_v, tau_bulk_out) |
| Compute bulk stress tensor components (diagonal only). | |
| subroutine | s_finalize_riemann_solver (flux_vf, flux_src_vf, flux_gsrc_vf, norm_dir) |
| Deallocation and/or disassociation procedures that are needed to finalize the selected Riemann problem solver. | |
Variables | |
| real(wp), dimension(:,:,:,:), allocatable | vel_src_rsx_vf |
| real(wp), dimension(:,:,:,:), allocatable | mom_sp_rsx_vf |
| real(wp), dimension(:,:,:,:), allocatable | re_avg_rsx_vf |
| real(wp), dimension(:), allocatable | gs_rs |
| real(wp), dimension(:,:), allocatable | res_gs |
| real(wp), dimension(:,:,:,:), allocatable | flux_rsx_vf |
| The cell-boundary values of the fluxes (src - source) that are computed through the chosen Riemann problem solver, and the direct evaluation of source terms, by using the left and right states given in qK_prim_rs_vf, dqK_prim_ds_vf where ds = dx, dy or dz. | |
| real(wp), dimension(:,:,:,:), allocatable | flux_src_rsx_vf |
| real(wp), dimension(:,:,:,:), allocatable | flux_gsrc_rsx_vf |
| The cell-boundary values of the geometrical source flux that are computed through the chosen Riemann problem solver by using the left and right states given in qK_prim_rs_vf. Currently 2D axisymmetric for inviscid only. | |
Indical bounds in the s1-, s2- and s3-directions | |
| type(int_bounds_info) | is1 |
| type(int_bounds_info) | is2 |
| type(int_bounds_info) | is3 |
| type(int_bounds_info) | isx |
| type(int_bounds_info) | isy |
| type(int_bounds_info) | isz |
Shared Riemann-solver module state and the per-sweep setup, state-buffer population, viscous source flux, and finalization helpers.
| subroutine m_riemann_state::s_calculate_bulk_stress_tensor | ( | real(wp), intent(in) | re_bulk, |
| real(wp), intent(in) | divergence_v, | ||
| real(wp), dimension(num_dims, num_dims), intent(out) | tau_bulk_out ) |
Compute bulk stress tensor components (diagonal only).
Definition at line 2490 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_calculate_shear_stress_tensor | ( | real(wp), dimension(num_dims, num_dims), intent(in) | vel_grad_avg, |
| real(wp), intent(in) | re_shear, | ||
| real(wp), intent(in) | divergence_v, | ||
| real(wp), dimension(num_dims, num_dims), intent(out) | tau_shear_out ) |
Compute shear stress tensor components.
Definition at line 2447 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_compute_cartesian_viscous_source_flux | ( | type(scalar_field), dimension(num_dims), intent(in) | dvell_dx_vf, |
| type(scalar_field), dimension(num_dims), intent(in) | dvell_dy_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvell_dz_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvelr_dx_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvelr_dy_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvelr_dz_vf, | ||
| type(scalar_field), dimension(sys_size), intent(inout) | flux_src_vf, | ||
| type(scalar_field), dimension(sys_size), intent(in) | q_prim_vf, | ||
| integer, intent(in) | norm_dir ) |
Compute Cartesian viscous source flux contributions for momentum and energy.
Definition at line 2249 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_compute_cylindrical_viscous_source_flux | ( | type(scalar_field), dimension(num_dims), intent(in) | vell_vf, |
| type(scalar_field), dimension(num_dims), intent(in) | dvell_dx_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvell_dy_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvell_dz_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | velr_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvelr_dx_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvelr_dy_vf, | ||
| type(scalar_field), dimension(num_dims), intent(in) | dvelr_dz_vf, | ||
| type(scalar_field), dimension(sys_size), intent(inout) | flux_src_vf, | ||
| type(scalar_field), dimension(sys_size), intent(in) | q_prim_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 ) |
Compute cylindrical viscous source flux contributions for momentum and energy.
Interface velocity ( \(v_1,v_2,v_3\)) (grid directions) for viscous work.
Shear stress vector ( \(\sigma_{N1}, \sigma_{N2}, \sigma_{N3}\)) on N-face (grid directions).
Definition at line 1986 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_compute_viscous_source_flux | ( | type(scalar_field), dimension(num_vels), intent(in) | vell_vf, |
| type(scalar_field), dimension(num_vels), intent(in) | dvell_dx_vf, | ||
| type(scalar_field), dimension(num_vels), intent(in) | dvell_dy_vf, | ||
| type(scalar_field), dimension(num_vels), intent(in) | dvell_dz_vf, | ||
| type(scalar_field), dimension(num_vels), intent(in) | velr_vf, | ||
| type(scalar_field), dimension(num_vels), intent(in) | dvelr_dx_vf, | ||
| type(scalar_field), dimension(num_vels), intent(in) | dvelr_dy_vf, | ||
| type(scalar_field), dimension(num_vels), intent(in) | dvelr_dz_vf, | ||
| type(scalar_field), dimension(sys_size), intent(inout) | flux_src_vf, | ||
| type(scalar_field), dimension(sys_size), intent(in) | q_prim_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 ) |
Dispatch to the subroutines that are utilized to compute the viscous source fluxes for either Cartesian or cylindrical geometries. For more information please refer to: 1) s_compute_cartesian_viscous_source_flux 2) s_compute_cylindrical_viscous_source_flux.
Definition at line 463 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_finalize_riemann_solver | ( | type(scalar_field), dimension(sys_size), intent(inout) | flux_vf, |
| type(scalar_field), dimension(sys_size), intent(inout) | flux_src_vf, | ||
| type(scalar_field), dimension(sys_size), intent(inout) | flux_gsrc_vf, | ||
| integer, intent(in) | norm_dir ) |
Deallocation and/or disassociation procedures that are needed to finalize the selected Riemann problem solver.
Definition at line 2526 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_initialize_riemann_solver | ( | type(scalar_field), dimension(sys_size), intent(inout) | flux_src_vf, |
| integer, intent(in) | norm_dir ) |
Set up the chosen Riemann solver algorithm for the current direction.
Definition at line 1567 of file m_riemann_state.fpp.f90.
| subroutine m_riemann_state::s_populate_riemann_states_variables_buffers | ( | real(wp), dimension(idwbuff(1)%beg:,idwbuff(2)%beg:,idwbuff(3)%beg:,1:), intent(inout) | ql_prim_rsx_vf, |
| type(scalar_field), dimension(:), intent(inout), allocatable | dql_prim_dx_vf, | ||
| type(scalar_field), dimension(:), intent(inout), allocatable | dql_prim_dy_vf, | ||
| type(scalar_field), dimension(:), intent(inout), allocatable | dql_prim_dz_vf, | ||
| real(wp), dimension(idwbuff(1)%beg:,idwbuff(2)%beg:,idwbuff(3)%beg:,1:), intent(inout) | qr_prim_rsx_vf, | ||
| type(scalar_field), dimension(:), intent(inout), allocatable | dqr_prim_dx_vf, | ||
| type(scalar_field), dimension(:), intent(inout), allocatable | dqr_prim_dy_vf, | ||
| type(scalar_field), dimension(:), intent(inout), allocatable | dqr_prim_dz_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 ) |
Populate the left and right Riemann state variable buffers based on boundary conditions.
Definition at line 486 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:,:,:), allocatable m_riemann_state::flux_gsrc_rsx_vf |
The cell-boundary values of the geometrical source flux that are computed through the chosen Riemann problem solver by using the left and right states given in qK_prim_rs_vf. Currently 2D axisymmetric for inviscid only.
Definition at line 359 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:,:,:), allocatable m_riemann_state::flux_rsx_vf |
The cell-boundary values of the fluxes (src - source) that are computed through the chosen Riemann problem solver, and the direct evaluation of source terms, by using the left and right states given in qK_prim_rs_vf, dqK_prim_ds_vf where ds = dx, dy or dz.
Definition at line 342 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:,:,:), allocatable m_riemann_state::flux_src_rsx_vf |
Definition at line 342 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:), allocatable m_riemann_state::gs_rs |
Definition at line 432 of file m_riemann_state.fpp.f90.
| type(int_bounds_info) m_riemann_state::is1 |
Definition at line 416 of file m_riemann_state.fpp.f90.
| type(int_bounds_info) m_riemann_state::is2 |
Definition at line 416 of file m_riemann_state.fpp.f90.
| type(int_bounds_info) m_riemann_state::is3 |
Definition at line 416 of file m_riemann_state.fpp.f90.
| type(int_bounds_info) m_riemann_state::isx |
Definition at line 417 of file m_riemann_state.fpp.f90.
| type(int_bounds_info) m_riemann_state::isy |
Definition at line 417 of file m_riemann_state.fpp.f90.
| type(int_bounds_info) m_riemann_state::isz |
Definition at line 417 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:,:,:), allocatable m_riemann_state::mom_sp_rsx_vf |
Definition at line 388 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:,:,:), allocatable m_riemann_state::re_avg_rsx_vf |
Definition at line 401 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:), allocatable m_riemann_state::res_gs |
Definition at line 445 of file m_riemann_state.fpp.f90.
| real(wp), dimension(:,:,:,:), allocatable m_riemann_state::vel_src_rsx_vf |
Definition at line 375 of file m_riemann_state.fpp.f90.