mod_atmos_phy_bl_driver.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
#include "scalelib.h"
module mod_atmos_phy_bl_driver
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_io
  use scale_prof
  use scale_atmos_grid_cartesc_index
  use scale_tracer
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: atmos_phy_bl_driver_tracer_setup
  public :: atmos_phy_bl_driver_setup
  public :: atmos_phy_bl_driver_finalize
  public :: atmos_phy_bl_driver_mkinit
  public :: atmos_phy_bl_driver_calc_tendency
 
  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  !-----------------------------------------------------------------------------
contains
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_bl_driver_tracer_setup
    use scale_prc, only: &
       prc_abort
    use scale_tracer, only: &
       tracer_regist
    use scale_atmos_phy_bl_mynn, only: &
       atmos_phy_bl_mynn_tracer_setup, &
       atmos_phy_bl_mynn_ntracer, &
       atmos_phy_bl_mynn_name, &
       atmos_phy_bl_mynn_desc, &
       atmos_phy_bl_mynn_units
    use scale_atmos_phy_bl_mynn_jmapplib, only: &
       atmos_phy_bl_mynn_jmapplib_ntracer, &
       atmos_phy_bl_mynn_jmapplib_name, &
       atmos_phy_bl_mynn_jmapplib_desc, &
       atmos_phy_bl_mynn_jmapplib_units
    use mod_atmos_admin, only: &
       atmos_phy_bl_type, &
       atmos_sw_phy_bl
    use mod_atmos_phy_bl_vars, only: &
       qs, qe
    implicit none
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("ATMOS_PHY_BL_driver_tracer_setup",*) 'Setup'
 
    if ( atmos_sw_phy_bl ) then
       select case ( atmos_phy_bl_type )
       case ( 'MYNN' )
          call atmos_phy_bl_mynn_tracer_setup
          call tracer_regist( &
               qs, &
               atmos_phy_bl_mynn_ntracer, &
               atmos_phy_bl_mynn_name,    &
               atmos_phy_bl_mynn_desc,    &
               atmos_phy_bl_mynn_units    )
          qe = qs + atmos_phy_bl_mynn_ntracer - 1
       case ( 'MYNN-JMAPPLIB' )
          call tracer_regist( &
               qs, &
               atmos_phy_bl_mynn_jmapplib_ntracer, &
               atmos_phy_bl_mynn_jmapplib_name,    &
               atmos_phy_bl_mynn_jmapplib_desc,    &
               atmos_phy_bl_mynn_jmapplib_units    )
          qe = qs + atmos_phy_bl_mynn_jmapplib_ntracer - 1
       case default
          log_error("ATMOS_PHY_BL_driver_tracer_setup",*) 'ATMOS_PHY_BL_TYPE is invalid: ', trim(atmos_phy_bl_type)
          call prc_abort
       end select
    end if
 
    return
  end subroutine atmos_phy_bl_driver_tracer_setup
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_bl_driver_setup
    use scale_atmos_phy_bl_mynn, only: &
       atmos_phy_bl_mynn_setup
    use scale_atmos_phy_bl_mynn_jmapplib, only: &
       atmos_phy_bl_mynn_jmapplib_setup
    use scale_atmos_grid_cartesc, only: &
       cz => atmos_grid_cartesc_cz, &
       dx
    use mod_atmos_admin, only: &
       atmos_phy_bl_type, &
       atmos_sw_phy_bl
    use mod_atmos_phy_bl_vars, only: &
       atmos_phy_bl_zi, &
       atmos_phy_bl_sflx_buoy
    use scale_bulkflux, only: &
       bulkflux_type
    use scale_time, only: &
       dt_bl => time_dtsec_atmos_phy_bl
    implicit none
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("ATMOS_PHY_BL_driver_setup",*) 'Setup'
 
    if ( atmos_sw_phy_bl ) then
       select case ( atmos_phy_bl_type )
       case ( 'MYNN' )
          call atmos_phy_bl_mynn_setup( &
               bulkflux_type, & ! (in)
               dx = dx        ) ! (in)
       case ( 'MYNN-JMAPPLIB' )
          call atmos_phy_bl_mynn_jmapplib_setup( &
               ka, ks, ke, &
               cz(:), &
               dt_bl )
       end select
    else
       log_info("ATMOS_PHY_BL_driver_setup",*) 'this component is never called.'
       !$acc kernels
       atmos_phy_bl_zi(:,:) = 0.0_rp
       atmos_phy_bl_sflx_buoy(:,:) = 0.0_rp
       !$acc end kernels
    endif
 
    return
  end subroutine atmos_phy_bl_driver_setup
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_bl_driver_finalize
    use scale_prc, only: &
       prc_abort
    use scale_atmos_phy_bl_mynn, only: &
       atmos_phy_bl_mynn_finalize
    use scale_atmos_phy_bl_mynn_jmapplib, only: &
       atmos_phy_bl_mynn_jmapplib_finalize
    use mod_atmos_admin, only: &
       atmos_phy_bl_type, &
       atmos_sw_phy_bl
    implicit none
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("ATMOS_PHY_BL_driver_finalize",*) 'Finalize'
 
    if ( atmos_sw_phy_bl ) then
       select case ( atmos_phy_bl_type )
       case ( 'MYNN' )
          call atmos_phy_bl_mynn_finalize
       case ( 'MYNN-JMAPPLIB' )
          call atmos_phy_bl_mynn_jmapplib_finalize
       case default
          log_error("ATMOS_PHY_BL_driver_finalize",*) 'ATMOS_PHY_BL_TYPE is invalid: ', trim(atmos_phy_bl_type)
          call prc_abort
       end select
    end if
 
    return
  end subroutine atmos_phy_bl_driver_finalize
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_bl_driver_mkinit( TKE_CONST )
    use scale_const, only: &
       undef => const_undef
    use scale_time, only: &
       dt_bl => time_dtsec_atmos_phy_bl
    use scale_atmos_phy_bl_mynn, only: &
       atmos_phy_bl_mynn_mkinit
    use scale_atmos_grid_cartesc_real, only: &
       z1 => atmos_grid_cartesc_real_z1, &
       cz => atmos_grid_cartesc_real_cz, &
       fz => atmos_grid_cartesc_real_fz, &
       f2h => atmos_grid_cartesc_real_f2h
    use scale_bulkflux, only: &
       bulkflux_type
    use scale_landuse, only: &
       frac_land => landuse_frac_land
    use scale_topography, only: &
       tansl_x => topography_tansl_x, &
       tansl_y => topography_tansl_y
    use scale_atmos_bottom, only: &
       bottom_estimate => atmos_bottom_estimate
    use scale_atmos_phy_sf_bulk, only: &
       atmos_phy_sf_bulk_flux
    use scale_comm_cartesc, only: &
       comm_vars8, &
       comm_wait
    use mod_atmos_phy_bl_vars, only: &
       qs, qe, &
       zi => atmos_phy_bl_zi
    use mod_atmos_phy_sf_vars, only: &
       sfc_temp => atmos_phy_sf_sfc_temp,  &
       sfc_dens => atmos_phy_sf_sfc_dens,  &
       sfc_pres => atmos_phy_sf_sfc_pres,  &
       sfc_z0m   => atmos_phy_sf_sfc_z0m,   &
       sfc_z0h   => atmos_phy_sf_sfc_z0h,   &
       sfc_z0e   => atmos_phy_sf_sfc_z0e,   &
       sflx_mu  => atmos_phy_sf_sflx_mu,   &
       sflx_mv  => atmos_phy_sf_sflx_mv,   &
       sflx_mw  => atmos_phy_sf_sflx_mv,   &
       sflx_sh  => atmos_phy_sf_sflx_sh,   &
       sflx_lh  => atmos_phy_sf_sflx_lh,   &
       sflx_qv  => atmos_phy_sf_sflx_qv,   &
       ustar    => atmos_phy_sf_ustar,     &
       tstar    => atmos_phy_sf_tstar,     &
       qstar    => atmos_phy_sf_qstar,     &
       wstar    => atmos_phy_sf_wstar,     &
       rlmo     => atmos_phy_sf_rlmo,      &
       u10       => atmos_phy_sf_u10,      &
       v10       => atmos_phy_sf_v10,      &
       t2        => atmos_phy_sf_t2,       &
       q2        => atmos_phy_sf_q2
    use mod_atmos_admin, only: &
       atmos_phy_bl_type, &
       atmos_sw_phy_bl
    use mod_atmos_vars, only: &
       dens,  &
       momz,  &
       momx,  &
       momy,  &
       rhot,  &
       qtrc,  &
       u,     &
       v,     &
       w,     &
       temp,  &
       pott,  &
       pres,  &
       exner, &
       qdry,  &
       qv,    &
       qc,    &
       qi,    &
       atmos_vars_calc_diagnostics, &
       atmos_vars_get_diagnostic
    implicit none
    real(rp), intent(in) :: tke_const
 
    real(rp) :: n2  (ka,ia,ja)
    real(rp) :: potl(ka,ia,ja)
    real(rp) :: potv(ka,ia,ja)
 
    real(rp) :: atm_temp(ia,ja)
    real(rp) :: atm_pres(ia,ja)
    real(rp) :: atm_u   (ia,ja)
    real(rp) :: atm_v   (ia,ja)
    real(rp) :: atm_w   (ia,ja)
    real(rp) :: atm_qv  (ia,ja)
 
    real(rp) :: qw(ka,ia,ja)
 
    integer  :: k, i, j, iq
    !---------------------------------------------------------------------------
 
    if ( .not. atmos_sw_phy_bl ) return
 
    select case ( atmos_phy_bl_type )
    case ( 'MYNN' )
       if ( qtrc(ks,is,je,qs) == undef ) then
 
          !$acc data create(N2,POTL,POTV,ATM_TEMP,ATM_PRES,ATM_U,ATM_V,ATM_W,ATM_QV,QW)
 
          call comm_vars8(dens, 1)
          call comm_vars8(momz, 2)
          call comm_vars8(momx, 3)
          call comm_vars8(momy, 4)
          call comm_vars8(rhot, 5)
          call comm_vars8(qv, 6)
          call comm_vars8(qc, 7)
          call comm_vars8(qi, 8)
          call comm_wait(dens, 1)
          call comm_wait(momz, 2)
          call comm_wait(momx, 3)
          call comm_wait(momy, 4)
          call comm_wait(rhot, 5)
          call comm_wait(qv, 6)
          call comm_wait(qc, 7)
          call comm_wait(qi, 8)
 
          call atmos_vars_calc_diagnostics
          call atmos_vars_get_diagnostic( "N2",   n2   )
          call atmos_vars_get_diagnostic( "POTL", potl )
          call atmos_vars_get_diagnostic( "POTV", potv )
 
          call bottom_estimate( ka, ks, ke, ia, isb, ieb, ja, jsb, jeb, &
                                dens(:,:,:), pres(:,:,:), qv(:,:,:), & ! [IN]
                                sfc_temp(:,:),                       & ! [IN]
                                fz(:,:,:),                           & ! [IN]
                                sfc_dens(:,:), sfc_pres(:,:)         ) ! [OUT]
 
          !$omp parallel do
          !$acc kernels
          do j = jsb, jeb
          do i = isb, ieb
             atm_temp(i,j) = temp(ks,i,j)
             atm_pres(i,j) = pres(ks,i,j)
             atm_u(i,j) = u(ks,i,j)
             atm_v(i,j) = v(ks,i,j)
             atm_w(i,j) = atm_u(i,j) * tansl_x(i,j) + atm_v(i,j) * tansl_y(i,j)
             atm_qv(i,j) = qv(ks,i,j)
          enddo
          enddo
          !$acc end kernels
 
          call atmos_phy_sf_bulk_flux( ia, isb, ieb, ja, jsb, jeb, &
                                       atm_w(:,:), atm_u(:,:), atm_v(:,:),          & ! [IN]
                                       atm_temp(:,:), atm_pres(:,:), atm_qv(:,:),   & ! [IN]
                                       sfc_dens(:,:), sfc_temp(:,:), sfc_pres(:,:), & ! [IN]
                                       sfc_z0m(:,:), sfc_z0h(:,:), sfc_z0e(:,:),    & ! [IN]
                                       zi(:,:), z1(:,:),                            & ! [IN]
                                       sflx_mw(:,:), sflx_mu(:,:), sflx_mv(:,:),    & ! [OUT]
                                       sflx_sh(:,:), sflx_lh(:,:), sflx_qv(:,:),    & ! [OUT]
                                       ustar(:,:), tstar(:,:), qstar(:,:),          & ! [OUT]
                                       wstar(:,:),                                  & ! [OUT]
                                       rlmo(:,:),                                   & ! [OUT]
                                       u10(:,:), v10(:,:), t2(:,:), q2(:,:)         ) ! [OUT]
 
          !$acc kernels
          do j = jsb, jeb
          do i = isb, ieb
          do k = ks, ke
             qw(k,i,j) = qv(k,i,j) + qc(k,i,j) + qi(k,i,j)
          end do
          end do
          end do
          !$acc end kernels
 
          call atmos_phy_bl_mynn_mkinit( &
               ka, ks, ke, ia, isb, ieb, ja, jsb, jeb, &
               qtrc(:,:,:,qs:qe),                                      & ! (out)
               dens(:,:,:), u(:,:,:), v(:,:,:), w(:,:,:), pott(:,:,:), & ! (in)
               pres(:,:,:), exner(:,:,:), n2(:,:,:),                   & ! (in)
               qdry(:,:,:), qv(:,:,:), qw(:,:,:),                      & ! (in)
               potl(:,:,:), potv(:,:,:),                               & ! (in)
               sfc_dens(:,:),                                          & ! (in)
               sflx_mu(:,:), sflx_mv(:,:), sflx_sh(:,:), sflx_qv(:,:), & ! (in)
               ustar(:,:), tstar(:,:), qstar(:,:), rlmo(:,:),          & ! (in)
               frac_land(:,:),                                         & ! (in)
               cz(:,:,:), fz(:,:,:), f2h(:,:,:,:),                     & ! (in)
               bulkflux_type                                           ) ! (in)
 
          !$acc end data
 
       end if
 
    case ( 'MYNN-JMAPPLIB' )
       !$acc kernels
       !$acc loop independent collapse(3)
       do j = jsb, jeb
       do i = isb, ieb
       do k = ks, ke
          if ( qtrc(k,i,j,qs) == undef ) then
             qtrc(k,i,j,qs) = tke_const
          end if
       end do
       end do
       end do
       !$acc end kernels
       !$acc kernels
       !$acc loop independent collapse(4)
       do iq = qs+1, qe
       do j = jsb, jeb
       do i = isb, ieb
       do k = ks, ke
          if ( qtrc(k,i,j,iq) == undef ) then
             qtrc(k,i,j,iq) = 0.0_rp
          end if
       end do
       end do
       end do
       end do
       !$acc end kernels
    end select
 
    return
  end subroutine atmos_phy_bl_driver_mkinit
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_bl_driver_calc_tendency( update_flag )
    use scale_statistics, only: &
       statistics_checktotal, &
       statistics_total
    use scale_file_history, only: &
       file_history_in
    use scale_time, only: &
       dt_bl => time_dtsec_atmos_phy_bl
    use scale_atmos_phy_bl_common, only: &
       atmos_phy_bl_tendency_tracer
    use scale_atmos_phy_bl_mynn, only: &
       atmos_phy_bl_mynn_tendency
    use scale_atmos_phy_bl_mynn_jmapplib, only: &
       atmos_phy_bl_mynn_jmapplib_tendency
    use scale_atmos_grid_cartesc_real, only: &
       cz => atmos_grid_cartesc_real_cz, &
       fz => atmos_grid_cartesc_real_fz, &
       f2h => atmos_grid_cartesc_real_f2h, &
       atmos_grid_cartesc_real_vol, &
       atmos_grid_cartesc_real_totvol
    use scale_atmos_hydrometeor, only: &
       i_qv
    use scale_bulkflux, only: &
       bulkflux_type
    use scale_landuse, only: &
       frac_land => landuse_frac_land
    use mod_atmos_admin, only: &
       atmos_phy_bl_type
    use mod_atmos_vars, only: &
       dens => dens_av, &
       qtrc => qtrc_av, &
       u,     &
       v,     &
       w,     &
       pott,  &
       pres,  &
       exner, &
       qdry,  &
       qv,    &
       qc,    &
       qi,    &
       rhou_t => rhou_tp, &
       rhov_t => rhov_tp, &
       rhot_t => rhot_tp, &
       rhoq_t => rhoq_tp, &
       atmos_vars_get_diagnostic
    use mod_atmos_phy_bl_vars, only: &
       qs, qe, &
       atmos_phy_bl_mix_tracers, &
       rhou_t_bl => atmos_phy_bl_rhou_t,    &
       rhov_t_bl => atmos_phy_bl_rhov_t,    &
       rhot_t_bl => atmos_phy_bl_rhot_t,    &
       rhoq_t_bl => atmos_phy_bl_rhoq_t,    &
       zi        => atmos_phy_bl_zi,        &
       sflx_buoy => atmos_phy_bl_sflx_buoy, &
       ql        => atmos_phy_bl_ql,        &
       cldfrac   => atmos_phy_bl_cldfrac
    use mod_atmos_phy_sf_vars, only: &
       sfc_dens => atmos_phy_sf_sfc_dens,  &
       sfc_pres => atmos_phy_sf_sfc_pres,  &
       sflx_mu  => atmos_phy_sf_sflx_mu,   &
       sflx_mv  => atmos_phy_sf_sflx_mv,   &
       sflx_sh  => atmos_phy_sf_sflx_sh,   &
       sflx_q   => atmos_phy_sf_sflx_qtrc, &
       sflx_qv  => atmos_phy_sf_sflx_qv,   &
       ustar    => atmos_phy_sf_ustar,     &
       tstar    => atmos_phy_sf_tstar,     &
       qstar    => atmos_phy_sf_qstar,     &
       rlmo     => atmos_phy_sf_rlmo
    implicit none
 
    logical, intent(in) :: update_flag
 
    real(rp) :: nu(ka,ia,ja)
    real(rp) :: kh(ka,ia,ja)
 
    real(rp) :: qw(ka,ia,ja)
 
    real(rp) :: n2  (ka,ia,ja)
    real(rp) :: potl(ka,ia,ja)
    real(rp) :: potv(ka,ia,ja)
 
    real(rp), pointer :: rhoqv_t(:,:,:)
 
    integer  :: k, i, j, iq
    !---------------------------------------------------------------------------
 
    if ( update_flag ) then
 
       !$acc data create(Nu,Kh,QW,N2,POTL,POTV)
 
       !$acc kernels
       rhoq_t_bl(:,:,:,:) = 0.0_rp
       !$acc end kernels
 
       select case ( atmos_phy_bl_type )
       case ( 'MYNN' )
          call atmos_vars_get_diagnostic( "N2",   n2   )
          call atmos_vars_get_diagnostic( "POTL", potl )
          call atmos_vars_get_diagnostic( "POTV", potv )
          !$acc kernels
          do j = jsb, jeb
          do i = isb, ieb
          do k = ks, ke
             qw(k,i,j) = qv(k,i,j) + qc(k,i,j) + qi(k,i,j)
          end do
          end do
          end do
          !$acc end kernels
          if ( i_qv > 0 ) then
             rhoqv_t => rhoq_t_bl(:,:,:,i_qv)
          else
             allocate( rhoqv_t(ka,ia,ja) )
             !$acc enter data create(RHOQV_t)
          end if
          call atmos_phy_bl_mynn_tendency( &
               ka, ks, ke, ia, is, ie, ja, js, je, &
               dens(:,:,:), u(:,:,:), v(:,:,:), w(:,:,:),              & ! (in)
               pott(:,:,:), qtrc(:,:,:,qs:qe),                         & ! (in)
               pres(:,:,:), exner(:,:,:), n2(:,:,:),                   & ! (in)
               qdry(:,:,:), qv(:,:,:), qw(:,:,:),                      & ! (in)
               potl(:,:,:), potv(:,:,:),                               & ! (in)
               sfc_dens(:,:),                                          & ! (in)
               sflx_mu(:,:), sflx_mv(:,:), sflx_sh(:,:), sflx_qv(:,:), & ! (in)
               ustar(:,:), tstar(:,:), qstar(:,:), rlmo(:,:),          & ! (in)
               frac_land(:,:),                                         & ! (in)
               cz(:,:,:), fz(:,:,:), f2h(:,:,:,:), dt_bl,              & ! (in)
               bulkflux_type,                                          & ! (in)
               rhou_t_bl(:,:,:), rhov_t_bl(:,:,:), rhot_t_bl(:,:,:),   & ! (out)
               rhoqv_t(:,:,:), rhoq_t_bl(:,:,:,qs:qe),                 & ! (out)
               nu(:,:,:), kh(:,:,:),                                   & ! (out)
               ql(:,:,:), cldfrac(:,:,:), zi(:,:), sflx_buoy(:,:)      ) ! (out)
          if ( i_qv <= 0 ) then
             !$acc exit data delete(RHOQV_t)
             deallocate( rhoqv_t )
          end if
       case ( 'MYNN-JMAPPLIB' )
          if ( i_qv > 0 ) then
             rhoqv_t => rhoq_t_bl(:,:,:,i_qv)
          else
             allocate( rhoqv_t(ka,ia,ja) )
             !$acc enter data create(RHOQV_t)
          end if
          !$acc update host(DENS,U,V,POTT,QTRC(:,:,:,QS:QE),PRES,EXNER,QDRY,QV,QC,QI,SFC_DENS,SFC_PRES,SFLX_MU,SFLX_MV,SFLX_SH,SFLX_QV,Ustar,RLmo)
          call atmos_phy_bl_mynn_jmapplib_tendency( &
               ka, ks, ke, ia, is, ie, ja, js, je, &
               dens(:,:,:), u(:,:,:), v(:,:,:),                        & ! (in)
               pott(:,:,:), qtrc(:,:,:,qs:qe),                         & ! (in)
               pres(:,:,:), exner(:,:,:),                              & ! (in)
               qdry(:,:,:), qv(:,:,:), qc(:,:,:), qi(:,:,:),           & ! (in)
               sfc_dens(:,:), sfc_pres(:,:),                           & ! (in)
               sflx_mu(:,:), sflx_mv(:,:), sflx_sh(:,:), sflx_qv(:,:), & ! (in)
               ustar(:,:), rlmo(:,:),                                  & ! (in)
               cz(:,:,:), fz(:,:,:), f2h(:,:,:,:), dt_bl,              & ! (in)
               rhou_t_bl(:,:,:), rhov_t_bl(:,:,:), rhot_t_bl(:,:,:),   & ! (out)
               rhoqv_t(:,:,:), rhoq_t_bl(:,:,:,qs:qe),                 & ! (out)
               nu(:,:,:), kh(:,:,:),                                   & ! (out)
               zi = zi(:,:), sflx_buoy = sflx_buoy(:,:)                ) ! (out)
          !$acc update device(RHOU_t_BL,RHOV_t_BL,RHOT_t_BL,RHOQV_t,RHOQ_t_BL,Nu,Kh,Zi,SFLX_BUOY)
          if ( i_qv <= 0 ) then
             !$acc exit data delete(RHOQV_t)
             deallocate( rhoqv_t )
          end if
       end select
 
       if ( atmos_phy_bl_mix_tracers ) then
          do iq = 1, qa
             if ( ( .not. tracer_advc(iq) ) .or. iq==i_qv .or. (iq>=qs .and. iq<=qe) ) cycle
             call atmos_phy_bl_tendency_tracer( &
                  ka, ks, ke, ia, is, ie, ja, js, je, &
                  dens(:,:,:), qtrc(:,:,:,iq),        & ! (in)
                  sflx_q(:,:,iq),                     & ! (in)
                  kh(:,:,:),                          & ! (in)
                  tracer_mass(iq),                    & ! (in)
                  cz(:,:,:), fz(:,:,:), f2h(:,:,:,:), & ! (in)
                  dt_bl, tracer_name(iq),             & ! (in)
                  rhoq_t_bl(:,:,:,iq)                 ) ! (out)
          end do
       end if
 
       call file_history_in( nu(:,:,:),     'Nu_BL',     'eddy viscosity',                         'm2/s',      fill_halo=.true., dim_type="ZHXY" )
       call file_history_in( kh(:,:,:),     'Kh_BL',     'eddy diffusion',                         'm2/s',      fill_halo=.true., dim_type="ZHXY" )
 
       call file_history_in( ql(:,:,:),    'QL_BL',      'liquid water content in partial condensation', 'kg/kg', fill_halo=.true. )
       call file_history_in( cldfrac(:,:,:),   'cldfrac_BL', 'cloud fraction in partial condensation',       '1',     fill_halo=.true. )
 
       call file_history_in( zi(:,:),     'Zi_BL',      'depth of the boundary layer', 'm', fill_halo=.true. )
 
       call file_history_in( rhou_t_bl(:,:,:), 'RHOU_t_BL', 'MOMX tendency (BL)',                     'kg/m2/s2',  fill_halo=.true. )
       call file_history_in( rhov_t_bl(:,:,:), 'RHOV_t_BL', 'MOMY tendency (BL)',                     'kg/m2/s2',  fill_halo=.true. )
       call file_history_in( rhot_t_bl(:,:,:), 'RHOT_t_BL', 'RHOT tendency (BL)',                     'K.kg/m3/s', fill_halo=.true. )
 
       do iq = 1, qa
          if ( .not. tracer_advc(iq) ) cycle
          call file_history_in( rhoq_t_bl(:,:,:,iq), trim(tracer_name(iq))//'_t_BL',                      &
                        'RHO*'//trim(tracer_name(iq))//' tendency (BL)', 'kg/m3/s', fill_halo=.true. )
       enddo
 
       if ( statistics_checktotal ) then
          call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                                 rhou_t_bl(:,:,:), 'RHOU_t_BL',      &
                                 atmos_grid_cartesc_real_vol(:,:,:), &
                                 atmos_grid_cartesc_real_totvol      )
          call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                                 rhov_t_bl(:,:,:), 'RHOV_t_BL',      &
                                 atmos_grid_cartesc_real_vol(:,:,:), &
                                 atmos_grid_cartesc_real_totvol      )
          call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                                 rhot_t_bl(:,:,:), 'RHOT_t_BL',      &
                                 atmos_grid_cartesc_real_vol(:,:,:), &
                                 atmos_grid_cartesc_real_totvol      )
          call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                                 nu(:,:,:),        'Nu_BL',          &
                                 atmos_grid_cartesc_real_vol(:,:,:), &
                                 atmos_grid_cartesc_real_totvol      )
          call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                                 kh(:,:,:),         'Kh_BL',         &
                                 atmos_grid_cartesc_real_vol(:,:,:), &
                                 atmos_grid_cartesc_real_totvol      )
 
          do iq = 1, qa
             if ( .not. tracer_advc(iq) ) cycle
             call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                                    rhoq_t_bl(:,:,:,iq), trim(tracer_name(iq))//'_t_BL', &
                                    atmos_grid_cartesc_real_vol(:,:,:),                  &
                                    atmos_grid_cartesc_real_totvol                       )
          enddo
       endif
 
       !$acc end data
 
    endif
 
    !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
    !$omp shared(JS,JE,IS,IE,KS,KE,RHOU_t,RHOU_t_BL,RHOV_t,RHOV_t_BL,RHOT_t,RHOT_t_BL)
    !$acc kernels
    do j = js, je
    do i = is, ie
    do k = ks, ke
       rhou_t(k,i,j) = rhou_t(k,i,j) + rhou_t_bl(k,i,j)
       rhov_t(k,i,j) = rhov_t(k,i,j) + rhov_t_bl(k,i,j)
       rhot_t(k,i,j) = rhot_t(k,i,j) + rhot_t_bl(k,i,j)
    enddo
    enddo
    enddo
    !$acc end kernels
 
    !$omp parallel private(iq,i,j,k)
    !$acc kernels
    do iq = 1,  qa
#ifndef _OPENACC
       if ( .not. tracer_advc(iq) ) cycle
#endif
       !$omp do OMP_SCHEDULE_ collapse(2)
       do j  = js, je
       do i  = is, ie
       do k  = ks, ke
#ifdef _OPENACC
          if ( tracer_advc(iq) ) &
#endif
          rhoq_t(k,i,j,iq) = rhoq_t(k,i,j,iq) + rhoq_t_bl(k,i,j,iq)
       enddo
       enddo
       enddo
       !$omp end do nowait
    enddo
    !$acc end kernels
    !$omp end parallel
 
    return
  end subroutine atmos_phy_bl_driver_calc_tendency
 
end module mod_atmos_phy_bl_driver