Compute the time step to be used at each cell center
Local time stepping can be used to get the solution advance towards steady state faster. If only the steady state solution is required, i.e., transients are irrelevant, use local time stepping.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in), | dimension(-2:dims%imx+2, -2:dims%jmx+2, -2:dims%kmx+2, 1:dims%n_var), target | :: | qp | Store primitive variable at cell center |
|
| real(kind=wp), | intent(inout), | dimension(1:dims%imx-1, 1:dims%jmx-1, 1:dims%kmx-1) | :: | delta_t | Local time increment value at each cell center |
|
| type(celltype), | intent(in), | dimension(-2:dims%imx+2,-2:dims%jmx+2,-2:dims%kmx+2) | :: | cells | Input cell quantities: volume |
|
| type(facetype), | intent(in), | dimension(-2:dims%imx+3,-2:dims%jmx+2,-2:dims%kmx+2) | :: | Ifaces | Store face quantites for I faces |
|
| type(facetype), | intent(in), | dimension(-2:dims%imx+2,-2:dims%jmx+3,-2:dims%kmx+2) | :: | Jfaces | Store face quantites for J faces |
|
| type(facetype), | intent(in), | dimension(-2:dims%imx+2,-2:dims%jmx+2,-2:dims%kmx+3) | :: | Kfaces | Store face quantites for K faces |
|
| real(kind=wp), | intent(in) | :: | CFL | CFL number |
||
| type(schemetype), | intent(in) | :: | scheme | finite-volume Schemes |
||
| type(flowtype), | intent(in) | :: | flow | Information about fluid flow: freestream-speed, ref-viscosity,etc. |
||
| type(extent), | intent(in) | :: | dims | Extent of the domain:imx,jmx,kmx |
subroutine compute_local_time_step(qp, delta_t, cells, Ifaces, Jfaces, Kfaces, CFL, scheme, flow, dims)
!< Compute the time step to be used at each cell center
!<
!< Local time stepping can be used to get the solution
!< advance towards steady state faster. If only the steady
!< state solution is required, i.e., transients are
!< irrelevant, use local time stepping.
!-----------------------------------------------------------
implicit none
real(wp), intent(in) :: CFL
!< CFL number
type(schemetype), intent(in) :: scheme
!< finite-volume Schemes
type(flowtype), intent(in) :: flow
!< Information about fluid flow: freestream-speed, ref-viscosity,etc.
type(extent), intent(in) :: dims
!< Extent of the domain:imx,jmx,kmx
real(wp), dimension(-2:dims%imx+2, -2:dims%jmx+2, -2:dims%kmx+2, 1:dims%n_var), intent(in), target :: qp
!< Store primitive variable at cell center
real(wp) , dimension(1:dims%imx-1, 1:dims%jmx-1, 1:dims%kmx-1), intent(inout) :: delta_t
!< Local time increment value at each cell center
type(celltype), dimension(-2:dims%imx+2,-2:dims%jmx+2,-2:dims%kmx+2), intent(in) :: cells
!< Input cell quantities: volume
type(facetype), dimension(-2:dims%imx+3,-2:dims%jmx+2,-2:dims%kmx+2), intent(in) :: Ifaces
!< Store face quantites for I faces
type(facetype), dimension(-2:dims%imx+2,-2:dims%jmx+3,-2:dims%kmx+2), intent(in) :: Jfaces
!< Store face quantites for J faces
type(facetype), dimension(-2:dims%imx+2,-2:dims%jmx+2,-2:dims%kmx+3), intent(in) :: Kfaces
!< Store face quantites for K faces
real(wp) :: lmx1, lmx2, lmx3, lmx4, lmx5, lmx6, lmxsum
real(wp) :: x_sound_speed_avg, y_sound_speed_avg, z_sound_speed_avg
integer :: i, j, k
DebugCall('compute_local_time_step')
do k = 1, kmx - 1
do j = 1, jmx - 1
do i = 1, imx - 1
! For orientation, refer to the report. The standard i,j,k
! direction are marked. All orientation notations are w.r.t
! to the perspective shown in the image.
! Faces with lower index
x_sound_speed_avg = 0.5 * (sqrt(flow%gm * x_qp_left(i, j, k, 5) / &
x_qp_left(i, j, k, 1)) + &
sqrt(flow%gm * x_qp_right(i, j, k, 5) / &
x_qp_right(i, j, k, 1)) )
y_sound_speed_avg = 0.5 * (sqrt(flow%gm * y_qp_left(i, j, k, 5) / &
y_qp_left(i, j, k, 1)) + &
sqrt(flow%gm * y_qp_right(i, j, k, 5) / &
y_qp_right(i, j, k, 1)) )
z_sound_speed_avg = 0.5 * (sqrt(flow%gm * z_qp_left(i, j, k, 5) / &
z_qp_left(i, j, k, 1)) + &
sqrt(flow%gm * z_qp_right(i, j, k, 5) / &
z_qp_right(i, j, k, 1)) )
! For left face: i.e., lower index face along xi direction
lmx1 = abs( &
(qp(i, j, k,2) * Ifaces(i, j, k)%nx) + &
(qp(i, j, k,3) * Ifaces(i, j, k)%ny) + &
(qp(i, j, k,4) * Ifaces(i, j, k)%nz)) + &
x_sound_speed_avg
! For front face, i.e., lower index face along eta direction
lmx2 = abs( &
(qp(i, j, k,2) * Jfaces(i, j, k)%nx) + &
(qp(i, j, k,3) * Jfaces(i, j, k)%ny) + &
(qp(i, j, k,4) * Jfaces(i, j, k)%nz)) + &
y_sound_speed_avg
! For bottom face, i.e., lower index face along zeta direction
lmx3 = abs( &
(qp(i, j, k,2) * Kfaces(i, j, k)%nx) + &
(qp(i, j, k,3) * Kfaces(i, j, k)%ny) + &
(qp(i, j, k,4) * Kfaces(i, j, k)%nz)) + &
z_sound_speed_avg
! Faces with higher index
x_sound_speed_avg = 0.5 * (sqrt(flow%gm * x_qp_left(i+1,j,k,5) / x_qp_left(i+1,j,k,1)) + &
sqrt(flow%gm * x_qp_right(i+1,j,k,5) / x_qp_right(i+1,j,k,1)) )
y_sound_speed_avg = 0.5 * (sqrt(flow%gm * y_qp_left(i,j+1,k,5) / y_qp_left(i,j+1,k,1)) + &
sqrt(flow%gm * y_qp_right(i,j+1,k,5) / y_qp_right(i,j+1,k,1)) )
z_sound_speed_avg = 0.5 * (sqrt(flow%gm * z_qp_left(i,j,k+1,5) / z_qp_left(i,j,k+1,1)) + &
sqrt(flow%gm * z_qp_right(i,j,k+1,5) / z_qp_right(i,j,k+1,1)) )
! For right face, i.e., higher index face along xi direction
lmx4 = abs( &
(qp(i+1, j, k,2) * Ifaces(i+1, j, k)%nx) + & !x_speed*xnx
(qp(i+1, j, k,3) * Ifaces(i+1, j, k)%ny) + & !y_speed*xny
(qp(i+1, j, k,4) * Ifaces(i+1, j, k)%nz)) + & !z_speed*xnz
x_sound_speed_avg
! For back face, i.e., higher index face along eta direction
lmx5 = abs( &
(qp(i, j+1, k,2) * Jfaces(i, j+1, k)%nx) + &
(qp(i, j+1, k,3) * Jfaces(i, j+1, k)%ny) + &
(qp(i, j+1, k,4) * Jfaces(i, j+1, k)%nz)) + &
y_sound_speed_avg
! For top face, i.e., higher index face along zeta direction
lmx6 = abs( &
(qp(i, j, k+1,2) * Kfaces(i, j, k+1)%nx) + &
(qp(i, j, k+1,3) * Kfaces(i, j, k+1)%ny) + &
(qp(i, j, k+1,4) * Kfaces(i, j, k+1)%nz)) + &
z_sound_speed_avg
lmxsum = (Ifaces(i, j, k)%A * lmx1) + &
(Jfaces(i, j, k)%A * lmx2) + &
(Kfaces(i, j, k)%A * lmx3) + &
(Ifaces(i+1, j, k)%A * lmx4) + &
(Jfaces(i, j+1, k)%A * lmx5) + &
(Kfaces(i, j, k+1)%A * lmx6)
delta_t(i, j, k) = 1. / lmxsum
delta_t(i, j, k) = delta_t(i, j, k) * cells(i, j, k)%volume * CFL
end do
end do
end do
if(flow%mu_ref/=0.0) then
call add_viscous_time(qp, delta_t, cells, Ifaces, Jfaces, Kfaces, CFL, flow, dims)
end if
if(flow%mu_ref/=0 .and. trim(scheme%turbulence)/='none')then
call add_turbulent_time(qp, delta_t, cells, Ifaces, Jfaces, Kfaces, CFL, flow, dims)
end if
end subroutine compute_local_time_step