mod_atmos_phy_cp_vars.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
#include "scalelib.h"
module mod_atmos_phy_cp_vars
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_io
  use scale_prof
  use scale_debug
  use scale_atmos_grid_cartesc_index
  use scale_tracer
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: atmos_phy_cp_vars_setup
  public :: atmos_phy_cp_vars_finalize
  public :: atmos_phy_cp_vars_fillhalo
  public :: atmos_phy_cp_vars_restart_read
  public :: atmos_phy_cp_vars_restart_write
 
  public :: atmos_phy_cp_vars_restart_create
  public :: atmos_phy_cp_vars_restart_open
  public :: atmos_phy_cp_vars_restart_def_var
  public :: atmos_phy_cp_vars_restart_enddef
  public :: atmos_phy_cp_vars_restart_close
 
  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  logical,               public :: atmos_phy_cp_restart_output                 = .false.                
 
  character(len=H_LONG),  public :: atmos_phy_cp_restart_in_basename           = ''
  logical,                public :: atmos_phy_cp_restart_in_postfix_timelabel  = .false.                
  logical,                public :: atmos_phy_cp_restart_in_aggregate
  character(len=H_LONG),  public :: atmos_phy_cp_restart_out_basename          = ''
  logical,                public :: atmos_phy_cp_restart_out_aggregate
  logical,                public :: atmos_phy_cp_restart_out_postfix_timelabel = .true.                 
  character(len=H_MID),   public :: atmos_phy_cp_restart_out_title             = 'ATMOS_PHY_CP restart'
  character(len=H_SHORT), public :: atmos_phy_cp_restart_out_dtype             = 'DEFAULT'
 
  ! restart variables
  real(rp), public, allocatable :: atmos_phy_cp_dens_t   (:,:,:)   ! tendency DENS     [kg/m3/s]
  real(rp), public, allocatable :: atmos_phy_cp_momz_t   (:,:,:)   ! tendency MOMZ     [kg/m2/s2]
  real(rp), public, allocatable :: atmos_phy_cp_rhot_t   (:,:,:)   ! tendency RHOT     [K*kg/m3/s]
  real(rp), public, allocatable :: atmos_phy_cp_rhoqv_t  (:,:,:)   ! tendency rho*QV   [kg/kg/s]
  real(rp), public, allocatable :: atmos_phy_cp_rhohyd_t (:,:,:,:) ! tendency rho*QHYD [kg/kg/s]
  real(rp), public, allocatable :: atmos_phy_cp_sflx_rain(:,:)     ! convective rain   [kg/m2/s]
  real(rp), public, allocatable :: atmos_phy_cp_sflx_snow(:,:)     ! convective snow   [kg/m2/s]
  real(rp), public, allocatable :: atmos_phy_cp_sflx_engi(:,:)     ! internal energy   [J/m2/s]
 
  ! only for K-F scheme
  real(rp), public, allocatable :: atmos_phy_cp_w0mean        (:,:,:) ! running mean vertical wind velocity [m/s]
  real(rp), public, allocatable :: atmos_phy_cp_kf_nca        (:,:)   ! advection/cumulus convection timescale for KF[sec]
 
 
  ! diagnostic variables
  real(rp), public, allocatable :: atmos_phy_cp_mflx_cloudbase(:,:)   ! cloud base mass flux [kg/m2/s]
  real(rp), public, allocatable :: atmos_phy_cp_cloudtop      (:,:)   ! cloud top  height [m]
  real(rp), public, allocatable :: atmos_phy_cp_cloudbase     (:,:)   ! cloud base height [m]
  real(rp), public, allocatable :: atmos_phy_cp_cldfrac_dp    (:,:,:) ! cloud fraction (deep    convection) (0-1)
  real(rp), public, allocatable :: atmos_phy_cp_cldfrac_sh    (:,:,:) ! cloud fraction (shallow convection) (0-1)
 
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  integer,                private, parameter :: vmax = 2
  integer,                private, parameter :: i_w0mean = 1
  integer,                private, parameter :: i_kf_nca = 2
 
 
  character(len=H_SHORT), private            :: var_name(vmax)
  character(len=H_MID),   private            :: var_desc(vmax)
  character(len=H_SHORT), private            :: var_unit(vmax)
  character(len=H_SHORT), private            :: var_dim (vmax)
  integer,                private            :: var_id  (vmax)
  integer,                private            :: restart_fid = -1  ! file ID
 
  data var_name / 'w0mean', &
                  'kf_nca'  /
 
  data var_desc / 'running mean vertical velocity',                &
                  'advection/cumulus convection timescale for KF'  /
 
  data var_unit / 'm/s', &
                  'sec'  /
 
  data var_dim  / 'ZXY',  &
                  'XY'    /
 
  ! surface flux and tendency names
  integer,                private              :: vmax_t
  integer,                private              :: i_cp_rain_sf = 1
  integer,                private              :: i_cp_snow_sf = 2
  integer,                private              :: i_cp_engi_sf = 3
  integer,                private              :: i_cp_dens_t  = 4
  integer,                private              :: i_cp_rhot_t  = 5
  integer,                private              :: i_cp_qv_t    = 6
 
  character(len=H_SHORT), private, allocatable :: var_t_name(:)
  character(len=H_MID),   private, allocatable :: var_t_desc(:)
  character(len=H_SHORT), private, allocatable :: var_t_unit(:)
  integer,                private, allocatable :: var_t_id  (:)
 
  !-----------------------------------------------------------------------------
