!= 簡単雲モデル
!
!= Simple cloud model
!
! Authors::   Yoshiyuki O. Takahashi
! Version::   $Id: cloud_simple.f90,v 1.4 2013-02-27 04:35:34 yot Exp $
! Tag Name::  $Name: dcpam5-20140204-3 $
! Copyright:: Copyright (C) GFD Dennou Club, 2008. All rights reserved.
! License::   See COPYRIGHT[link:../../../COPYRIGHT]
!

module cloud_simple
  !
  != 簡単雲モデル
  !
  != Simple cloud model
  !
  ! <b>Note that Japanese and English are described in parallel.</b>
  !
  ! 簡単雲モデルによる雲の計算.
  !
  ! In this module, the amount of cloud is calculated by use of a simple
  ! cloud model.
  !
  !== Procedures List
  !
!!$  ! RadiationFluxDennouAGCM :: 放射フラックスの計算
!!$  ! ------------            :: ------------
!!$  ! RadiationFluxDennouAGCM :: Calculate radiation flux
  !
  !== NAMELIST
  !
  ! NAMELIST#cloud_simple_nml
  !

  ! モジュール引用 ; USE statements

  !
  ! Kind type parameter
  !
  use dc_types, only: DP, &      ! Double precision.
    &                 STRING, &  ! Strings.
    &                 TOKEN      ! Keywords.

  ! メッセージ出力
  ! Message output
  !
  use dc_message, only: MessageNotify

  ! 格子点設定
  ! Grid points settings
  !
  use gridset, only: imax, & ! 経度格子点数.
                             ! Number of grid points in longitude
    &                jmax, & ! 緯度格子点数.
                             ! Number of grid points in latitude
    &                kmax    ! 鉛直層数.
                             ! Number of vertical level

  implicit none

  private


  ! 公開手続き
  ! Public procedure
  !
  public :: CloudSimple
  public :: CloudSimpleWithIce
  public :: CloudSimpleCalcCloudCover
  public :: CloudSimpleInit

  public :: CloudSimpleDivideWatAndIce
  public :: CloudSimpleCalcPRCPWithPC


  ! 公開変数
  ! Public variables
  !


  ! 非公開変数
  ! Private variables
  !
  logical , save        :: FlagSnow
                           ! A flag for snow

  integer , save        :: IDCloudCoverMethod
  integer , parameter   :: IDCloudCoverMethodConst    = 1
  integer , parameter   :: IDCloudCoverMethodRH       = 2

  real(DP), save        :: CloudLifeTime
  real(DP), save        :: CloudWatLifeTime
  real(DP), save        :: CloudIceLifeTime

  real(DP), save        :: CloudCover

  real(DP), save        :: RHCrtl


  logical, save :: cloud_simple_inited = .false.
                              ! 初期設定フラグ.
                              ! Initialization flag

  character(*), parameter:: module_name = 'cloud_simple'
                              ! モジュールの名称.
                              ! Module name
  character(*), parameter:: version = &
    & '$Name: dcpam5-20140204-3 $' // &
    & '$Id: cloud_simple.f90,v 1.4 2013-02-27 04:35:34 yot Exp $'
                              ! モジュールのバージョン
                              ! Module version

  !--------------------------------------------------------------------------------------

