scale_atmos_phy_tb_hybrid.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
module scale_atmos_phy_tb_hybrid
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_stdio
  use scale_prof
  use scale_grid_index
  use scale_tracer

#ifdef DEBUG
  use scale_debug, only: &
     check
  use scale_const, only: &
     undef  => const_undef, &
     iundef => const_undef2
#endif
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: atmos_phy_tb_hybrid_setup
  public :: atmos_phy_tb_hybrid

  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  abstract interface
     subroutine tb( &
       qflx_sgs_momz, qflx_sgs_momx, qflx_sgs_momy, & ! (out)
       qflx_sgs_rhot, qflx_sgs_rhoq,                & ! (out)
       tke,                                         & ! (inout)
       tke_t, nu_c, ri, pr, n2,                     & ! (out) diagnostic variables
       momz, momx, momy, rhot, dens, qtrc,          & ! (in)
       sflx_mw, sflx_mu, sflx_mv, sflx_sh, sflx_qv, & ! (in)
       gsqrt, j13g, j23g, j33g, mapf, dt            ) ! (in)
       use scale_precision
       use scale_grid_index
       use scale_tracer
       implicit none
       real(RP), intent(out) :: qflx_sgs_momz(ka,ia,ja,3)
       real(RP), intent(out) :: qflx_sgs_momx(ka,ia,ja,3)
       real(RP), intent(out) :: qflx_sgs_momy(ka,ia,ja,3)
       real(RP), intent(out) :: qflx_sgs_rhot(ka,ia,ja,3)
       real(RP), intent(out) :: qflx_sgs_rhoq(ka,ia,ja,3,qa)

       real(RP), intent(inout) :: tke (ka,ia,ja) ! TKE
       real(RP), intent(out)   :: tke_t(ka,ia,ja) ! TKE
       real(RP), intent(out)   :: nu_c(ka,ia,ja) ! eddy viscosity (center)
       real(RP), intent(out)   :: ri  (ka,ia,ja) ! Richardson number
       real(RP), intent(out)   :: pr  (ka,ia,ja) ! Prantle number
       real(RP), intent(out)   :: n2  (ka,ia,ja) ! squared Brunt-Vaisala frequency

       real(RP), intent(in)  :: momz(ka,ia,ja)
       real(RP), intent(in)  :: momx(ka,ia,ja)
       real(RP), intent(in)  :: momy(ka,ia,ja)
       real(RP), intent(in)  :: rhot(ka,ia,ja)
       real(RP), intent(in)  :: dens(ka,ia,ja)
       real(RP), intent(in)  :: qtrc(ka,ia,ja,qa)

       real(RP), intent(in)  :: sflx_mw(ia,ja)
       real(RP), intent(in)  :: sflx_mu(ia,ja)
       real(RP), intent(in)  :: sflx_mv(ia,ja)
       real(RP), intent(in)  :: sflx_sh(ia,ja)
       real(RP), intent(in)  :: sflx_qv(ia,ja)

       real(RP), intent(in)  :: gsqrt   (ka,ia,ja,7)
       real(RP), intent(in)  :: j13g    (ka,ia,ja,7)
       real(RP), intent(in)  :: j23g    (ka,ia,ja,7)
       real(RP), intent(in)  :: j33g
       real(RP), intent(in)  :: mapf    (ia,ja,2,4)
       real(DP), intent(in)  :: dt
     end subroutine tb
  end interface
  procedure(tb), pointer :: sgs_tb => null()
  procedure(tb), pointer :: pbl_tb => null()
  real(RP), allocatable :: frac_sgs(:,:)
  real(RP), allocatable :: frac_pbl(:,:)
  real(RP), allocatable :: frac_sgs_tke(:,:)
  real(RP), allocatable :: frac_pbl_tke(:,:)
  !-----------------------------------------------------------------------------