contains
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_setup
    use scale_prc, only: &
       prc_abort
    use scale_const, only: &
       undef => const_undef
    use scale_atmos_hydrometeor, only: &
       n_hyd, &
       hyd_name
    implicit none
 
    namelist / param_atmos_phy_cp_vars / &
       atmos_phy_cp_restart_in_basename,           &
       atmos_phy_cp_restart_in_aggregate,          &
       atmos_phy_cp_restart_in_postfix_timelabel,  &
       atmos_phy_cp_restart_output,                &
       atmos_phy_cp_restart_out_basename,          &
       atmos_phy_cp_restart_out_aggregate,         &
       atmos_phy_cp_restart_out_postfix_timelabel, &
       atmos_phy_cp_restart_out_title,             &
       atmos_phy_cp_restart_out_dtype
 
    integer :: ierr
    integer :: iv
    integer :: iq
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("ATMOS_PHY_CP_vars_setup",*) 'Setup'
 
    allocate( atmos_phy_cp_dens_t(ka,ia,ja)       )
    allocate( atmos_phy_cp_momz_t(ka,ia,ja)       )
    allocate( atmos_phy_cp_rhot_t(ka,ia,ja)       )
    allocate( atmos_phy_cp_rhoqv_t(ka,ia,ja)       )
    allocate( atmos_phy_cp_rhohyd_t(ka,ia,ja,n_hyd) )
    atmos_phy_cp_dens_t(:,:,:)   = 0.0_rp
    atmos_phy_cp_momz_t(:,:,:)   = undef
    atmos_phy_cp_rhot_t(:,:,:)   = 0.0_rp
    atmos_phy_cp_rhoqv_t(:,:,:)   = 0.0_rp
    atmos_phy_cp_rhohyd_t(:,:,:,:) = 0.0_rp
    !$acc enter data copyin(ATMOS_PHY_CP_DENS_t,ATMOS_PHY_CP_MOMZ_t,ATMOS_PHY_CP_RHOT_t,ATMOS_PHY_CP_RHOQV_t,ATMOS_PHY_CP_RHOHYD_t)
 
    allocate( atmos_phy_cp_w0mean(ka,ia,ja) )
    allocate( atmos_phy_cp_kf_nca(ia,ja)    )
    atmos_phy_cp_w0mean(:,:,:) =    0.0_rp
    atmos_phy_cp_kf_nca(:,:)   = -100.0_rp
    !$acc enter data copyin(ATMOS_PHY_CP_w0mean,ATMOS_PHY_CP_kf_nca)
 
    ! for surface flux and tendency restart
    vmax_t = 6 + n_hyd
    allocate( var_t_name(vmax_t) )
    allocate( var_t_desc(vmax_t) )
    allocate( var_t_unit(vmax_t) )
    allocate( var_t_id(vmax_t) )
 
    var_t_name(i_cp_rain_sf) = 'SFLX_RAIN_CP'
    var_t_desc(i_cp_rain_sf) = 'surface rain flux in CP'
    var_t_unit(i_cp_rain_sf) = 'kg/m2/s'
    var_t_name(i_cp_snow_sf) = 'SFLX_SNOW_CP'
    var_t_desc(i_cp_snow_sf) = 'surface snow flux in CP'
    var_t_unit(i_cp_snow_sf) = 'kg/m2/s'
    var_t_name(i_cp_engi_sf) = 'SFLX_ENGI_CP'
    var_t_desc(i_cp_engi_sf) = 'surface internal energy flux in CP'
    var_t_unit(i_cp_engi_sf) = 'J/m2/s'
 
    var_t_name(i_cp_dens_t) = 'DENS_t_CP'
    var_t_desc(i_cp_dens_t) = 'tendency DENS in CP'
    var_t_unit(i_cp_dens_t) = 'kg/m3/s'
    var_t_name(i_cp_rhot_t) = 'RHOT_t_CP'
    var_t_desc(i_cp_rhot_t) = 'tendency RHOT in CP'
    var_t_unit(i_cp_rhot_t) = 'K*kg/m3/s'
 
    var_t_name(i_cp_qv_t) = 'QV_t_CP'
    var_t_desc(i_cp_qv_t) = 'tendency rho*QV in CP'
    var_t_unit(i_cp_qv_t) = 'kg/m3/s'
    do iq = 1, n_hyd
       var_t_name(6+iq) = trim(hyd_name(iq))//'_t_CP'
       var_t_desc(6+iq) = 'tendency rho*'//trim(hyd_name(iq))//' in CP'
       var_t_unit(6+iq) = 'kg/m3/s'
    enddo
 
 
    allocate( atmos_phy_cp_mflx_cloudbase(ia,ja)    )
    allocate( atmos_phy_cp_sflx_rain(ia,ja)    )
    allocate( atmos_phy_cp_sflx_snow(ia,ja)    )
    allocate( atmos_phy_cp_sflx_engi(ia,ja)    )
    allocate( atmos_phy_cp_cloudtop(ia,ja)    )
    allocate( atmos_phy_cp_cloudbase(ia,ja)    )
    allocate( atmos_phy_cp_cldfrac_dp(ka,ia,ja) )
    allocate( atmos_phy_cp_cldfrac_sh(ka,ia,ja) )
    atmos_phy_cp_mflx_cloudbase(:,:)   =    0.0_rp
    atmos_phy_cp_sflx_rain(:,:)   =    0.0_rp
    atmos_phy_cp_sflx_snow(:,:)   =    0.0_rp
    atmos_phy_cp_sflx_engi(:,:)   =    0.0_rp
    atmos_phy_cp_cloudtop(:,:)   =    0.0_rp
    atmos_phy_cp_cloudbase(:,:)   =    0.0_rp
    atmos_phy_cp_cldfrac_dp(:,:,:) =    0.0_rp
    atmos_phy_cp_cldfrac_sh(:,:,:) =    0.0_rp
    !$acc enter data copyin(ATMOS_PHY_CP_MFLX_cloudbase,ATMOS_PHY_CP_SFLX_rain,ATMOS_PHY_CP_SFLX_snow,ATMOS_PHY_CP_SFLX_ENGI,ATMOS_PHY_CP_cloudtop,ATMOS_PHY_CP_cloudbase,ATMOS_PHY_CP_cldfrac_dp,ATMOS_PHY_CP_cldfrac_sh)
 
    !--- read namelist
    rewind(io_fid_conf)
    read(io_fid_conf,nml=param_atmos_phy_cp_vars,iostat=ierr)
    if( ierr < 0 ) then !--- missing
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Not found namelist. Default used.'
    elseif( ierr > 0 ) then !--- fatal error
       log_error("ATMOS_PHY_CP_vars_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_PHY_CP_VARS. Check!'
       call prc_abort
    endif
    log_nml(param_atmos_phy_cp_vars)
 
    log_newline
    log_info("ATMOS_PHY_CP_vars_setup",*) '[ATMOS_PHY_CP] prognostic/diagnostic variables'
    log_info_cont('(1x,A,A24,A,A48,A,A12,A)') &
               '      |', 'VARNAME                 ','|', &
               'DESCRIPTION                                     ', '[', 'UNIT        ', ']'
    do iv = 1, vmax
       log_info_cont('(1x,A,I3,A,A24,A,A48,A,A12,A)') &
                  'NO.',iv,'|',var_name(iv),'|',var_desc(iv),'[',var_unit(iv),']'
    enddo
 
    ! tendency
    do iv = 1, vmax_t
       log_info_cont('(1x,A,I3,A,A24,A,A48,A,A12,A)') &
                  'NO.',iv+vmax,'|',var_t_name(iv),'|',var_t_desc(iv),'[',var_t_unit(iv),']'
    enddo
 
    log_newline
    if ( atmos_phy_cp_restart_in_basename /= '' ) then
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Restart input?  : YES, file = ', trim(atmos_phy_cp_restart_in_basename)
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Add timelabel?  : ', atmos_phy_cp_restart_in_postfix_timelabel
    else
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Restart input?  : NO'
    endif
    if (       atmos_phy_cp_restart_output             &
         .AND. atmos_phy_cp_restart_out_basename /= '' ) then
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Restart output? : YES, file = ', trim(atmos_phy_cp_restart_out_basename)
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Add timelabel?  : ', atmos_phy_cp_restart_out_postfix_timelabel
    else
       log_info("ATMOS_PHY_CP_vars_setup",*) 'Restart output? : NO'
       atmos_phy_cp_restart_output = .false.
    endif
 
    return
  end subroutine atmos_phy_cp_vars_setup
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_finalize
    implicit none
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("ATMOS_PHY_CP_vars_finalize",*) 'Finalize'
 
    !$acc exit data delete(ATMOS_PHY_CP_DENS_t,ATMOS_PHY_CP_MOMZ_t,ATMOS_PHY_CP_RHOT_t,ATMOS_PHY_CP_RHOQV_t,ATMOS_PHY_CP_RHOHYD_t)
    deallocate( atmos_phy_cp_dens_t     )
    deallocate( atmos_phy_cp_momz_t     )
    deallocate( atmos_phy_cp_rhot_t     )
    deallocate( atmos_phy_cp_rhoqv_t      )
    deallocate( atmos_phy_cp_rhohyd_t )
 
    !$acc exit data delete(ATMOS_PHY_CP_w0mean,ATMOS_PHY_CP_kf_nca)
    deallocate( atmos_phy_cp_w0mean        )
    deallocate( atmos_phy_cp_kf_nca        )
 
    ! for tendency restart
    deallocate( var_t_name )
    deallocate( var_t_desc )
    deallocate( var_t_unit )
    deallocate( var_t_id   )
 
    !$acc exit data delete(ATMOS_PHY_CP_MFLX_cloudbase,ATMOS_PHY_CP_SFLX_rain,ATMOS_PHY_CP_SFLX_snow,ATMOS_PHY_CP_SFLX_ENGI,ATMOS_PHY_CP_cloudtop,ATMOS_PHY_CP_cloudbase,ATMOS_PHY_CP_cldfrac_dp,ATMOS_PHY_CP_cldfrac_sh)
    deallocate( atmos_phy_cp_mflx_cloudbase    )
    deallocate( atmos_phy_cp_sflx_rain         )
    deallocate( atmos_phy_cp_sflx_snow         )
    deallocate( atmos_phy_cp_sflx_engi         )
    deallocate( atmos_phy_cp_cloudtop          )
    deallocate( atmos_phy_cp_cloudbase         )
    deallocate( atmos_phy_cp_cldfrac_dp    )
    deallocate( atmos_phy_cp_cldfrac_sh    )
 
    return
  end subroutine atmos_phy_cp_vars_finalize
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_fillhalo
    use scale_comm_cartesc, only: &
       comm_vars8, &
       comm_wait
    use scale_atmos_hydrometeor, only: &
       n_hyd
    implicit none
 
    integer :: i, j
    integer :: iq
    !---------------------------------------------------------------------------
 
    !$omp parallel do
    !$acc parallel vector_length(32)
    !$acc loop collapse(2) independent
    do j = js, je
    do i = is, ie
       atmos_phy_cp_w0mean(   1:ks-1,i,j) = atmos_phy_cp_w0mean(ks,i,j)
       atmos_phy_cp_w0mean(ke+1:ka  ,i,j) = atmos_phy_cp_w0mean(ke,i,j)
       atmos_phy_cp_dens_t(   1:ks-1,i,j) = atmos_phy_cp_dens_t(ks,i,j)
       atmos_phy_cp_dens_t(ke+1:ka  ,i,j) = atmos_phy_cp_dens_t(ke,i,j)
       atmos_phy_cp_rhot_t(   1:ks-1,i,j) = atmos_phy_cp_rhot_t(ks,i,j)
       atmos_phy_cp_rhot_t(ke+1:ka  ,i,j) = atmos_phy_cp_rhot_t(ke,i,j)
       atmos_phy_cp_rhoqv_t(   1:ks-1,i,j) = atmos_phy_cp_rhoqv_t(ks,i,j)
       atmos_phy_cp_rhoqv_t(ke+1:ka  ,i,j) = atmos_phy_cp_rhoqv_t(ke,i,j)
    end do
    end do
    !$acc end parallel
    !$omp parallel do collapse(2)
    !$acc parallel vector_length(32)
    !$acc loop collapse(3) independent
    do iq = 1, n_hyd
       do j = js, je
       do i = is, ie
          atmos_phy_cp_rhohyd_t(   1:ks-1,i,j,iq) = atmos_phy_cp_rhohyd_t(ks,i,j,iq)
          atmos_phy_cp_rhohyd_t(ke+1:ka  ,i,j,iq) = atmos_phy_cp_rhohyd_t(ke,i,j,iq)
       enddo
       enddo
    end do
    !$acc end parallel
 
    call comm_vars8( atmos_phy_cp_w0mean(:,:,:), 1 )
    call comm_vars8( atmos_phy_cp_kf_nca(:,:),   2 )
 
    ! tendency
    call comm_vars8( atmos_phy_cp_dens_t(:,:,:), vmax+1 )
    call comm_vars8( atmos_phy_cp_rhot_t(:,:,:), vmax+2 )
    call comm_vars8( atmos_phy_cp_rhoqv_t(:,:,:), vmax+3 )
    do iq = 1, n_hyd
       call comm_vars8( atmos_phy_cp_rhohyd_t(:,:,:,iq), vmax+3+iq )
    enddo
 
    call comm_wait ( atmos_phy_cp_w0mean(:,:,:), 1 )
    call comm_wait ( atmos_phy_cp_kf_nca(:,:),   2 )
 
    call comm_wait ( atmos_phy_cp_dens_t(:,:,:), vmax+1 )
    call comm_wait ( atmos_phy_cp_rhot_t(:,:,:), vmax+2 )
    call comm_wait ( atmos_phy_cp_rhoqv_t(:,:,:), vmax+3 )
    do iq = 1, n_hyd
       call comm_wait ( atmos_phy_cp_rhohyd_t(:,:,:,iq), vmax+3+iq )
    enddo
 
    return
  end subroutine atmos_phy_cp_vars_fillhalo
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_open
    use scale_time, only: &
       time_gettimelabel
    use scale_file_cartesc, only: &
       file_cartesc_open
    implicit none
 
    character(len=19)     :: timelabel
    character(len=H_LONG) :: basename
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("ATMOS_PHY_CP_vars_restart_open",*) 'Open restart file (ATMOS_PHY_CP) '
 
    if ( atmos_phy_cp_restart_in_basename /= '' ) then
 
       if ( atmos_phy_cp_restart_in_postfix_timelabel ) then
          call time_gettimelabel( timelabel )
          basename = trim(atmos_phy_cp_restart_in_basename)//'_'//trim(timelabel)
       else
          basename = trim(atmos_phy_cp_restart_in_basename)
       endif
 
       log_info("ATMOS_PHY_CP_vars_restart_open",*) 'basename: ', trim(basename)
 
       call file_cartesc_open( basename, restart_fid, aggregate=atmos_phy_cp_restart_in_aggregate )
    else
       log_info("ATMOS_PHY_CP_vars_restart_open",*) 'restart file for ATMOS_PHY_CP is not specified.'
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_open
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_read
    use scale_file, only: &
       file_get_aggregate
    use scale_file_cartesc, only: &
       file_cartesc_read, &
       file_cartesc_flush
    use scale_atmos_hydrometeor, only: &
       n_hyd
    implicit none
 
    integer  :: i, j, iq
    !---------------------------------------------------------------------------
 
    if ( restart_fid /= -1 ) then
       log_newline
       log_info("ATMOS_PHY_CP_vars_restart_read",*) 'Read from restart file (ATMOS_PHY_CP) '
 
       call file_cartesc_read( restart_fid, var_name(1), 'ZXY', & ! [IN]
                               atmos_phy_cp_w0mean(:,:,:)       ) ! [OUT]
       call file_cartesc_read( restart_fid, var_name(2), 'XY',  & ! [IN]
                               atmos_phy_cp_kf_nca(:,:)         ) ! [OUT]
 
       ! surface flux
       call file_cartesc_read( restart_fid, var_t_name(1), 'XY', & ! [IN]
                               atmos_phy_cp_sflx_rain(:,:)       ) ! [OUT]
       call file_cartesc_read( restart_fid, var_t_name(2), 'XY', & ! [IN]
                               atmos_phy_cp_sflx_snow(:,:)       ) ! [OUT]
       call file_cartesc_read( restart_fid, var_t_name(3), 'XY', & ! [IN]
                               atmos_phy_cp_sflx_engi(:,:)       ) ! [OUT]
 
       ! tendency
       call file_cartesc_read( restart_fid, var_t_name(4), 'ZXY', & ! [IN]
                               atmos_phy_cp_dens_t(:,:,:)         ) ! [OUT]
       call file_cartesc_read( restart_fid, var_t_name(5), 'ZXY', & ! [IN]
                               atmos_phy_cp_rhot_t(:,:,:)         ) ! [OUT]
       call file_cartesc_read( restart_fid, var_t_name(6), 'ZXY', & ! [IN]
                               atmos_phy_cp_rhoqv_t(:,:,:)        ) ! [OUT]
       do iq = 1, n_hyd
          call file_cartesc_read( restart_fid, var_t_name(6+iq), 'ZXY', & ! [IN]
                                  atmos_phy_cp_rhohyd_t(:,:,:,iq)       ) ! [OUT]
       enddo
 
       if ( file_get_aggregate(restart_fid) ) then
          call file_cartesc_flush( restart_fid ) ! X/Y halos have been read from file
          !$acc update device(ATMOS_PHY_CP_w0mean,ATMOS_PHY_CP_kf_nca,ATMOS_PHY_CP_SFLX_rain,ATMOS_PHY_CP_SFLX_snow,ATMOS_PHY_CP_SFLX_ENGI,ATMOS_PHY_CP_DENS_t,ATMOS_PHY_CP_RHOT_t,ATMOS_PHY_CP_RHOQV_t,ATMOS_PHY_CP_RHOHYD_t)
 
          ! fill K halos
          !$omp parallel do
          !$acc parallel vector_length(32)
          !$acc loop collapse(2) independent
          do j  = 1, ja
          do i  = 1, ia
             atmos_phy_cp_w0mean(   1:ks-1,i,j) = atmos_phy_cp_w0mean(ks,i,j)
             atmos_phy_cp_w0mean(ke+1:ka  ,i,j) = atmos_phy_cp_w0mean(ke,i,j)
             atmos_phy_cp_dens_t(   1:ks-1,i,j) = atmos_phy_cp_dens_t(ks,i,j)
             atmos_phy_cp_dens_t(ke+1:ka  ,i,j) = atmos_phy_cp_dens_t(ke,i,j)
             atmos_phy_cp_rhot_t(   1:ks-1,i,j) = atmos_phy_cp_rhot_t(ks,i,j)
             atmos_phy_cp_rhot_t(ke+1:ka  ,i,j) = atmos_phy_cp_rhot_t(ke,i,j)
             atmos_phy_cp_rhoqv_t(   1:ks-1,i,j) = atmos_phy_cp_rhoqv_t(ks,i,j)
             atmos_phy_cp_rhoqv_t(ke+1:ka  ,i,j) = atmos_phy_cp_rhoqv_t(ke,i,j)
          end do
          end do
          !$acc end parallel
          !$omp parallel do collapse(2)
          !$acc parallel vector_length(32)
          !$acc loop collapse(3) independent
          do iq = 1, n_hyd
             do j  = 1, ja
             do i  = 1, ia
                atmos_phy_cp_rhohyd_t(   1:ks-1,i,j,iq) = atmos_phy_cp_rhohyd_t(ks,i,j,iq)
                atmos_phy_cp_rhohyd_t(ke+1:ka,  i,j,iq) = atmos_phy_cp_rhohyd_t(ke,i,j,iq)
             enddo
             enddo
          enddo
          !$acc end parallel
       else
          call atmos_phy_cp_vars_fillhalo
       end if
 
       call atmos_phy_cp_vars_check
 
    else
       log_info("ATMOS_PHY_CP_vars_restart_read",*) 'invalid restart file ID for ATMOS_PHY_CP.'
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_read
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_create
    use scale_time, only: &
       time_gettimelabel
    use scale_file_cartesc, only: &
       file_cartesc_create
    implicit none
 
    character(len=19)     :: timelabel
    character(len=H_LONG) :: basename
    !---------------------------------------------------------------------------
 
    if ( atmos_phy_cp_restart_out_basename /= '' ) then
 
       log_newline
       log_info("ATMOS_PHY_CP_vars_restart_create",*) 'Create restart file (ATMOS_PHY_AE) '
 
       if ( atmos_phy_cp_restart_out_postfix_timelabel ) then
          call time_gettimelabel( timelabel )
          basename = trim(atmos_phy_cp_restart_out_basename)//'_'//trim(timelabel)
       else
          basename = trim(atmos_phy_cp_restart_out_basename)
       endif
 
       log_info("ATMOS_PHY_CP_vars_restart_create",*) 'basename: ', trim(basename)
 
       call file_cartesc_create( &
            basename, atmos_phy_cp_restart_out_title, atmos_phy_cp_restart_out_dtype, & ! [IN]
            restart_fid,                                                              & ! [OUT]
            aggregate=atmos_phy_cp_restart_out_aggregate                              ) ! [IN]
 
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_create
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_enddef
    use scale_file_cartesc, only: &
       file_cartesc_enddef
    implicit none
 
    if ( restart_fid /= -1 ) then
       call file_cartesc_enddef( restart_fid ) ! [IN]
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_enddef
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_close
    use scale_file_cartesc, only: &
       file_cartesc_close
    implicit none
    !---------------------------------------------------------------------------
 
    if ( restart_fid /= -1 ) then
       log_newline
       log_info("ATMOS_PHY_CP_vars_restart_close",*) 'Close restart file (ATMOS_PHY_CP) '
 
       call file_cartesc_close( restart_fid ) ! [IN]
 
       restart_fid = -1
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_close
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_def_var
    use scale_file_cartesc, only: &
       file_cartesc_def_var
    use scale_atmos_hydrometeor, only: &
       n_hyd
    implicit none
 
    integer :: i, iq
    !---------------------------------------------------------------------------
 
    if ( restart_fid /= -1 ) then
 
       do i = 1, vmax
          call file_cartesc_def_var( restart_fid,          & ! [IN]
               var_name(i), var_desc(i), var_unit(i),      & ! [IN]
               var_dim(i), atmos_phy_cp_restart_out_dtype, & ! [IN]
               var_id(i)                                   ) ! [OUT]
       end do
 
       do i = 1, 3
          call file_cartesc_def_var( restart_fid,           & ! [IN]
               var_t_name(i), var_t_desc(i), var_t_unit(i), & ! [IN]
               'XY',  atmos_phy_cp_restart_out_dtype,       & ! [IN]
               var_t_id(i)                                  ) ! [OUT]
       end do
 
       do i = 4, 6
          call file_cartesc_def_var( restart_fid,           & ! [IN]
               var_t_name(i), var_t_desc(i), var_t_unit(i), & ! [IN]
               'ZXY',  atmos_phy_cp_restart_out_dtype,      & ! [IN]
               var_t_id(i)                                  ) ! [OUT]
       end do
 
       do iq = 1, n_hyd
          call file_cartesc_def_var( restart_fid,                    & ! [IN]
               var_t_name(6+iq), var_t_desc(6+iq), var_t_unit(6+iq), & ! [IN]
               'ZXY',  atmos_phy_cp_restart_out_dtype,               & ! [IN]
               var_t_id(6+iq)                                        ) ! [OUT]
       enddo
 
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_def_var
 
  !-----------------------------------------------------------------------------
  subroutine atmos_phy_cp_vars_restart_write
    use scale_file_cartesc, only: &
       file_cartesc_write => file_cartesc_write_var
    use scale_atmos_hydrometeor, only: &
       n_hyd
    implicit none
 
    integer  :: iq
    !---------------------------------------------------------------------------
 
    if ( restart_fid /= -1 ) then
 
       call atmos_phy_cp_vars_fillhalo
 
       call atmos_phy_cp_vars_check
 
       call file_cartesc_write( restart_fid, var_id(1), atmos_phy_cp_w0mean(:,:,:), & ! [IN]
                                var_name(1), 'ZXY' ) ! [IN]
       call file_cartesc_write( restart_fid, var_id(2), atmos_phy_cp_kf_nca(:,:),   & ! [IN]
                                var_name(2), 'XY'  ) ! [IN]
 
       ! surface flux
       call file_cartesc_write( restart_fid, var_t_id(1), atmos_phy_cp_sflx_rain(:,:), & ! [IN]
                                var_t_name(1), 'XY' ) ! [IN]
       call file_cartesc_write( restart_fid, var_t_id(2), atmos_phy_cp_sflx_snow(:,:), & ! [IN]
                                var_t_name(2), 'XY' ) ! [IN]
       call file_cartesc_write( restart_fid, var_t_id(3), atmos_phy_cp_sflx_engi(:,:), & ! [IN]
                                var_t_name(3), 'XY' ) ! [IN]
 
       ! tendency
       call file_cartesc_write( restart_fid, var_t_id(4), atmos_phy_cp_dens_t(:,:,:), & ! [IN]
                                var_t_name(4), 'ZXY' ) ! [IN]
       call file_cartesc_write( restart_fid, var_t_id(5), atmos_phy_cp_rhot_t(:,:,:), & ! [IN]
                                var_t_name(5), 'ZXY' ) ! [IN]
       call file_cartesc_write( restart_fid, var_t_id(6), atmos_phy_cp_rhoqv_t(:,:,:), & ! [IN]
                                var_t_name(6), 'ZXY' ) ! [IN]
       do iq = 1, n_hyd
          call file_cartesc_write( restart_fid, var_t_id(6+iq), atmos_phy_cp_rhohyd_t(:,:,:,iq), & ! [IN]
                             var_t_name(6+iq), 'ZXY' ) ! [IN]
       enddo
 
    endif
 
    return
  end subroutine atmos_phy_cp_vars_restart_write
 
  subroutine atmos_phy_cp_vars_check
    use scale_statistics, only: &
       statistics_checktotal, &
       statistics_total
    use scale_atmos_grid_cartesc_real, only: &
       atmos_grid_cartesc_real_area,    &
       atmos_grid_cartesc_real_totarea, &
       atmos_grid_cartesc_real_vol,     &
       atmos_grid_cartesc_real_totvol
    use scale_atmos_hydrometeor, only: &
       n_hyd
    implicit none
    integer :: iq
    !---------------------------------------------------------------------------
 
    call valcheck( ka, ks, ke, ia, is, ie, ja, js, je, &
                   atmos_phy_cp_w0mean(:,:,:),         & ! (in)
                   -100.0_rp, 100.0_rp, var_name(1),   & ! (in)
                   __file__, __line__                  ) ! (in)
    call valcheck( ia, is, ie, ja, js, je, &
                   atmos_phy_cp_kf_nca(:,:),           & ! (in)
                   -100.0_rp, 1.0e5_rp, var_name(2),   & ! (in)
                   __file__, __line__                  ) ! (in)
    ! surface flux
    call valcheck( ia, is, ie, ja, js, je, &
                   atmos_phy_cp_sflx_rain(:,:),        & ! (in)
                   0.0e0_rp, 1.0e0_rp, var_t_name(1),  & ! (in)
                   __file__, __line__                  ) ! (in)
    call valcheck( ia, is, ie, ja, js, je, &
                   atmos_phy_cp_sflx_snow(:,:),        & ! (in)
                   0.0e0_rp, 1.0e0_rp, var_t_name(2),  & ! (in)
                   __file__, __line__                  ) ! (in)
    call valcheck( ia, is, ie, ja, js, je, &
                   atmos_phy_cp_sflx_engi(:,:),        & ! (in)
                   0.0e0_rp, 5.0e3_rp, var_t_name(3),  & ! (in)
                   __file__, __line__                  ) ! (in)
    ! tendency
    call valcheck( ka, ks, ke, ia, is, ie, ja, js, je, &
                   atmos_phy_cp_dens_t(:,:,:),         & ! (in)
                   -1.0e0_rp, 1.0e0_rp, var_t_name(4), & ! (in)
                   __file__, __line__                  ) ! (in)
    call valcheck( ka, ks, ke, ia, is, ie, ja, js, je, &
                   atmos_phy_cp_rhot_t(:,:,:),         & ! (in)
                   -1.0e3_rp, 1.0e3_rp, var_t_name(5), & ! (in)
                   __file__, __line__                  ) ! (in)
    call valcheck( ka, ks, ke, ia, is, ie, ja, js, je, &
                   atmos_phy_cp_rhoqv_t(:,:,:),        & ! (in)
                   -1.0e0_rp, 1.0e0_rp, var_t_name(6), & ! (in)
                   __file__, __line__                  ) ! (in)
    do iq = 1, n_hyd
       call valcheck( ka, ks, ke, ia, is, ie, ja, js, je, &
                      atmos_phy_cp_rhohyd_t(:,:,:,iq),       & ! (in)
                      -1.0e0_rp, 1.0e0_rp, var_t_name(6+iq), & ! (in)
                      __file__, __line__                     ) ! (in)
    end do
 
    call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                           atmos_phy_cp_w0mean(:,:,:), var_name(1), & ! (in)
                           atmos_grid_cartesc_real_vol(:,:,:),              & ! (in)
                           atmos_grid_cartesc_real_totvol                   ) ! (in)
    call statistics_total( ia, is, ie, ja, js, je, &
                           atmos_phy_cp_kf_nca(:,:),  var_name(2), & ! (in)
                           atmos_grid_cartesc_real_area(:,:),               & ! (in)
                           atmos_grid_cartesc_real_totarea                  ) ! (in)
    ! surface flux
    call statistics_total( ia, is, ie, ja, js, je, &
                           atmos_phy_cp_sflx_rain(:,:), var_t_name(1), & ! (in)
                           atmos_grid_cartesc_real_area(:,:),                & ! (in)
                           atmos_grid_cartesc_real_totarea                   ) ! (in)
    call statistics_total( ia, is, ie, ja, js, je, &
                           atmos_phy_cp_sflx_snow(:,:), var_t_name(2), & ! (in)
                           atmos_grid_cartesc_real_area(:,:),                & ! (in)
                           atmos_grid_cartesc_real_totarea                   ) ! (in)
    call statistics_total( ia, is, ie, ja, js, je, &
                           atmos_phy_cp_sflx_engi(:,:), var_t_name(3), & ! (in)
                           atmos_grid_cartesc_real_area(:,:),                & ! (in)
                           atmos_grid_cartesc_real_totarea                   ) ! (in)
    ! tendency
    call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                           atmos_phy_cp_dens_t(:,:,:), var_t_name(4), & ! (in)
                           atmos_grid_cartesc_real_vol(:,:,:),                & ! (in)
                           atmos_grid_cartesc_real_totvol                     ) ! (in)
    call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                           atmos_phy_cp_rhot_t(:,:,:), var_t_name(5), & ! (in)
                           atmos_grid_cartesc_real_vol(:,:,:),                & ! (in)
                           atmos_grid_cartesc_real_totvol                     ) ! (in)
    call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                           atmos_phy_cp_rhoqv_t(:,:,:), var_t_name(6), & ! (in)
                           atmos_grid_cartesc_real_vol(:,:,:),                & ! (in)
                           atmos_grid_cartesc_real_totvol                     ) ! (in)
    do iq = 1, n_hyd
       call statistics_total( ka, ks, ke, ia, is, ie, ja, js, je, &
                              atmos_phy_cp_rhohyd_t(:,:,:,iq), var_t_name(6+iq), & ! (in)
                              atmos_grid_cartesc_real_vol(:,:,:),                & ! (in)
                              atmos_grid_cartesc_real_totvol                     ) ! (in)
    enddo
 
    return
  end subroutine atmos_phy_cp_vars_check
 
end module mod_atmos_phy_cp_vars