contains

  !--------------------------------------------------------------------------------------

  subroutine CloudSimple(                             &
    & xyr_Press, xyz_Temp,                            & ! (in)
    & xyz_DQH2OLiqDtCum, xyz_DQH2OLiqDtLSC,           & ! (in)
    & xyz_QH2OLiq, xyz_DTempDtCond,                   & ! (inout)
    & xy_RainCum, xy_SnowCum, xy_RainLsc, xy_SnowLsc  & ! (out)
    & )

    ! USE statements
    !

    ! 時刻管理
    ! Time control
    !
    use timeset, only: &
      & DelTime            ! $ \Delta t $ [s]

    ! 雲関系ルーチン
    ! Cloud-related routines
    !
    use cloud_utils, only : CloudUtilsCalcPRCPKeyLLTemp3D


    real(DP), intent(in   ) :: xyr_Press        ( 0:imax-1, 1:jmax, 0:kmax )
    real(DP), intent(in   ) :: xyz_Temp         ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in   ) :: xyz_DQH2OLiqDtCum( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in   ) :: xyz_DQH2OLiqDtLSC( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_QH2OLiq      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_DTempDtCond  ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(out  ) :: xy_RainCum       ( 0:imax-1, 1:jmax )
    real(DP), intent(out  ) :: xy_SnowCum       ( 0:imax-1, 1:jmax )
    real(DP), intent(out  ) :: xy_RainLsc       ( 0:imax-1, 1:jmax )
    real(DP), intent(out  ) :: xy_SnowLsc       ( 0:imax-1, 1:jmax )


    real(DP) :: xyz_QH2OLiqB( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQRainDt( 0:imax-1, 1:jmax, 1:kmax )

!!$    real(DP) :: xyz_DTempDtPrcpPCCum( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xyz_DTempDtPrcpPCLsc( 0:imax-1, 1:jmax, 1:kmax )

    ! 実行文 ; Executable statement
    !

    ! 初期化確認
    ! Initialization check
    !
    if ( .not. cloud_simple_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    ! Cloud optical depth
    !

    ! Numerical solution

!!$      xyz_DQCloudWaterDt = xyz_DQCloudWaterDt &
!!$        & - xyz_QCloudWater / ( CloudLifeTime + 1.0d-100 )


!      ( X_{t+1} - X_{t-1} ) / ( 2 \Delta t ) = Q - X_{t+1} / \tau
!
!      X_{t+1} / ( 2 \Delta t )  + X_{t+1} / \tau = X_{t-1} / ( 2 \Delta t ) + Q
!      ( 1 / ( 2 \Delta t )  + 1 / \tau ) X_{t+1} = X_{t-1} / ( 2 \Delta t ) + Q
!      X_{t+1} = ( X_{t-1} / ( 2 \Delta t ) + Q ) / ( 1 / ( 2 \Delta t )  + 1 / \tau ) 

!!$    xyz_QH2OLiq =                                                           &
!!$      &   (                                                                 &
!!$      &       xyz_QH2OLiq / ( 2.0_DP * DelTime )                            &
!!$      &     + xyz_DQH2OLiqDtCum + xyz_DQH2OLiqDtLSC                         &
!!$      &   )                                                                 &
!!$      & / ( 1.0_DP / ( 2.0_DP * DelTime ) + 1.0_DP / ( CloudLifeTime + 1.0d-100 ) )
!!$
!!$    call CloudUtilsCalcPRCPKeyLLTemp3D(          &
!!$      & xyr_Press, xyz_Temp, xyz_DQH2OLiqDtCum,  &  ! (in )
!!$      & xy_RainCum, xy_SnowCum                   &  ! (out)
!!$      & )
!!$    call CloudUtilsCalcPRCPKeyLLTemp3D(          &
!!$      & xyr_Press, xyz_Temp, xyz_DQH2OLiqDtLsc,  &  ! (in )
!!$      & xy_RainLsc, xy_SnowLsc                   &  ! (out)
!!$      & )
    !-----


    ! Analytical solution

    ! save cloud water before adjustment
    xyz_QH2OLiqB = xyz_QH2OLiq

    xyz_QH2OLiq =                                                                 &
      &   xyz_QH2OLiq * exp( - 2.0_DP * DelTime / ( CloudLifeTime + 1.0d-100 ) )  &
      & + ( xyz_DQH2OLiqDtCum + xyz_DQH2OLiqDtLSC ) * CloudLifeTime               &
      &   * ( 1.0_DP - exp( - 2.0_DP * DelTime / ( CloudLifeTime + 1.0d-100 ) ) )

    xyz_DQRainDt =                                                     &
      &   xyz_QH2OLiqB                                                 &
      & + ( xyz_DQH2OLiqDtCum + xyz_DQH2OLiqDtLSC ) * 2.0_DP * DelTime &
      & - xyz_QH2OLiq
    xyz_DQRainDt = xyz_DQRainDt / ( 2.0_DP * DelTime )

    call CloudUtilsCalcPRCPKeyLLTemp3D(    &
      & xyr_Press, xyz_Temp, xyz_DQRainDt, &  ! (in )
      & xy_RainCum, xy_SnowCum             &  ! (out)
      & )
    xy_RainLsc = 0.0_DP
    xy_SnowLsc = 0.0_DP


!!$        call CloudSimpleCalcPRCPWithPC(                &
!!$          & xyr_Press, xyz_TempA, xyz_DQH2OLiqDtCum,   & ! (in )
!!$          & xy_RainCumulus, xy_SnowCumulus,            & ! (out)
!!$          & xyz_DTempDtPrcpPCCum                       & ! (out)
!!$          & )
!!$        call CloudSimpleCalcPRCPWithPC(                &
!!$          & xyr_Press, xyz_TempA, xyz_DQH2OLiqDtLsc,   & ! (in )
!!$          & xy_RainLsc, xy_SnowLsc,                    & ! (out)
!!$          & xyz_DTempDtPrcpPCLsc                       & ! (out)
!!$          & )
!!$        xyz_TempA = xyz_TempA                                 &
!!$          & + ( xyz_DTempDtPrcpPCCum + xyz_DTempDtPrcpPCLsc ) &
!!$          &     * 2.0_DP * DelTime

!!$        xyz_DTempDtCond = xyz_DTempDtCond                  &
!!$          & + xyz_DTempDtPrcpPCCum + xyz_DTempDtPrcpPCLsc


  end subroutine CloudSimple

  !--------------------------------------------------------------------------------------

  subroutine CloudSimpleWithIce(                      &
    & xyr_Press,                                      & ! (in)
    & xyz_DQH2OLiqDtCum, xyz_DQH2OLiqDtLSC,           & ! (in)
    & xyz_Temp,                                       & ! (inout)
    & xyz_QH2OLiq, xyz_QH2OSol, xyz_DTempDtCond,      & ! (inout)
    & xy_RainCum, xy_SnowCum, xy_RainLsc, xy_SnowLsc  & ! (out)
    & )

    ! USE statements
    !

    ! 時刻管理
    ! Time control
    !
    use timeset, only: &
      & DelTime            ! $ \Delta t $ [s]

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only:  &
      & CpDry,            &
                              ! $ C_p $ [J kg-1 K-1].
                              ! 乾燥大気の定圧比熱.
                              ! Specific heat of air at constant pressure
      & Grav,             &
                              ! $ g $ [m s-2].
                              ! 重力加速度.
                              ! Gravitational acceleration
      & LatentHeatFusion
                              ! $ L $ [J kg-1] .
                              ! 融解の潜熱.
                              ! Latent heat of fusion

    ! 雪と海氷の定数の設定
    ! Setting constants of snow and sea ice
    !
    use constants_snowseaice, only: TempCondWater

    ! 雲関系ルーチン
    ! Cloud-related routines
    !
    use cloud_utils, only : CloudUtilsWatFraction

    real(DP), intent(in   ) :: xyr_Press        ( 0:imax-1, 1:jmax, 0:kmax )
    real(DP), intent(in   ) :: xyz_DQH2OLiqDtCum( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in   ) :: xyz_DQH2OLiqDtLSC( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_Temp         ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_QH2OLiq      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_QH2OSol      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_DTempDtCond  ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(out  ) :: xy_RainCum       ( 0:imax-1, 1:jmax )
    real(DP), intent(out  ) :: xy_SnowCum       ( 0:imax-1, 1:jmax )
    real(DP), intent(out  ) :: xy_RainLsc       ( 0:imax-1, 1:jmax )
    real(DP), intent(out  ) :: xy_SnowLsc       ( 0:imax-1, 1:jmax )


    real(DP) :: xyz_TempB     ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_QH2OLiqB  ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_QH2OSolB  ( 0:imax-1, 1:jmax, 1:kmax )

    real(DP) :: xyz_DelMass      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DTempDtFreeze( 0:imax-1, 1:jmax, 1:kmax )

    real(DP) :: xyz_WatFrac         ( 0:imax-1, 1:jmax, 1:kmax )

    real(DP) :: xyz_DQH2OLiqDt      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQH2OSolDt      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQH2OLiqDtCumLSC( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQH2OSolDtCumLSC( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQH2OLiqDtI2W   ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQH2OSolDtW2I   ( 0:imax-1, 1:jmax, 1:kmax )

    real(DP) :: xyz_DQRainDt  ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xyz_DQSnowDt  ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xy_Rain       ( 0:imax-1, 1:jmax )
    real(DP) :: xy_Snow       ( 0:imax-1, 1:jmax )
    real(DP) :: TempTentative
    real(DP) :: Melt

    integer  :: i
    integer  :: j
    integer  :: k


    ! 実行文 ; Executable statement
    !

    ! 初期化確認
    ! Initialization check
    !
    if ( .not. cloud_simple_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    ! save temperature, cloud water and ice before adjustment
    xyz_TempB    = xyz_Temp
    xyz_QH2OLiqB = xyz_QH2OLiq
    xyz_QH2OSolB = xyz_QH2OSol

    ! mass of an atmospheric layer
    do k = 1, kmax
      xyz_DelMass(:,:,k) = ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
    end do

    ! fraction of cloud water
    call CloudUtilsWatFraction(      &
      & xyz_Temp,                    & ! (in )
      & xyz_WatFrac                  & ! (out)
      & )

    ! partition cloud water production rate into water and ice
    xyz_DQH2OLiqDtCumLSC = &
      & ( xyz_DQH2OLiqDtCum + xyz_DQH2OLiqDtLSC ) * xyz_WatFrac
    xyz_DQH2OSolDtCumLSC = &
      & ( xyz_DQH2OLiqDtCum + xyz_DQH2OLiqDtLSC ) - xyz_DQH2OLiqDtCumLSC

    ! conversion rate from water to ice, and vice versa
    xyz_DQH2OLiqDtI2W = &
      & ( xyz_QH2OLiq + xyz_QH2OSol ) * xyz_WatFrac - xyz_QH2OLiq
    xyz_DQH2OLiqDtI2W = xyz_DQH2OLiqDtI2W / ( 2.0_DP * DelTime )
    xyz_DQH2OSolDtW2I = - xyz_DQH2OLiqDtI2W

    xyz_DQH2OLiqDt = xyz_DQH2OLiqDtCumLSC + xyz_DQH2OLiqDtI2W
    xyz_DQH2OSolDt = xyz_DQH2OSolDtCumLSC + xyz_DQH2OSolDtW2I

    ! temperature increase due to freezing
    xyz_DTempDtFreeze = LatentHeatFusion * xyz_DQH2OSolDt / CpDry
    xyz_Temp = xyz_Temp + xyz_DTempDtFreeze * 2.0_DP * DelTime

    xyz_QH2OLiq =                                                                 &
      &   xyz_QH2OLiq * exp( - 2.0_DP * DelTime / ( CloudWatLifeTime + 1.0d-100 ) )  &
      & + xyz_DQH2OLiqDt * CloudWatLifeTime                                          &
      &   * ( 1.0_DP - exp( - 2.0_DP * DelTime / ( CloudWatLifeTime + 1.0d-100 ) ) )

    xyz_DQRainDt =                                         &
      &   xyz_QH2OLiqB + xyz_DQH2OLiqDt * 2.0_DP * DelTime &
      & - xyz_QH2OLiq
    xyz_DQRainDt = xyz_DQRainDt / ( 2.0_DP * DelTime )


    xyz_QH2OSol =                                                                 &
      &   xyz_QH2OSol * exp( - 2.0_DP * DelTime / ( CloudIceLifeTime + 1.0d-100 ) )  &
      & + xyz_DQH2OSolDt * CloudIceLifeTime                                          &
      &   * ( 1.0_DP - exp( - 2.0_DP * DelTime / ( CloudIceLifeTime + 1.0d-100 ) ) )

    xyz_DQSnowDt =                                         &
      &   xyz_QH2OSolB + xyz_DQH2OSolDt * 2.0_DP * DelTime &
      & - xyz_QH2OSol
    xyz_DQSnowDt = xyz_DQSnowDt / ( 2.0_DP * DelTime )



    xy_Rain = 0.0_DP
    xy_Snow = 0.0_DP
    do k = kmax, 1, -1
      do j = 1, jmax
        do i = 0, imax-1
          xy_Rain(i,j) = xy_Rain(i,j)                               &
            & + xyz_DQRainDt(i,j,k) * xyz_DelMass(i,j,k)
          xy_Snow(i,j) = xy_Snow(i,j)                               &
            & + xyz_DQSnowDt(i,j,k) * xyz_DelMass(i,j,k)

          if ( xyz_Temp(i,j,k) > TempCondWater ) then
            Melt = xy_Snow(i,j)
            TempTentative = xyz_Temp(i,j,k)                  &
              & - LatentHeatFusion * Melt * 2.0_DP * DelTime &
              &   / ( CpDry * xyz_DelMass(i,j,k) )
            if ( TempTentative < TempCondWater ) then
              TempTentative = TempCondWater
              Melt =                                          &
                & CpDry * xyz_DelMass(i,j,k)                  &
                &   * ( xyz_Temp(i,j,k) - TempTentative )     &
                &   / ( LatentHeatFusion * 2.0_DP * DelTime )
            end if
            xyz_Temp(i,j,k) = TempTentative
            xy_Rain(i,j) = xy_Rain(i,j) + Melt
            xy_Snow(i,j) = xy_Snow(i,j) - Melt
          end if

        end do
      end do
    end do

    xy_RainCum = xy_Rain
    xy_SnowCum = xy_Snow
    xy_RainLsc = 0.0_DP
    xy_SnowLsc = 0.0_DP


    xyz_DTempDtCond = xyz_DTempDtCond &
      & + ( xyz_Temp - xyz_TempB ) / ( 2.0_DP * DelTime )


  end subroutine CloudSimpleWithIce

  !--------------------------------------------------------------------------------------

  subroutine CloudSimpleCalcCloudCover(  &
    & xyz_Press, xyz_Temp, xyz_QH2OTot,  & ! (in )
    & xyz_CloudCover                     & ! (out)
    & )

    ! USE statements
    !

    ! 飽和比湿の算出
    ! Evaluate saturation specific humidity
    !
    use saturate, only: xyz_CalcQVapSat

    real(DP), intent(in ) :: xyz_Press     ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in ) :: xyz_Temp      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in ) :: xyz_QH2OTot   ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(out) :: xyz_CloudCover( 0:imax-1, 1:jmax, 1:kmax )


    real(DP) :: xyz_RH(0:imax-1, 1:jmax, 1:kmax)


    ! 実行文 ; Executable statement
    !

    ! 初期化確認
    ! Initialization check
    !
    if ( .not. cloud_simple_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    select case ( IDCloudCoverMethod )
    case ( IDCloudCoverMethodConst )

      xyz_CloudCover = CloudCover

    case ( IDCloudCoverMethodRH )

      ! see Sundqvist et al. (1989), Del Genio et al. (1996)
      xyz_RH = xyz_QH2OTot / xyz_CalcQVapSat( xyz_Temp, xyz_Press )
      xyz_RH = min( xyz_RH, 1.0_DP )

      xyz_CloudCover = 1.0d0 - sqrt( ( 1.0d0 - xyz_RH ) / ( 1.0d0 - RHCrtl ) )

      xyz_CloudCover = max( xyz_CloudCover, 0.0_DP )
      xyz_CloudCover = min( xyz_CloudCover, 1.0_DP )

    end select


  end subroutine CloudSimpleCalcCloudCover

  !--------------------------------------------------------------------------------------

  subroutine CloudSimpleDivideWatAndIce(  &
    & xyz_Temp,                           & ! (in )
    & xyz_QH2OWatAndIce,                  & ! (in )
    & xyz_QH2OWat, xyz_QH2OIce            & ! (out)
    & )

    ! USE statements
    !

    ! 雲関系ルーチン
    ! Cloud-related routines
    !
    use cloud_utils, only : CloudUtilsWatFraction


    real(DP), intent(in ) :: xyz_Temp         (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in ) :: xyz_QH2OWatAndIce(0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(out) :: xyz_QH2OWat      (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(out) :: xyz_QH2OIce      (0:imax-1, 1:jmax, 1:kmax)


    real(DP) :: xyz_WatFrac(0:imax-1, 1:jmax, 1:kmax)

    ! 実行文 ; Executable statement
    !

    ! 初期化確認
    ! Initialization check
    !
    if ( .not. cloud_simple_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    call CloudUtilsWatFraction(      &
      & xyz_Temp,                    & ! (in )
      & xyz_WatFrac                  & ! (out)
      & )

    xyz_QH2OWat = xyz_QH2OWatAndIce * xyz_WatFrac
    xyz_QH2OIce = xyz_QH2OWatAndIce * ( 1.0_DP - xyz_WatFrac )


  end subroutine CloudSimpleDivideWatAndIce

  !--------------------------------------------------------------------------------------

  subroutine CloudSimpleCalcPRCPWithPC(    &
    & xyr_Press, xyz_Temp, xyz_DQH2OLiqDt, & ! (in )
    & xy_SurfRainFlux, xy_SurfSnowFlux,    & ! (out)
    & xyz_DTempDt                          & ! (out)
    & )

    ! 時刻管理
    ! Time control
    !
    use timeset, only: &
      & DelTime            ! $ \Delta t $ [s]

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only:  &
      & CpDry,            &
                              ! $ C_p $ [J kg-1 K-1].
                              ! 乾燥大気の定圧比熱.
                              ! Specific heat of air at constant pressure
      & Grav,             &
                              ! $ g $ [m s-2].
                              ! 重力加速度.
                              ! Gravitational acceleration
      & LatentHeatFusion
                              ! $ L $ [J kg-1] .
                              ! 融解の潜熱.
                              ! Latent heat of fusion

    ! 雲関系ルーチン
    ! Cloud-related routines
    !
    use cloud_utils, only : CloudUtilsWatFraction


    real(DP), intent(in ) :: xyr_Press      ( 0:imax-1, 1:jmax, 0:kmax )
    real(DP), intent(in ) :: xyz_Temp       ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in ) :: xyz_DQH2OLiqDt ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(out) :: xy_SurfRainFlux( 0:imax-1, 1:jmax )
    real(DP), intent(out) :: xy_SurfSnowFlux( 0:imax-1, 1:jmax )
    real(DP), intent(out) :: xyz_DTempDt    ( 0:imax-1, 1:jmax, 1:kmax )


    ! 作業変数
    ! Work variables
    !
    real(DP) :: xyz_WatFrac( 0:imax-1, 1:jmax, 1:kmax )
    real(DP) :: xy_Rain    ( 0:imax-1, 1:jmax )
    real(DP) :: xy_Snow    ( 0:imax-1, 1:jmax )
    real(DP) :: xy_Tot     ( 0:imax-1, 1:jmax )

    integer:: k


    ! 初期化確認
    ! Initialization check
    !
    if ( .not. cloud_simple_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    call CloudUtilsWatFraction(   &
      & xyz_Temp,                 & ! (in )
      & xyz_WatFrac               & ! (out)
      & )


    xy_SurfRainFlux = 0.0_DP
    xy_SurfSnowFlux = 0.0_DP

    do k = kmax, 1, -1
      xy_Rain = xy_SurfRainFlux   &
        & + xyz_DQH2OLiqDt(:,:,k) &
        &     * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
      xy_Snow = xy_SurfSnowFlux
      xy_Tot  = xy_Rain + xy_Snow

      xyz_DTempDt(:,:,k) =                                                &
        & + ( xy_Rain - xyz_WatFrac(:,:,k) * xy_Tot ) * 2.0_DP * DelTime  &
        &     / ( ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav )      &
        &     * LatentHeatFusion / CpDry                                  &
        & / ( 2.0_DP * DelTime )

      xy_SurfRainFlux = xyz_WatFrac(:,:,k)              * xy_Tot
      xy_SurfSnowFlux = ( 1.0_DP - xyz_WatFrac(:,:,k) ) * xy_Tot
    end do


  end subroutine CloudSimpleCalcPRCPWithPC

  !--------------------------------------------------------------------------------------

  subroutine CloudSimpleInit( &
    & ArgFlagSnow             &
    & )

    ! ファイル入出力補助
    ! File I/O support
    !
    use dc_iounit, only: FileOpen

    ! NAMELIST ファイル入力に関するユーティリティ
    ! Utilities for NAMELIST file input
    !
    use namelist_util, only: namelist_filename, NmlutilMsg, NmlutilAryValid

    ! ヒストリデータ出力
    ! History data output
    !
    use gtool_historyauto, only: HistoryAutoAddVariable

    ! 大規模凝結 (非対流性凝結)
    ! Large scale condensation (non-convective condensation)
    !
    use lscond, only : LScaleCondInit

    ! 雲関系ルーチン
    ! Cloud-related routines
    !
    use cloud_utils, only : CloudUtilsInit


    ! 宣言文 ; Declaration statements
    !

    logical, intent(in) :: ArgFlagSnow


    character(STRING) :: CloudCoverMethod

    integer:: unit_nml        ! NAMELIST ファイルオープン用装置番号.
                              ! Unit number for NAMELIST file open
    integer:: iostat_nml      ! NAMELIST 読み込み時の IOSTAT.
                              ! IOSTAT of NAMELIST read

    ! NAMELIST 変数群
    ! NAMELIST group name
    !
    namelist /cloud_simple_nml/ &
      & CloudLifeTime,       &
      & CloudWatLifeTime,    &
      & CloudIceLifeTime,    &
      & CloudCoverMethod,    &
      & RHCrtl,              &
      & CloudCover
          !
          ! デフォルト値については初期化手続 "cloud_simple#CloudSimpleInit"
          ! のソースコードを参照のこと.
          !
          ! Refer to source codes in the initialization procedure
          ! "cloud_simple#CloudSimpleInit" for the default values.
          !

    ! 実行文 ; Executable statement
    !

    if ( cloud_simple_inited ) return


    FlagSnow = ArgFlagSnow


    ! デフォルト値の設定
    ! Default values settings
    !
    CloudLifeTime       = 3600.0_DP
    CloudWatLifeTime    = 3600.0_DP
    CloudIceLifeTime    = 3600.0_DP

    CloudCoverMethod    = 'Const'
!!$    CloudCoverMethod    = 'RH'

    RHCrtl              = 0.8_DP

    CloudCover          = 1.0_DP


    ! NAMELIST の読み込み
    ! NAMELIST is input
    !
    if ( trim(namelist_filename) /= '' ) then
      call FileOpen( unit_nml, &          ! (out)
        & namelist_filename, mode = 'r' ) ! (in)

      rewind( unit_nml )
      read( unit_nml,                     & ! (in)
        & nml = cloud_simple_nml,         & ! (out)
        & iostat = iostat_nml )             ! (out)
      close( unit_nml )

      call NmlutilMsg( iostat_nml, module_name ) ! (in)
    end if


    select case ( CloudCoverMethod )
    case ( 'Const' )
      IDCloudCoverMethod = IDCloudCoverMethodConst
    case ( 'RH' )
      IDCloudCoverMethod = IDCloudCoverMethodRH
    case default
      call MessageNotify( 'E', module_name,         &
        & 'CloudCoverMethod=<%c> is not supported.', &
        & c1 = trim(CloudCoverMethod) )
    end select


    ! Initialization of modules used in this module
    !

    ! 大規模凝結 (非対流性凝結) (Manabe, 1965)
    ! Large scale condensation (non-convective condensation) (Le Treut and Li, 1991)
    !
    call LScaleCondInit

    ! 雲関系ルーチン
    ! Cloud-related routines
    !
    call CloudUtilsInit( &
      & FlagSnow         &
      & )


    ! ヒストリデータ出力のためのへの変数登録
    ! Register of variables for history data output
    !
!!$    call HistoryAutoAddVariable( 'EffCloudCover', &
!!$      & (/ 'lon ', 'lat ', 'time' /), &
!!$      & 'effective cloud cover', '1' )



    ! 印字 ; Print
    !
    call MessageNotify( 'M', module_name, '----- Initialization Messages -----' )
    call MessageNotify( 'M', module_name, 'CloudLifeTime       = %f', d = (/ CloudLifeTime /) )
    call MessageNotify( 'M', module_name, 'CloudWatLifeTime    = %f', d = (/ CloudWatLifeTime /) )
    call MessageNotify( 'M', module_name, 'CloudIceLifeTime    = %f', d = (/ CloudIceLifeTime /) )
    call MessageNotify( 'M', module_name, 'CloudCoverMethod    = %c', c1 = trim(CloudCoverMethod) )
    call MessageNotify( 'M', module_name, 'RHCrtl              = %f', d = (/ RHCrtl /) )
    call MessageNotify( 'M', module_name, 'CloudCover          = %f', d = (/ CloudCover /) )
    call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )


    cloud_simple_inited = .true.

  end subroutine CloudSimpleInit

  !--------------------------------------------------------------------------------------

end module cloud_simple