contains
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_tb_hybrid_setup( &
       TB_TYPE,       &
       CDZ, CDX, CDY, &
       CZ             )
    use scale_process, only: &
       prc_mpistop
    use scale_atmos_phy_tb_smg, only: &
       atmos_phy_tb_smg_setup, &
       atmos_phy_tb_smg
    use scale_atmos_phy_tb_mynn, only: &
       atmos_phy_tb_mynn_setup, &
       atmos_phy_tb_mynn
    implicit none

    character(len=*), intent(in) :: TB_TYPE

    real(RP), intent(in) :: CDZ(ka)
    real(RP), intent(in) :: CDX(ia)
    real(RP), intent(in) :: CDY(ja)
    real(RP), intent(in) :: CZ (ka,ia,ja)

    real(RP) :: ATMOS_PHY_TB_HYBRID_SGS_DX = 100.0_rp                      
    real(RP) :: ATMOS_PHY_TB_HYBRID_PBL_DX = 500.0_rp                      
    character(len=H_SHORT) :: ATMOS_PHY_TB_HYBRID_SGS_TYPE = 'SMAGORINSKY'
    character(len=H_SHORT) :: ATMOS_PHY_TB_HYBRID_PBL_TYPE = 'MYNN'

    character(len=H_SHORT) :: ATMOS_PHY_TB_HYBRID_TKE_TYPE = 'PBL'

    namelist / param_atmos_phy_tb_hybrid / &
         atmos_phy_tb_hybrid_sgs_dx, &
         atmos_phy_tb_hybrid_pbl_dx, &
         atmos_phy_tb_hybrid_sgs_type, &
         atmos_phy_tb_hybrid_pbl_type, &
         atmos_phy_tb_hybrid_tke_type


    real(RP) :: dxy

    integer :: i, j
    integer :: ierr
    !---------------------------------------------------------------------------

    if( io_l ) write(io_fid_log,*)
    if( io_l ) write(io_fid_log,*) '++++++ Module[TURBULENCE] / Categ[ATMOS PHYSICS] / Origin[SCALElib]'
    if( io_l ) write(io_fid_log,*) '+++ SGS-parameterization hybrid Model'

    if ( tb_type /= 'HYBRID' ) then
       write(*,*) 'xxx ATMOS_PHY_TB_TYPE is not HYBRID. Check!'
       call prc_mpistop
    endif


    !--- read namelist
    rewind(io_fid_conf)
    read(io_fid_conf,nml=param_atmos_phy_tb_hybrid,iostat=ierr)
    if( ierr < 0 ) then !--- missing
       if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
    elseif( ierr > 0 ) then !--- fatal error
       write(*,*) 'xxx Not appropriate names in namelist PARAM_ATMOS_PHY_TB_HYBRID. Check!'
       call prc_mpistop
    endif
    if( io_lnml ) write(io_fid_log,nml=param_atmos_phy_tb_hybrid)

    select case ( atmos_phy_tb_hybrid_sgs_type )
    case ('SMAGORINSKY')
       call atmos_phy_tb_smg_setup( &
            atmos_phy_tb_hybrid_sgs_type, &
            cdz, cdx, cdy, &
            cz             )
       sgs_tb => atmos_phy_tb_smg
    case default
       write(*,*) 'xxx ATMOS_PHY_TB_HYBRID_SGS_TYPE is invalid'
       call prc_mpistop
    end select

    select case ( atmos_phy_tb_hybrid_pbl_type )
    case ('MYNN')
       call atmos_phy_tb_mynn_setup( &
            atmos_phy_tb_hybrid_pbl_type, &
            cdz, cdx, cdy, &
            cz             )
       pbl_tb => atmos_phy_tb_mynn
    case default
       write(*,*) 'xxx ATMOS_PHY_TB_HYBRID_PBL_TYPE is invalid'
       call prc_mpistop
    end select

    allocate( frac_sgs(ia,ja) )
    allocate( frac_pbl(ia,ja) )
    allocate( frac_sgs_tke(ia,ja) )
    allocate( frac_pbl_tke(ia,ja) )

    do j = 1, ja
    do i = 1, ia
       dxy = sqrt( ( cdx(i)**2 + cdy(j)**2 )*0.5_rp )
       frac_pbl(i,j) = &
               min( 1.0_rp, &
               max( 0.0_rp, &
                    ( dxy - atmos_phy_tb_hybrid_sgs_dx ) &
                  / ( atmos_phy_tb_hybrid_pbl_dx - atmos_phy_tb_hybrid_sgs_dx ) ) )
       frac_sgs(i,j) = 1.0_rp - frac_pbl(i,j)
    end do
    end do

    select case ( atmos_phy_tb_hybrid_tke_type )
    case ('SGS')
       frac_pbl_tke(:,:) = 0.0_rp
       frac_sgs_tke(:,:) = 1.0_rp
    case ('PBL')
       frac_pbl_tke(:,:) = 1.0_rp
       frac_sgs_tke(:,:) = 0.0_rp
    case ('MIXED')
       frac_pbl_tke(:,:) = frac_pbl(:,:)
       frac_sgs_tke(:,:) = frac_sgs(:,:)
    case default
       write(*,*) 'xxx ATMOS_PHY_TB_HYBRID_TKE_TYPE is invalid'
       call prc_mpistop
    end select

    return
  end subroutine atmos_phy_tb_hybrid_setup

  !-----------------------------------------------------------------------------
  subroutine atmos_phy_tb_hybrid( &
       qflx_sgs_momz, qflx_sgs_momx, qflx_sgs_momy, & ! (out)
       qflx_sgs_rhot, qflx_sgs_rhoq,                & ! (out)
       tke,                                         & ! (inout)
       tke_t, nu, ri, pr, n2,                       & ! (out) diagnostic variables
       momz, momx, momy, rhot, dens, qtrc,          & ! (in)
       sflx_mw, sflx_mu, sflx_mv, sflx_sh, sflx_qv, & ! (in)
       gsqrt, j13g, j23g, j33g, mapf, dt            ) ! (in)
    use scale_precision
    use scale_grid_index
    use scale_tracer
    use scale_const, only: &
       grav => const_grav
    implicit none

    real(RP), intent(out) :: qflx_sgs_momz(ka,ia,ja,3)
    real(RP), intent(out) :: qflx_sgs_momx(ka,ia,ja,3)
    real(RP), intent(out) :: qflx_sgs_momy(ka,ia,ja,3)
    real(RP), intent(out) :: qflx_sgs_rhot(ka,ia,ja,3)
    real(RP), intent(out) :: qflx_sgs_rhoq(ka,ia,ja,3,qa)

    real(RP), intent(inout) :: tke (ka,ia,ja) ! TKE
    real(RP), intent(out) :: tke_t(ka,ia,ja) ! tendency TKE
    real(RP), intent(out) :: Nu(ka,ia,ja) ! eddy viscosity (center)
    real(RP), intent(out) :: Pr(ka,ia,ja) ! Plandtle number
    real(RP), intent(out) :: Ri(ka,ia,ja) ! Richardson number
    real(RP), intent(out) :: N2(ka,ia,ja) ! squared Brunt-Vaisala frequency

    real(RP), intent(in)  :: MOMZ(ka,ia,ja)
    real(RP), intent(in)  :: MOMX(ka,ia,ja)
    real(RP), intent(in)  :: MOMY(ka,ia,ja)
    real(RP), intent(in)  :: RHOT(ka,ia,ja)
    real(RP), intent(in)  :: DENS(ka,ia,ja)
    real(RP), intent(in)  :: QTRC(ka,ia,ja,qa)

    real(RP), intent(in)  :: SFLX_MW(ia,ja)
    real(RP), intent(in)  :: SFLX_MU(ia,ja)
    real(RP), intent(in)  :: SFLX_MV(ia,ja)
    real(RP), intent(in)  :: SFLX_SH(ia,ja)
    real(RP), intent(in)  :: SFLX_QV(ia,ja)

    real(RP), intent(in)  :: GSQRT   (ka,ia,ja,7)
    real(RP), intent(in)  :: J13G    (ka,ia,ja,7)
    real(RP), intent(in)  :: J23G    (ka,ia,ja,7)
    real(RP), intent(in)  :: J33G
    real(RP), intent(in)  :: MAPF    (ia,ja,2,4)
    real(DP), intent(in)  :: dt

    real(RP) :: w_qflx_sgs_momz(ka,ia,ja,3,2)
    real(RP) :: w_qflx_sgs_momx(ka,ia,ja,3,2)
    real(RP) :: w_qflx_sgs_momy(ka,ia,ja,3,2)
    real(RP) :: w_qflx_sgs_rhot(ka,ia,ja,3,2)
    real(RP) :: w_qflx_sgs_rhoq(ka,ia,ja,3,qa,2)

    real(RP) :: w_tke(ka,ia,ja,2)
    real(RP) :: w_tke_t(ka,ia,ja,2)
    real(RP) :: w_Nu(ka,ia,ja,2)
    real(RP) :: w_Ri(ka,ia,ja,2)
    real(RP) :: w_Pr(ka,ia,ja,2)
    real(RP) :: w_N2(ka,ia,ja,2)

    integer :: k, i, j, iq

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       w_tke(k,i,j,1) = tke(k,i,j)
       w_tke(k,i,j,2) = tke(k,i,j)
    end do
    end do
    end do

    call sgs_tb( &
         w_qflx_sgs_momz(:,:,:,:,1), w_qflx_sgs_momx(:,:,:,:,1), & ! (out)
         w_qflx_sgs_momy(:,:,:,:,1), w_qflx_sgs_rhot(:,:,:,:,1), & ! (out
         w_qflx_sgs_rhoq(:,:,:,:,:,1),                           & ! (out)
         w_tke(:,:,:,1),                                         & ! (inout)
         w_tke_t(:,:,:,1),                                       & ! (out)
         w_nu(:,:,:,1), w_ri(:,:,:,1), w_pr(:,:,:,1), w_n2(:,:,:,1), & ! (out)
         momz, momx, momy, rhot, dens, qtrc,                     & ! (in)
         sflx_mw, sflx_mu, sflx_mv, sflx_sh, sflx_qv,            & ! (in)
         gsqrt, j13g, j23g, j33g, mapf, dt                       ) ! (in)

    call pbl_tb( &
         w_qflx_sgs_momz(:,:,:,:,2), w_qflx_sgs_momx(:,:,:,:,2), & ! (out)
         w_qflx_sgs_momy(:,:,:,:,2), w_qflx_sgs_rhot(:,:,:,:,2), & ! (out
         w_qflx_sgs_rhoq(:,:,:,:,:,2),                           & ! (out)
         w_tke(:,:,:,2),                                         & ! (inout)
         w_tke_t(:,:,:,2),                                       & ! (out)
         w_nu(:,:,:,2), w_ri(:,:,:,2), w_pr(:,:,:,2), w_n2(:,:,:,2), & ! (out)
         momz, momx, momy, rhot, dens, qtrc,                     & ! (in)
         sflx_mw, sflx_mu, sflx_mv, sflx_sh, sflx_qv,            & ! (in)
         gsqrt, j13g, j23g, j33g, mapf, dt                       ) ! (in)

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       qflx_sgs_momz(k,i,j,zdir) = w_qflx_sgs_momz(k,i,j,zdir,1) * frac_sgs(i,j) &
                                 + w_qflx_sgs_momz(k,i,j,zdir,2) * frac_pbl(i,j)
       qflx_sgs_momz(k,i,j,xdir) = w_qflx_sgs_momz(k,i,j,xdir,1)
       qflx_sgs_momz(k,i,j,ydir) = w_qflx_sgs_momz(k,i,j,ydir,1)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       qflx_sgs_momx(k,i,j,zdir) = w_qflx_sgs_momx(k,i,j,zdir,1) * frac_sgs(i,j) &
                                 + w_qflx_sgs_momx(k,i,j,zdir,2) * frac_pbl(i,j)
       qflx_sgs_momx(k,i,j,xdir) = w_qflx_sgs_momx(k,i,j,xdir,1)
       qflx_sgs_momx(k,i,j,ydir) = w_qflx_sgs_momx(k,i,j,ydir,1)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       qflx_sgs_momy(k,i,j,zdir) = w_qflx_sgs_momy(k,i,j,zdir,1) * frac_sgs(i,j) &
                                 + w_qflx_sgs_momy(k,i,j,zdir,2) * frac_pbl(i,j)
       qflx_sgs_momy(k,i,j,xdir) = w_qflx_sgs_momy(k,i,j,xdir,1)
       qflx_sgs_momy(k,i,j,ydir) = w_qflx_sgs_momy(k,i,j,ydir,1)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       qflx_sgs_rhot(k,i,j,zdir) = w_qflx_sgs_rhot(k,i,j,zdir,1) * frac_sgs(i,j) &
                                 + w_qflx_sgs_rhot(k,i,j,zdir,2) * frac_pbl(i,j)
       qflx_sgs_rhot(k,i,j,xdir) = w_qflx_sgs_rhot(k,i,j,xdir,1)
       qflx_sgs_rhot(k,i,j,ydir) = w_qflx_sgs_rhot(k,i,j,ydir,1)
    end do
    end do
    end do

    do iq = 1, qa
    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       qflx_sgs_rhoq(k,i,j,zdir,iq) = w_qflx_sgs_rhoq(k,i,j,zdir,iq,1) * frac_sgs(i,j) &
                                    + w_qflx_sgs_rhoq(k,i,j,zdir,iq,2) * frac_pbl(i,j)
       qflx_sgs_rhoq(k,i,j,xdir,iq) = w_qflx_sgs_rhoq(k,i,j,xdir,iq,1)
       qflx_sgs_rhoq(k,i,j,ydir,iq) = w_qflx_sgs_rhoq(k,i,j,ydir,iq,1)
    end do
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       tke(k,i,j) = w_tke(k,i,j,1) * frac_sgs_tke(i,j) &
                  + w_tke(k,i,j,2) * frac_pbl_tke(i,j)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       tke_t(k,i,j) = w_tke_t(k,i,j,1) * frac_sgs_tke(i,j) &
                    + w_tke_t(k,i,j,2) * frac_pbl_tke(i,j)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       nu(k,i,j) = w_nu(k,i,j,1) * frac_sgs(i,j) &
                 + w_nu(k,i,j,2) * frac_pbl(i,j)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       ri(k,i,j) = w_ri(k,i,j,1) * frac_sgs(i,j) &
                 + w_ri(k,i,j,2) * frac_pbl(i,j)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       pr(k,i,j) = w_pr(k,i,j,1) * frac_sgs(i,j) &
                 + w_pr(k,i,j,2) * frac_pbl(i,j)
    end do
    end do
    end do

    do j = 1, ja
    do i = 1, ia
    do k = ks, ke
       n2(k,i,j) = w_n2(k,i,j,1) * frac_sgs(i,j) &
                 + w_n2(k,i,j,2) * frac_pbl(i,j)
    end do
    end do
    end do

    return
  end subroutine atmos_phy_tb_hybrid

end module scale_atmos_phy_tb_hybrid