Addition to local time step due to turbulence
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(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 add_turbulent_time(qp,delta_t,cells,Ifaces,Jfaces,Kfaces,CFL,flow,dims)
!< Addition to local time step due to turbulence
implicit none
real(wp), intent(in) :: CFL
!< CFL number
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
integer :: i, j, k
DebugCall('add_viscous_time_step')
do k = 1, kmx - 1
do j = 1, jmx - 1
do i = 1, imx - 1
! Faces with lower index
! For left face: i.e., lower index face along xi direction
lmx1 = mu_t(i,j,k)/(qp(i,j,k,1)*abs( &
((cells(i-1,j,k)%centerx - cells(i,j,k)%centerx) * Ifaces(i, j, k)%nx) + &
((cells(i-1,j,k)%centery - cells(i,j,k)%centery) * Ifaces(i, j, k)%ny) + &
((cells(i-1,j,k)%centerz - cells(i,j,k)%centerz) * Ifaces(i, j, k)%nz)))
! For front face, i.e., lower index face along eta direction
lmx2 = mu_t(i,j,k)/(qp(i,j,k,1)*abs( &
((cells(i,j-1,k)%centerx - cells(i,j,k)%centerx) * Jfaces(i, j, k)%nx) + &
((cells(i,j-1,k)%centery - cells(i,j,k)%centery) * Jfaces(i, j, k)%ny) + &
((cells(i,j-1,k)%centerz - cells(i,j,k)%centerz) * Jfaces(i, j, k)%nz)))
! For bottom face, i.e., lower index face along zeta direction
lmx3 = mu_t(i,j,k)/(qp(i,j,k,1)*abs( &
((cells(i,j,k-1)%centerx - cells(i,j,k)%centerx) * Kfaces(i, j, k)%nx) + &
((cells(i,j,k-1)%centery - cells(i,j,k)%centery) * Kfaces(i, j, k)%ny) + &
((cells(i,j,k-1)%centerz - cells(i,j,k)%centerz) * Kfaces(i, j, k)%nz)))
! For right face, i.e., higher index face along xi direction
lmx4 = mu_t(i+1,j,k)/(qp(i+1,j,k,1)*abs( &
((cells(i,j,k)%centerx - cells(i+1,j,k)%centerx) * Ifaces(i+1, j, k)%nx) + &
((cells(i,j,k)%centery - cells(i+1,j,k)%centery) * Ifaces(i+1, j, k)%ny) + &
((cells(i,j,k)%centerz - cells(i+1,j,k)%centerz) * Ifaces(i+1, j, k)%nz)))
! For back face, i.e., higher index face along eta direction
lmx5 = mu_t(i,j+1,k)/(qp(i,j+1,k,1)*abs( &
((cells(i,j,k)%centerx - cells(i,j+1,k)%centerx) * Jfaces(i, j+1, k)%nx) + &
((cells(i,j,k)%centery - cells(i,j+1,k)%centery) * Jfaces(i, j+1, k)%ny) + &
((cells(i,j,k)%centerz - cells(i,j+1,k)%centerz) * Jfaces(i, j+1, k)%nz)))
! For top face, i.e., higher index face along zeta direction
lmx6 = mu_t(i,j,k+1)/(qp(i,j,k+1,1)*abs( &
((cells(i,j,k)%centerx - cells(i,j,k+1)%centerx) * Kfaces(i, j, k+1)%nx) + &
((cells(i,j,k)%centery - cells(i,j,k+1)%centery) * Kfaces(i, j, k+1)%ny) + &
((cells(i,j,k)%centerz - cells(i,j,k+1)%centerz) * Kfaces(i, j, k+1)%nz)))
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)
lmxsum = flow%gm*lmxsum/flow%tPr
lmxsum = 2./(lmxsum + (2.*CFL*cells(i,j,k)%volume/delta_t(i,j,k)))
delta_t(i, j, k) = CFL*( lmxsum * cells(i, j, k)%volume)
end do
end do
end do
end subroutine add_turbulent_time