Class set_cloud
In: radiation/set_cloud.f90

雲の設定

Set cloud

Note that Japanese and English are described in parallel.

雲の分布を設定.

In this module, the amount of cloud or cloud optical depth are set. This module is under development and is still a preliminary version.

Procedures List

!$ ! RadiationFluxDennouAGCM :放射フラックスの計算
!$ ! ———— :————
!$ ! RadiationFluxDennouAGCM :Calculate radiation flux

NAMELIST

NAMELIST#set_cloud_nml

Methods

Included Modules

dc_types dc_message gridset timeset constants saturate gtool_historyauto constants_snowseaice dc_iounit namelist_util

Public Instance methods

Subroutine :
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(in )
xyz_CloudCover( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(out)

[Source]

  subroutine SetCloudCalcCloudCover( xyz_Press, xyz_Temp, xyz_QH2OTot, xyz_CloudCover )

    ! 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. set_cloud_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    ! Cloud optical depth
    !
    if ( .not. FlagCloud ) then

      xyz_CloudCover = 0.0_DP

    else

      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 if


  end subroutine SetCloudCalcCloudCover
Subroutine :
xyz_TransCloudOneLayer(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyrr_OverlappedCloudTrans(0:imax-1, 1:jmax, 0:kmax, 0:kmax) :real(DP), intent(out)

[Source]

  subroutine SetCloudCalcOverlappedCloudTrans( xyz_TransCloudOneLayer, xyz_CloudCover, xyrr_OverlappedCloudTrans )

    ! USE statements
    !

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

    ! 時刻管理
    ! Time control
    !
    use timeset, only: TimeN, EndTime, TimesetClockStart, TimesetClockStop

!!$    use sort, only : SortQuick

    real(DP), intent(in ) :: xyz_TransCloudOneLayer   (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in ) :: xyz_CloudCover           (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(out) :: xyrr_OverlappedCloudTrans(0:imax-1, 1:jmax, 0:kmax, 0:kmax)


    real(DP) :: xyz_EffCloudCover           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_CloudCoverSorted        (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_EffCloudCoverSorted     (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_TransCloudOneLayerSorted(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: CloudCoverSortedCur
    real(DP) :: EffCloudCoverSortedCur
    real(DP) :: TransCloudOneLayerSortedCur
    integer  :: KInsPos
    integer  :: i
    integer  :: j
    integer  :: k
    integer  :: kk
    integer  :: kkk



    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth
    !
    if ( .not. FlagCloud ) then

      xyrr_OverlappedCloudTrans = 1.0_DP

    else

      select case ( IDCloudOverlapType )

        ! This is old version to be deleted. 

!!$      case ( IDCloudOverlapTypeFillGrid )
!!$        ! Cloud fills a grid box.
!!$
!!$        do k = 0, kmax
!!$          kk = k
!!$          xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0d0
!!$          do kk = k+1, kmax
!!$            xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,k,kk-1) &
!!$              & * xyz_TransCloudOneLayer(:,:,kk)
!!$          end do
!!$        end do
!!$        do k = 0, kmax
!!$          do kk = 0, k-1
!!$            xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,kk,k)
!!$          end do
!!$        end do

      case ( IDCloudOverlapTypeRandom )

        xyz_EffCloudCover = xyz_CloudCover * ( 1.0_DP - xyz_TransCloudOneLayer )

        do k = 0, kmax
          kk = k
          xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP
          do kk = k+1, kmax
            xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,k,kk-1) * ( 1.0_DP - xyz_EffCloudCover(:,:,kk) )
          end do
        end do

        do k = 0, kmax
          do kk = 0, k-1
            xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,kk,k)
          end do
        end do

      case ( IDCloudOverlapTypeMaxOverlap )

!!$        call MessageNotify( 'E', module_name, 'This has not been implemented, yet.' )


        ! see Chou et al. (2001)

        xyz_EffCloudCover = xyz_CloudCover * ( 1.0_DP - xyz_TransCloudOneLayer )


        ! Original method (computationally expensive, probably)
        !
!!$        do k = 0, kmax
!!$          kk = k
!!$          xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP
!!$          do kk = k+1, kmax
!!$
!!$            xyz_CloudCoverSorted         = xyz_CloudCover
!!$            xyz_EffCloudCoverSorted      = xyz_EffCloudCover
!!$            xyz_TransCloudOneLayerSorted = xyz_TransCloudOneLayer
!!$
!!$            call SortQuick( imax, jmax, kk-k,             &
!!$              & xyz_CloudCoverSorted        (:,:,k+1:kk), &
!!$              & xyz_EffCloudCoverSorted     (:,:,k+1:kk), &
!!$              & xyz_TransCloudOneLayerSorted(:,:,k+1:kk)  &
!!$              & )
!!$
!!$            xyrr_OverlappedCloudTrans(:,:,k,kk) = 0.0_DP
!!$            do kkk = k+1, kk
!!$              xyrr_OverlappedCloudTrans(:,:,k,kk) =         &
!!$                & xyz_EffCloudCoverSorted(:,:,kkk)          &
!!$                & + xyrr_OverlappedCloudTrans(:,:,k,kk)     &
!!$                &   * xyz_TransCloudOneLayerSorted(:,:,kkk)
!!$            end do
!!$            xyrr_OverlappedCloudTrans(:,:,k,kk) = &
!!$              & 1.0_DP - xyrr_OverlappedCloudTrans(:,:,k,kk)
!!$
!!$          end do
!!$        end do



        ! Economical method (probably)
        !
        do k = 0, kmax

!!$          do kkk = 1, kmax
!!$              xyz_CloudCoverSorted(:,:,kkk) = real( kmax-kkk ) / real(kmax)
!!$!              xyz_CloudCoverSorted(:,:,kkk) = abs( 0.55d0 - real( kmax-kkk ) / real(kmax) )
!!$          end do
!!$          ! debug output
!!$          if ( k == 0 ) then
!!$            kk = kmax
!!$            do kkk = k+1, kk
!!$              write( 6, * ) kkk, xyz_CloudCoverSorted(0,jmax/2+1,kkk)
!!$            end do
!!$          end if

          xyz_CloudCoverSorted         = xyz_CloudCover
          xyz_EffCloudCoverSorted      = xyz_EffCloudCover
          xyz_TransCloudOneLayerSorted = xyz_TransCloudOneLayer

          kk = k
          xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP
          do kk = k+1, kmax


            do j = 1, jmax
              do i = 0, imax-1

                ! xyz_CloudCoverSorted(i,j,kk) is inserved in an appropriate position.
                !
                KInsPos = kk
                loop : do kkk = k+1, kk-1

                  if ( xyz_CloudCoverSorted(i,j,kk) < xyz_CloudCoverSorted(i,j,kkk) ) then
                    KInsPos = kkk
                    exit loop
                  end if

                end do loop

                ! values are saved
                CloudCoverSortedCur         = xyz_CloudCoverSorted        (i,j,kk)
                EffCloudCoverSortedCur      = xyz_EffCloudCoverSorted     (i,j,kk)
                TransCloudOneLayerSortedCur = xyz_TransCloudOneLayerSorted(i,j,kk)

                ! values are shifted upward to empty an array at insert position
                do kkk = kk, KInsPos+1, -1
                  xyz_CloudCoverSorted        (i,j,kkk) = xyz_CloudCoverSorted        (i,j,kkk-1)
                  xyz_EffCloudCoverSorted     (i,j,kkk) = xyz_EffCloudCoverSorted     (i,j,kkk-1)
                  xyz_TransCloudOneLayerSorted(i,j,kkk) = xyz_TransCloudOneLayerSorted(i,j,kkk-1)
                end do
                kkk = KInsPos
                xyz_CloudCoverSorted        (i,j,kkk) = CloudCoverSortedCur
                xyz_EffCloudCoverSorted     (i,j,kkk) = EffCloudCoverSortedCur
                xyz_TransCloudOneLayerSorted(i,j,kkk) = TransCloudOneLayerSortedCur

              end do
            end do


!!$            xyz_CloudCoverSorted         = xyz_CloudCover
!!$            do kkk = 1, kmax
!!$              xyz_CloudCoverSorted(:,:,kkk) = real( kmax-kkk ) / real(kmax)
!!$            end do
!!$            xyz_EffCloudCoverSorted      = xyz_EffCloudCover
!!$            xyz_TransCloudOneLayerSorted = xyz_TransCloudOneLayer
!!$
!!$            call SortQuick( imax, jmax, kk-k,             &
!!$              & xyz_CloudCoverSorted        (:,:,k+1:kk), &
!!$              & xyz_EffCloudCoverSorted     (:,:,k+1:kk), &
!!$              & xyz_TransCloudOneLayerSorted(:,:,k+1:kk)  &
!!$              & )


!!$            ! debug output
!!$            if ( ( k == 0 ) .and. ( kk == kmax-2 ) ) then
!!$              do kkk = k+1, kk
!!$                write( 6, * ) kkk, xyz_CloudCoverSorted(0,jmax/2+1,kkk)
!!$              end do
!!$              write( 6, * ) '-----'
!!$            end if


            xyrr_OverlappedCloudTrans(:,:,k,kk) = 0.0_DP
            do kkk = k+1, kk
              xyrr_OverlappedCloudTrans(:,:,k,kk) = xyz_EffCloudCoverSorted(:,:,kkk) + xyrr_OverlappedCloudTrans(:,:,k,kk) * xyz_TransCloudOneLayerSorted(:,:,kkk)
            end do
            xyrr_OverlappedCloudTrans(:,:,k,kk) = 1.0_DP - xyrr_OverlappedCloudTrans(:,:,k,kk)

          end do
        end do



        do k = 0, kmax
          do kk = 0, k-1
            xyrr_OverlappedCloudTrans(:,:,k,kk) = xyrr_OverlappedCloudTrans(:,:,kk,k)
          end do
        end do


      end select

    end if


    ! Output effective cloud cover
    !
!!$    call HistoryAutoPut( TimeN, 'EffCloudCover', &
!!$      & 1.0_DP - xyrr_OverlappedCloudTrans(:,:,0,kmax) )


  end subroutine SetCloudCalcOverlappedCloudTrans
Subroutine :
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(in )
xy_SurfRainFlux( 0:imax-1, 1:jmax ) :real(DP), intent(out)
xy_SurfSnowFlux( 0:imax-1, 1:jmax ) :real(DP), intent(out)

[Source]

  subroutine SetCloudCalcPRCP( xyr_Press, xyz_Temp, xyz_DQH2OLiqDt, xy_SurfRainFlux, xy_SurfSnowFlux )

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: CpDry, Grav, LatentHeatFusion
                              ! $ L $ [J kg-1] .
                              ! 融解の潜熱.
                              ! Latent heat of fusion

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


    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 )


    ! 作業変数
    ! Work variables
    !
    real(DP) :: xy_PRCP( 0:imax-1, 1:jmax )
    real(DP) :: xy_TempIncByFusion( 0:imax-1, 1:jmax )
                              ! Temperature increase by fusion

    integer:: i               ! 経度方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in longitude
    integer:: j               ! 緯度方向に回る DO ループ用作業変数
                              ! Work variables for DO loop in latitude
    integer:: k


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



    xy_PRCP = 0.0d0
    do k = kmax, 1, -1
      xy_PRCP = xy_PRCP + xyz_DQH2OLiqDt(:,:,k) * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
    end do


    ! Precipitation is added to soil moisture or surface snow
    !
    do j = 1, jmax
      do i = 0, imax-1
        if ( xyz_Temp(i,j,1) > TempCondWater ) then
          xy_SurfRainFlux(i,j) = xy_PRCP(i,j)
          xy_SurfSnowFlux(i,j) = 0.0_DP
        else
          xy_SurfRainFlux(i,j) = 0.0_DP
          xy_SurfSnowFlux(i,j) = xy_PRCP(i,j)

!!$            ! Heating by latent heat release by fusion.
!!$            !   It is assumed that all of latent heat for fusion is used to heat 
!!$            !   the lowest layer. 
!!$            !   It should be checked whether this causes unstable layer or not. 
!!$            !   (yot, 2010/08/17)
!!$            !   This is commented out. (yot, 2010/08/21)
!!$            !
!!$            xy_TempIncByFusion(i,j) =                                     &
!!$              & xy_SurfPRCPFlux(i,j) * LatentHeatFusion * 2.0d0 * DelTime &
!!$              &   / ( CpDry * z_DelSigma(1) * xy_Ps(i,j) / Grav )
!!$
!!$!            xyz_Temp(i,j,1) = xyz_Temp(i,j,1) + xy_TempIncByFusion(i,j)
!!$
!!$!            if ( xy_TempIncByFusion(i,j) > 0.0d0 ) then
!!$!              write( 6, * ) i, j, xyz_Temp(i,j,1), xy_TempIncByFusion(i,j)
!!$!            end if


        end if
      end do
    end do


  end subroutine SetCloudCalcPRCP
Subroutine :
xyz_DQCloudWaterDtCum( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(in )
xyz_DQCloudWaterDtLSC( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(in )
xyz_QCloudWater( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(inout)

[Source]

  subroutine SetCloudCloudWaterAdjust( xyz_DQCloudWaterDtCum, xyz_DQCloudWaterDtLSC, xyz_QCloudWater )

    ! USE statements
    !

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


    real(DP), intent(in   ) :: xyz_DQCloudWaterDtCum( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in   ) :: xyz_DQCloudWaterDtLSC( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_QCloudWater      ( 0:imax-1, 1:jmax, 1:kmax )


    real(DP) :: xyz_QCloudWaterB( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xyz_PRCP        ( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xy_PRCP         ( 0:imax-1, 1:jmax )

    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth
    !
    if ( FlagCloud ) then

      ! save cloud water before adjustment
      xyz_QCloudWaterB = xyz_QCloudWater

!!$      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_QCloudWater = ( xyz_QCloudWater / ( 2.0_DP * DelTime ) + xyz_DQCloudWaterDtCum + xyz_DQCloudWaterDtLSC ) / ( 1.0_DP / ( 2.0_DP * DelTime ) + 1.0_DP / ( CloudLifeTime + 1.0d-100 ) )


!!$      do k = 1, kmax
!!$        xyz_PRCP(:,:,k) =                                                       &
!!$          & (                                                                   &
!!$          &     xyz_QCloudWaterB(:,:,k)                                         &
!!$          &   + ( xyz_DQCloudWaterDtCum(:,:,k) + xyz_DQCloudWaterDtLSC(:,:,k) ) &
!!$          &     * ( 2.0_DP * DelTime )                                          &
!!$          &   - xyz_QCloudWater(:,:,k)                                          &
!!$          & )                                                                   &
!!$          & * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
!!$      end do
!!$      xyz_PRCP = xyz_PRCP / ( 2.0_DP * DelTime )
!!$
!!$      xy_PRCP = 0.0_DP
!!$      do k = kmax, 1, -1
!!$        xy_PRCP = xy_PRCP + xyz_PRCP(:,:,k)
!!$      end do


    end if


  end subroutine SetCloudCloudWaterAdjust
Subroutine :
xyz_Press(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyr_Press(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
xyz_VirTemp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DQCloudWaterDtCum(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_MoistConvDetTend(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_OMG(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_MoistConvSubsidMassFlux(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_Temp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_QH2OVap(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_QCloudWater(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_Rain(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out )

[Source]

  subroutine SetCloudCloudWaterAdjustWithCloudCoverV2( xyz_Press, xyr_Press, xyz_VirTemp, xyz_DQCloudWaterDtCum, xyz_MoistConvDetTend, xyz_OMG, xyz_MoistConvSubsidMassFlux, xyz_DTempDtPhy, xyz_Temp, xyz_QH2OVap, xyz_QCloudWater, xyz_CloudCover, xyz_Rain )

    ! USE statements
    !

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

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: Grav, CpDry, GasRDry, GasRWet, LatentHeat
                              ! $ L $ [J kg-1] .
                              ! 凝結の潜熱.
                              ! Latent heat of condensation

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

    real(DP), intent(in   ) :: xyz_Press                  (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyr_Press                  (0:imax-1, 1:jmax, 0:kmax)
    real(DP), intent(in   ) :: xyz_VirTemp                (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_DQCloudWaterDtCum      (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_MoistConvDetTend       (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_OMG                    (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_MoistConvSubsidMassFlux(0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_DTempDtPhy             (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_QH2OVap                (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_QCloudWater            (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_CloudCover             (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(out  ) :: xyz_Rain                   (0:imax-1, 1:jmax, 1:kmax)


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

    real(DP) :: xyz_QH2OVapSat       (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDPress(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDTemp (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDt    (0:imax-1, 1:jmax, 1:kmax)

    real(DP), parameter :: MixCoef = 1.0d-6

    real(DP) :: xyz_ZeroOneDQsDt    (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DelCloudCoverStr(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA1          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA2          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactB           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC1          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC2          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactD           (0:imax-1, 1:jmax, 1:kmax)

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

!!$    real(DP) :: xyz_PRCP        ( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xy_PRCP         ( 0:imax-1, 1:jmax )

    real(DP), parameter :: QCloudWaterThreshold = 1.0d-10


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


    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth
    !
    if ( FlagCloud ) then


      ! Save cloud water amount
      !
      xyz_QCloudWaterB = xyz_QCloudWater


      xyz_QH2OVapSat       = xyz_CalcQVapSat( xyz_Temp, xyz_Press )
      xyz_DQH2OVapSatDTemp = xyz_CalcDQVapSatDTemp( xyz_Temp, xyz_QH2OVapSat )

      xyz_DQH2OVapSatDPress = xyz_QH2OVapSat / xyz_Press * ( LatentHeat * GasRDry * xyz_VirTemp - CpDry * GasRWet * xyz_Temp**2 ) / ( LatentHeat**2 * xyz_QH2OVapSat     + CpDry * GasRWet * xyz_Temp**2 )


      xyz_DQH2OVapSatDt = xyz_DQH2OVapSatDPress * ( xyz_OMG + Grav * xyz_MoistConvSubsidMassFlux ) + xyz_DQH2OVapSatDTemp * xyz_DTempDtPhy


      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_DQH2OVapSatDt(i,j,k) < 0.0_DP ) then
              xyz_ZeroOneDQsDt(i,j,k) = 1.0_DP
            else
              xyz_ZeroOneDQsDt(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do


!!$      do k = 1, kmax
!!$        do j = 1, jmax
!!$          do i = 0, imax-1
!!$            if ( xyz_QH2OVap(i,j,k) >= xyz_QH2OVapSat(i,j,k) * ( 1.0_DP - 1.0d-10 ) ) then
!!$              xyz_QH2OVapLim(i,j,k) = xyz_QH2OVap(i,j,k) * ( 1.0_DP - 1.0d-10 )
!!$            else
!!$              xyz_QH2OVapLim(i,j,k) = xyz_QH2OVap(i,j,k)
!!$            end if
!!$          end do
!!$        end do
!!$      end do


      xyz_FactC1 = xyz_MoistConvDetTend
!!$      xyz_FactC2 =                                             &
!!$        & - ( 1.0_DP - xyz_CloudCover )                        &
!!$        &   / ( 2.0_DP * ( xyz_QH2OVapSat - xyz_QH2OVapLim ) ) &
!!$        &   * xyz_DQH2OVapSatDt * xyz_ZeroOneDQsDt
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVap(i,j,k) < 0.99_DP * xyz_QH2OVapSat(i,j,k) ) then
              xyz_FactC2(i,j,k) = - ( 1.0_DP - xyz_CloudCover(i,j,k) ) / ( 2.0_DP * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) ) * xyz_DQH2OVapSatDt(i,j,k) * xyz_ZeroOneDQsDt(i,j,k)
            else
              xyz_FactC2(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      xyz_FactC = xyz_FactC1 + xyz_FactC2
!!$      xyz_FactD = xyz_CloudCover * MixCoef * ( xyz_QH2OVapSat - xyz_QH2OVap ) &
!!$        &  / ( xyz_QCloudWater + 1.0d-100 )
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVapSat(i,j,k) > xyz_QH2OVap(i,j,k) ) then
              xyz_FactD(i,j,k) = xyz_CloudCover(i,j,k) * MixCoef * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) / ( xyz_QCloudWater(i,j,k) + 1.0d-100 )
            else
              xyz_FactD(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      !
      xyz_DelCloudCoverStr = xyz_FactC2 * 2.0_DP * DelTime - xyz_FactC2 / ( xyz_FactC + xyz_FactD + 1.0d-100 ) * ( xyz_FactC * 2.0_DP * DelTime + ( xyz_CloudCover - xyz_FactC / ( xyz_FactC + xyz_FactD + 1.0d-100 ) ) * ( 1.0_DP - exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) ) )
      !
      xyz_FactA1 = xyz_DQCloudWaterDtCum
      xyz_FactA2 = - xyz_CloudCover * xyz_DQH2OVapSatDt - xyz_DelCloudCoverStr * xyz_DQH2OVapSatDt * xyz_ZeroOneDQsDt - xyz_CloudCover * MixCoef * ( xyz_QH2OVapSat - xyz_QH2OVap )
      !    The value of xyz_FactA2 is checked, and is updated.
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_FactA2(i,j,k) * 2.0_DP * DelTime > xyz_QH2OVap(i,j,k) ) then
              xyz_FactA2(i,j,k) = xyz_QH2OVap(i,j,k) / ( 2.0_DP * DelTime )
!!$              xyz_FactA2(i,j,k) = ( xyz_QH2OVap(i,j,k) - 1.0d-100 ) / ( 2.0_DP * DelTime )
            end if
          end do
        end do
      end do
      xyz_FactA = xyz_FactA1 + xyz_FactA2
      !
      xyz_FactB = 1.0_DP / ( CloudLifeTime + 1.0d-100 )



      ! Values at next time step are calculated.
      !
      xyz_QCloudWater = xyz_QCloudWater * exp( - xyz_FactB * 2.0_DP * DelTime ) + xyz_FactA / xyz_FactB * ( 1.0_DP - exp( - xyz_FactB * 2.0_DP * DelTime ) )
      !   The value of cloud water amount is checked, and xyz_FactA is updated.
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < 0.0_DP ) then
              xyz_FactA2(i,j,k) = - xyz_FactA1(i,j,k) - xyz_FactB(i,j,k) * xyz_QCloudWaterB(i,j,k) * exp( - xyz_FactB(i,j,k) * 2.0_DP * DelTime ) / ( 1.0_DP - exp( - xyz_FactB(i,j,k) * 2.0_DP * DelTime ) )
              xyz_QCloudWater(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      xyz_FactA = xyz_FactA1 + xyz_FactA2
      !
      xyz_CloudCover = xyz_CloudCover * exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) + xyz_FactC / ( xyz_FactC + xyz_FactD + 1.0d-100 ) * ( 1.0_DP - exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) )
      !
      xyz_QH2OVap = xyz_QH2OVap - xyz_FactA2 * 2.0_DP * DelTime
      !
      xyz_Temp = xyz_Temp + xyz_FactA2 * 2.0_DP * DelTime * LatentHeat / CpDry


      ! Rain
      !
      xyz_Rain = xyz_FactA * 2.0_DP * DelTime + xyz_QCloudWaterB * ( 1.0_DP - exp( - xyz_FactB * 2.0_DP * DelTime ) ) - xyz_FactA / xyz_FactB * ( 1.0_DP - exp( - xyz_FactB * 2.0_DP * DelTime ) )
      xyz_Rain = xyz_Rain / ( 2.0_DP * DelTime )
!!$      do k = 1, kmax
!!$        xyz_Rain(:,:,k) = xyz_Rain(:,:,k) &
!!$          & * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
!!$      end do



      ! Evaporation
      !
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < QCloudWaterThreshold ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do

      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_CloudCover(i,j,k) > 1.0_DP ) then
              xyz_CloudCover(i,j,k) = 1.0_DP
            else if ( xyz_CloudCover(i,j,k) < 0.0_DP ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do


      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < 0.0_DP ) then
              write( 6, * ) i, j, k, xyz_QCloudWater(i,j,k)
            end if
          end do
        end do
      end do


    end if


  end subroutine SetCloudCloudWaterAdjustWithCloudCoverV2
Subroutine :
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(in )
xyz_QH2OWat(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)
xyz_QH2OIce(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out)

[Source]

  subroutine SetCloudDivideWatAndIce( xyz_Temp, xyz_QH2OWatAndIce, xyz_QH2OWat, xyz_QH2OIce )

    ! USE statements
    !

    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. set_cloud_inited ) then
      call MessageNotify( 'E', module_name, 'This module has not been initialized.' )
    end if


    xyz_WatFrac = ( 1.0_DP - 0.0_DP ) / ( TempWatLim - TempIceLim ) * ( xyz_Temp - TempIceLim )
    xyz_WatFrac = min( xyz_WatFrac, 1.0_DP )
    xyz_WatFrac = max( xyz_WatFrac, 0.0_DP )

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


  end subroutine SetCloudDivideWatAndIce
Subroutine :

This procedure input/output NAMELIST#set_cloud_nml .

[Source]

  subroutine SetCloudInit

    ! ファイル入出力補助
    ! 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


    ! 宣言文 ; Declaration statements
    !
    character(STRING) :: CloudOverlapType
    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 /set_cloud_nml/ FlagCloud, CloudLifeTime, CloudOverlapType, CloudCoverMethod, RHCrtl, CloudCover, TempWatLim, TempIceLim
          !
          ! デフォルト値については初期化手続 "set_cloud#setCloudInit"
          ! のソースコードを参照のこと.
          !
          ! Refer to source codes in the initialization procedure
          ! "set_cloud#SetCloudInit" for the default values.
          !

    ! 実行文 ; Executable statement
    !

    if ( set_cloud_inited ) return


    ! デフォルト値の設定
    ! Default values settings
    !

    FlagCloud           = .true.

    CloudLifeTime       = 3600.0_DP

    CloudOverlapType    = "Random"
!!$    CloudOverlapType    = "MaxOverlap"

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

    RHCrtl              = 0.8_DP

    CloudCover          = 1.0_DP

    TempWatLim          = 273.15_DP
    TempIceLim          = 273.15_DP - 40.0_DP


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

      rewind( unit_nml )
      read( unit_nml, nml = set_cloud_nml, iostat = iostat_nml )             ! (out)
      close( unit_nml )

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


    select case ( CloudOverlapType )
    case ( 'Random' )
      IDCloudOverlapType = IDCloudOverlapTypeRandom
    case ( 'MaxOverlap' )
      IDCloudOverlapType = IDCloudOverlapTypeMaxOverlap
    case default
      call MessageNotify( 'E', module_name, 'CloudOverlapType=<%c> is not supported.', c1 = trim(CloudOverlapType) )
    end select

    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


    ! ヒストリデータ出力のためのへの変数登録
    ! 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, 'FlagCloud           = %b', l = (/ FlagCloud /) )
    call MessageNotify( 'M', module_name, 'CloudLifeTime       = %f', d = (/ CloudLifeTime /) )
    call MessageNotify( 'M', module_name, 'CloudOverlapType    = %c', c1 = trim(CloudOverlapType) )
    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, 'TempWatLim          = %f', d = (/ TempWatLim /) )
    call MessageNotify( 'M', module_name, 'TempIceLim          = %f', d = (/ TempIceLim /) )
    call MessageNotify( 'M', module_name, '-- version = %c', c1 = trim(version) )


    set_cloud_inited = .true.

  end subroutine SetCloudInit
Subroutine :
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)

[Source]

  subroutine SetCloudLocalizeCloud( xyz_CloudCover, xyz_DelCloudOptDep )

    ! USE statements
    !

    real(DP), intent(in   ) :: xyz_CloudCover    (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax)


    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth is scaled by considering cloud cover less than 1. 

    xyz_DelCloudOptDep = xyz_DelCloudOptDep / max( xyz_CloudCover, 1.0d-3 )


  end subroutine SetCloudLocalizeCloud

Private Instance methods

CloudCover
Variable :
CloudCover :real(DP), save
CloudLifeTime
Variable :
CloudLifeTime :real(DP), save
FlagCloud
Variable :
FlagCloud :logical , save
: A flag for cloud set. If FlagCloud is true, the effect of cloud is considered.
IDCloudCoverMethod
Variable :
IDCloudCoverMethod :integer , save
IDCloudCoverMethodConst
Constant :
IDCloudCoverMethodConst = 1 :integer , parameter
IDCloudCoverMethodRH
Constant :
IDCloudCoverMethodRH = 2 :integer , parameter
IDCloudOverlapType
Variable :
IDCloudOverlapType :integer , save
IDCloudOverlapTypeMaxOverlap
Constant :
IDCloudOverlapTypeMaxOverlap = 2 :integer , parameter
IDCloudOverlapTypeRandom
Constant :
IDCloudOverlapTypeRandom = 1 :integer , parameter
RHCrtl
Variable :
RHCrtl :real(DP), save
Subroutine :
xyz_DQCloudWaterDtCum( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(in )
xyz_DQCloudWaterDtLSC( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(in )
xyz_MoistConvDetTend( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(in )
xyz_QCloudWater( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(inout)
xyz_CloudCover( 0:imax-1, 1:jmax, 1:kmax ) :real(DP), intent(inout)

[Source]

  subroutine SetCloudCloudWaterAdjustWithCloudCoverV1( xyz_DQCloudWaterDtCum, xyz_DQCloudWaterDtLSC, xyz_MoistConvDetTend, xyz_QCloudWater, xyz_CloudCover )

    ! USE statements
    !

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

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: GasRDry
                              ! $ R $ [J kg-1 K-1].
                              ! 乾燥大気の気体定数.
                              ! Gas constant of air

    real(DP), intent(in   ) :: xyz_DQCloudWaterDtCum( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in   ) :: xyz_DQCloudWaterDtLSC( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(in   ) :: xyz_MoistConvDetTend ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_QCloudWater      ( 0:imax-1, 1:jmax, 1:kmax )
    real(DP), intent(inout) :: xyz_CloudCover       ( 0:imax-1, 1:jmax, 1:kmax )


    real(DP) :: xyz_QCloudWaterB( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xyz_PRCP        ( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xy_PRCP         ( 0:imax-1, 1:jmax )

    real(DP), parameter :: QCloudWaterThreshold = 1.0d-10

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


    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth
    !
    if ( FlagCloud ) then

      ! save cloud water before adjustment
      xyz_QCloudWaterB = xyz_QCloudWater

!!$      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_QCloudWater = ( xyz_QCloudWater / ( 2.0_DP * DelTime ) + xyz_DQCloudWaterDtCum + xyz_DQCloudWaterDtLSC ) / ( 1.0_DP / ( 2.0_DP * DelTime ) + 1.0_DP / ( CloudLifeTime + 1.0d-100 ) )


!!$      do k = 1, kmax
!!$        xyz_PRCP(:,:,k) =                                                       &
!!$          & (                                                                   &
!!$          &     xyz_QCloudWaterB(:,:,k)                                         &
!!$          &   + ( xyz_DQCloudWaterDtCum(:,:,k) + xyz_DQCloudWaterDtLSC(:,:,k) ) &
!!$          &     * ( 2.0_DP * DelTime )                                          &
!!$          &   - xyz_QCloudWater(:,:,k)                                          &
!!$          & )                                                                   &
!!$          & * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
!!$      end do
!!$      xyz_PRCP = xyz_PRCP / ( 2.0_DP * DelTime )
!!$
!!$      xy_PRCP = 0.0_DP
!!$      do k = kmax, 1, -1
!!$        xy_PRCP = xy_PRCP + xyz_PRCP(:,:,k)
!!$      end do


      ! a2 = a1 + ( 1 - a2 ) * P * 2 * dt
      !    = a1 + P * 2 * dt - a2 * P * 2 * dt
      ! ( 1 + P * 2 * dt ) * a2 = a1 + P * 2 * dt
      ! a2 = ( a1 + P * 2 * dt ) / ( 1 + P * 2 * dt )

!!$      xyz_CloudCover =                                                 &
!!$        & ( xyz_CloudCover + xyz_MoistConvDetTend * 2.0_DP * DelTime ) &
!!$        & / ( 1.0_DP + xyz_MoistConvDetTend * 2.0_DP * DelTime )

      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_DQCloudWaterDtLSC(i,j,k) > 0.0_DP ) then
              xyz_CloudCover(i,j,k) = 1.0_DP
            else
              xyz_CloudCover(i,j,k) = (   xyz_CloudCover(i,j,k) + xyz_MoistConvDetTend(i,j,k) * 2.0_DP * DelTime ) / ( 1.0_DP + xyz_MoistConvDetTend(i,j,k) * 2.0_DP * DelTime )
            end if
          end do
        end do
      end do

      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < QCloudWaterThreshold ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do

      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_CloudCover(i,j,k) > 1.0_DP ) then
              xyz_CloudCover(i,j,k) = 1.0_DP
            end if
          end do
        end do
      end do

    end if


  end subroutine SetCloudCloudWaterAdjustWithCloudCoverV1
Subroutine :
xyz_Press(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyr_Press(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
xyz_VirTemp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DQCloudWaterDtCum(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_MoistConvDetTend(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_OMG(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_MoistConvSubsidMassFlux(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_Temp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_QH2OVap(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_QCloudWater(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_Rain(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out )

[Source]

  subroutine SetCloudCloudWaterAdjustWithCloudCoverV3( xyz_Press, xyr_Press, xyz_VirTemp, xyz_DQCloudWaterDtCum, xyz_MoistConvDetTend, xyz_OMG, xyz_MoistConvSubsidMassFlux, xyz_DTempDtPhy, xyz_Temp, xyz_QH2OVap, xyz_QCloudWater, xyz_CloudCover, xyz_Rain )

    ! USE statements
    !

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

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: Grav, CpDry, GasRDry, GasRWet, LatentHeat
                              ! $ L $ [J kg-1] .
                              ! 凝結の潜熱.
                              ! Latent heat of condensation

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

    real(DP), intent(in   ) :: xyz_Press                  (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyr_Press                  (0:imax-1, 1:jmax, 0:kmax)
    real(DP), intent(in   ) :: xyz_VirTemp                (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_DQCloudWaterDtCum      (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_MoistConvDetTend       (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_OMG                    (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_MoistConvSubsidMassFlux(0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_DTempDtPhy             (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_QH2OVap                (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_QCloudWater            (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_CloudCover             (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(out  ) :: xyz_Rain                   (0:imax-1, 1:jmax, 1:kmax)


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

    real(DP) :: xyz_QH2OVapSat       (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDPress(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDTemp (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDt    (0:imax-1, 1:jmax, 1:kmax)

    real(DP), parameter :: MixCoef = 1.0d-6

    real(DP) :: xyz_ZeroOneDQsDt    (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DelCloudCoverStr(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA1          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA2          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactB           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC1          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC2          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactD           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactE           (0:imax-1, 1:jmax, 1:kmax)

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

!!$    real(DP) :: xyz_PRCP        ( 0:imax-1, 1:jmax, 1:kmax )
!!$    real(DP) :: xy_PRCP         ( 0:imax-1, 1:jmax )

    real(DP), parameter :: QCloudWaterThreshold = 1.0d-10

    real(DP), parameter :: RHThreshold = 0.99_DP


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


    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth
    !
    if ( FlagCloud ) then


      ! Save cloud water amount
      !
      xyz_QCloudWaterB = xyz_QCloudWater


      xyz_QH2OVapSat       = xyz_CalcQVapSat( xyz_Temp, xyz_Press )
      xyz_DQH2OVapSatDTemp = xyz_CalcDQVapSatDTemp( xyz_Temp, xyz_QH2OVapSat )

!!$      xyz_DQH2OVapSatDPress = xyz_QH2OVapSat / xyz_Press                           &
!!$        & * ( LatentHeat * GasRDry * xyz_VirTemp - CpDry * GasRWet * xyz_Temp**2 ) &
!!$        & / ( LatentHeat**2 * xyz_QH2OVapSat     + CpDry * GasRWet * xyz_Temp**2 )
      xyz_DQH2OVapSatDPress = xyz_QH2OVapSat / xyz_Press * ( LatentHeat * GasRDry * xyz_VirTemp - CpDry * GasRWet * xyz_Temp**2 ) / ( xyz_CloudCover * LatentHeat**2 * xyz_QH2OVapSat + CpDry * GasRWet * xyz_Temp**2 )

!!$      xyz_DQH2OVapSatDt =                                                            &
!!$        &   xyz_DQH2OVapSatDPress * ( xyz_OMG + Grav * xyz_MoistConvSubsidMassFlux ) &
!!$        & + xyz_DQH2OVapSatDTemp * xyz_DTempDtPhy
      xyz_DQH2OVapSatDt = xyz_DQH2OVapSatDPress * ( xyz_OMG + Grav * xyz_MoistConvSubsidMassFlux ) + xyz_DQH2OVapSatDTemp / ( 1.0_DP + xyz_CloudCover * LatentHeat / CpDry * xyz_DQH2OVapSatDTemp ) * xyz_DTempDtPhy


      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_DQH2OVapSatDt(i,j,k) < 0.0_DP ) then
              xyz_ZeroOneDQsDt(i,j,k) = 1.0_DP
            else
              xyz_ZeroOneDQsDt(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do


!!$      do k = 1, kmax
!!$        do j = 1, jmax
!!$          do i = 0, imax-1
!!$            if ( xyz_QH2OVap(i,j,k) >= xyz_QH2OVapSat(i,j,k) * ( 1.0_DP - 1.0d-10 ) ) then
!!$              xyz_QH2OVapLim(i,j,k) = xyz_QH2OVap(i,j,k) * ( 1.0_DP - 1.0d-10 )
!!$            else
!!$              xyz_QH2OVapLim(i,j,k) = xyz_QH2OVap(i,j,k)
!!$            end if
!!$          end do
!!$        end do
!!$      end do


      xyz_FactC1 = xyz_MoistConvDetTend
!!$      xyz_FactC2 =                                             &
!!$        & - ( 1.0_DP - xyz_CloudCover )                        &
!!$        &   / ( 2.0_DP * ( xyz_QH2OVapSat - xyz_QH2OVapLim ) ) &
!!$        &   * xyz_DQH2OVapSatDt * xyz_ZeroOneDQsDt
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVap(i,j,k) < RHThreshold * xyz_QH2OVapSat(i,j,k) ) then
              xyz_FactC2(i,j,k) = - ( 1.0_DP - xyz_CloudCover(i,j,k) ) / ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) * xyz_DQH2OVapSatDt(i,j,k) * xyz_ZeroOneDQsDt(i,j,k)
            else
              xyz_FactC2(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      xyz_FactC = xyz_FactC1 + xyz_FactC2
!!$      xyz_FactD = xyz_CloudCover * MixCoef * ( xyz_QH2OVapSat - xyz_QH2OVap ) &
!!$        &  / ( xyz_QCloudWater + 1.0d-100 )
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVapSat(i,j,k) > xyz_QH2OVap(i,j,k) ) then
              xyz_FactD(i,j,k) = 2.0_DP * xyz_CloudCover(i,j,k) * MixCoef * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) / ( xyz_QCloudWater(i,j,k) + 1.0d-100 )
            else
              xyz_FactD(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      !
      xyz_DelCloudCoverStr = xyz_FactC2 * 2.0_DP * DelTime - xyz_FactC2 / ( xyz_FactC + xyz_FactD + 1.0d-100 ) * ( xyz_FactC * 2.0_DP * DelTime + ( xyz_CloudCover - xyz_FactC / ( xyz_FactC + xyz_FactD + 1.0d-100 ) ) * ( 1.0_DP - exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) ) )
      !
      xyz_FactA1 = xyz_DQCloudWaterDtCum
      xyz_FactA2 = - xyz_CloudCover * xyz_DQH2OVapSatDt - xyz_DelCloudCoverStr * xyz_DQH2OVapSatDt * xyz_ZeroOneDQsDt - xyz_CloudCover * MixCoef * ( xyz_QH2OVapSat - xyz_QH2OVap )
      !    The value of xyz_FactA2 is checked, and is updated.
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_FactA2(i,j,k) * 2.0_DP * DelTime > xyz_QH2OVap(i,j,k) ) then
              xyz_FactA2(i,j,k) = xyz_QH2OVap(i,j,k) / ( 2.0_DP * DelTime )
!!$              xyz_FactA2(i,j,k) = ( xyz_QH2OVap(i,j,k) - 1.0d-100 ) / ( 2.0_DP * DelTime )
            end if
          end do
        end do
      end do
      xyz_FactA = xyz_FactA1 + xyz_FactA2
      !
      xyz_FactB = 1.0_DP / ( CloudLifeTime + 1.0d-100 )
      !
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVap(i,j,k) < RHThreshold * xyz_QH2OVapSat(i,j,k) ) then
              xyz_FactE(i,j,k) = ( 1.0_DP - xyz_CloudCover(i,j,k) )**2 / ( 2.0_DP * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) ) * xyz_DQH2OVapSatDt(i,j,k) * xyz_ZeroOneDQsDt(i,j,k)
            else
              xyz_FactE(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVapSat(i,j,k) > xyz_QH2OVap(i,j,k) ) then
              xyz_FactE(i,j,k) = xyz_FactE(i,j,k) + xyz_CloudCover(i,j,k)**2 * MixCoef * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) / ( xyz_QCloudWater(i,j,k) + 1.0d-100 )
            end if
          end do
        end do
      end do


      ! Values at next time step are calculated.
      !
      xyz_QCloudWater = xyz_QCloudWater * exp( - xyz_FactB * 2.0_DP * DelTime ) + xyz_FactA / xyz_FactB * ( 1.0_DP - exp( - xyz_FactB * 2.0_DP * DelTime ) )
      !   The value of cloud water amount is checked, and xyz_FactA is updated.
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < 0.0_DP ) then
              xyz_FactA2(i,j,k) = - xyz_FactA1(i,j,k) - xyz_FactB(i,j,k) * xyz_QCloudWaterB(i,j,k) * exp( - xyz_FactB(i,j,k) * 2.0_DP * DelTime ) / ( 1.0_DP - exp( - xyz_FactB(i,j,k) * 2.0_DP * DelTime ) )
              xyz_QCloudWater(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      xyz_FactA = xyz_FactA1 + xyz_FactA2
      !
      xyz_CloudCover = xyz_CloudCover * exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) + ( xyz_FactC + xyz_FactE ) / ( xyz_FactC + xyz_FactD + 1.0d-100 ) * ( 1.0_DP - exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) )
      !
      xyz_QH2OVap = xyz_QH2OVap - xyz_FactA2 * 2.0_DP * DelTime
      !
      xyz_Temp = xyz_Temp + xyz_FactA2 * 2.0_DP * DelTime * LatentHeat / CpDry


      ! Rain
      !
      xyz_Rain = xyz_FactA * 2.0_DP * DelTime + xyz_QCloudWaterB * ( 1.0_DP - exp( - xyz_FactB * 2.0_DP * DelTime ) ) - xyz_FactA / xyz_FactB * ( 1.0_DP - exp( - xyz_FactB * 2.0_DP * DelTime ) )
      xyz_Rain = xyz_Rain / ( 2.0_DP * DelTime )
!!$      do k = 1, kmax
!!$        xyz_Rain(:,:,k) = xyz_Rain(:,:,k) &
!!$          & * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav
!!$      end do



      ! Evaporation
      !
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < QCloudWaterThreshold ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do

      ! Cloud cover is restricted.
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_CloudCover(i,j,k) > 1.0_DP ) then
              xyz_CloudCover(i,j,k) = 1.0_DP
            else if ( xyz_CloudCover(i,j,k) < 0.0_DP ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do


      ! Check values
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < 0.0_DP ) then
              write( 6, * ) i, j, k, xyz_QCloudWater(i,j,k)
            end if
          end do
        end do
      end do


    end if


  end subroutine SetCloudCloudWaterAdjustWithCloudCoverV3
Subroutine :
xyz_Press(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyr_Press(0:imax-1, 1:jmax, 0:kmax) :real(DP), intent(in )
xyz_VirTemp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DQCloudWaterDtCum(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_MoistConvDetTend(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_OMG(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_MoistConvSubsidMassFlux(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DTempDtPhy(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_Temp(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_QH2OVap(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_QCloudWater(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)
xyz_Rain(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(out )

[Source]

  subroutine SetCloudCloudWaterAdjustWithCloudCoverV4( xyz_Press, xyr_Press, xyz_VirTemp, xyz_DQCloudWaterDtCum, xyz_MoistConvDetTend, xyz_OMG, xyz_MoistConvSubsidMassFlux, xyz_DTempDtPhy, xyz_Temp, xyz_QH2OVap, xyz_QCloudWater, xyz_CloudCover, xyz_Rain )

    ! USE statements
    !

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

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: Grav, CpDry, GasRDry, GasRWet, LatentHeat
                              ! $ L $ [J kg-1] .
                              ! 凝結の潜熱.
                              ! Latent heat of condensation

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

    real(DP), intent(in   ) :: xyz_Press                  (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyr_Press                  (0:imax-1, 1:jmax, 0:kmax)
    real(DP), intent(in   ) :: xyz_VirTemp                (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_DQCloudWaterDtCum      (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_MoistConvDetTend       (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_OMG                    (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_MoistConvSubsidMassFlux(0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(in   ) :: xyz_DTempDtPhy             (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_QH2OVap                (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_QCloudWater            (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_CloudCover             (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(out  ) :: xyz_Rain                   (0:imax-1, 1:jmax, 1:kmax)


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

    real(DP) :: xyz_QH2OVapSat       (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDPress(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDTemp (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DQH2OVapSatDt    (0:imax-1, 1:jmax, 1:kmax)

    real(DP) :: xyz_ZeroOneDQsDt    (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_DelCloudCoverStr(0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA1          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactA2          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactB           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC1          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactC2          (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactD           (0:imax-1, 1:jmax, 1:kmax)
    real(DP) :: xyz_FactE           (0:imax-1, 1:jmax, 1:kmax)

    real(DP) :: xy_Rain                    (0:imax-1, 1:jmax)
    real(DP) :: xy_RainConvFactor          (0:imax-1, 1:jmax)
    real(DP) :: xy_FactCo                  (0:imax-1, 1:jmax)
    real(DP) :: xy_FactBF                  (0:imax-1, 1:jmax)
    real(DP) :: xy_QCloudWaterConvThreshold(0:imax-1, 1:jmax)


    real(DP), parameter :: MixCoef                  = 1.0d-6
    real(DP), parameter :: QCloudWaterEvapThreshold = 1.0d-10
    real(DP), parameter :: RHThreshold              = 0.99_DP
    ! Values below are obtained from Sundqvist et al. (1989), but some of 
    ! those are arbitrarily selected.
    real(DP), parameter :: C0                       = 1.0d-3
    real(DP), parameter :: C1                       = 300.0_DP
    real(DP), parameter :: C2                       = 0.5_DP
    real(DP), parameter :: QCloudWaterConvThreshold = 4.0d-4

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


    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth
    !
    if ( FlagCloud ) then


      ! Save cloud water amount
      !
      xyz_QCloudWaterB = xyz_QCloudWater


      xyz_QH2OVapSat       = xyz_CalcQVapSat( xyz_Temp, xyz_Press )
      xyz_DQH2OVapSatDTemp = xyz_CalcDQVapSatDTemp( xyz_Temp, xyz_QH2OVapSat )

      xyz_DQH2OVapSatDPress = xyz_QH2OVapSat / xyz_Press * ( LatentHeat * GasRDry * xyz_VirTemp - CpDry * GasRWet * xyz_Temp**2 ) / ( xyz_CloudCover * LatentHeat**2 * xyz_QH2OVapSat + CpDry * GasRWet * xyz_Temp**2 )

      xyz_DQH2OVapSatDt = xyz_DQH2OVapSatDPress * ( xyz_OMG + Grav * xyz_MoistConvSubsidMassFlux ) + xyz_DQH2OVapSatDTemp / ( 1.0_DP + xyz_CloudCover * LatentHeat / CpDry * xyz_DQH2OVapSatDTemp ) * xyz_DTempDtPhy


      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_DQH2OVapSatDt(i,j,k) < 0.0_DP ) then
              xyz_ZeroOneDQsDt(i,j,k) = 1.0_DP
            else
              xyz_ZeroOneDQsDt(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do



      xyz_FactC1 = xyz_MoistConvDetTend
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVap(i,j,k) < RHThreshold * xyz_QH2OVapSat(i,j,k) ) then
              xyz_FactC2(i,j,k) = - ( 1.0_DP - xyz_CloudCover(i,j,k) ) / ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) * xyz_DQH2OVapSatDt(i,j,k) * xyz_ZeroOneDQsDt(i,j,k)
            else
              xyz_FactC2(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      xyz_FactC = xyz_FactC1 + xyz_FactC2
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVapSat(i,j,k) > xyz_QH2OVap(i,j,k) ) then
              xyz_FactD(i,j,k) = 2.0_DP * xyz_CloudCover(i,j,k) * MixCoef * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) / ( xyz_QCloudWater(i,j,k) + 1.0d-100 )
            else
              xyz_FactD(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      !
      xyz_DelCloudCoverStr = xyz_FactC2 * 2.0_DP * DelTime - xyz_FactC2 / ( xyz_FactC + xyz_FactD + 1.0d-100 ) * ( xyz_FactC * 2.0_DP * DelTime + ( xyz_CloudCover - xyz_FactC / ( xyz_FactC + xyz_FactD + 1.0d-100 ) ) * ( 1.0_DP - exp( - ( xyz_FactC + xyz_FactD ) * 2.0_DP * DelTime ) ) )
      !
      xyz_FactA1 = xyz_DQCloudWaterDtCum
      xyz_FactA2 = - xyz_CloudCover * xyz_DQH2OVapSatDt - xyz_DelCloudCoverStr * xyz_DQH2OVapSatDt * xyz_ZeroOneDQsDt - xyz_CloudCover * MixCoef * ( xyz_QH2OVapSat - xyz_QH2OVap )
      !    The value of xyz_FactA2 is checked, and is updated.
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_FactA2(i,j,k) * 2.0_DP * DelTime > xyz_QH2OVap(i,j,k) ) then
              xyz_FactA2(i,j,k) = xyz_QH2OVap(i,j,k) / ( 2.0_DP * DelTime )
            end if
          end do
        end do
      end do
      xyz_FactA = xyz_FactA1 + xyz_FactA2
      !
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVap(i,j,k) < RHThreshold * xyz_QH2OVapSat(i,j,k) ) then
              xyz_FactE(i,j,k) = ( 1.0_DP - xyz_CloudCover(i,j,k) )**2 / ( 2.0_DP * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) ) * xyz_DQH2OVapSatDt(i,j,k) * xyz_ZeroOneDQsDt(i,j,k)
            else
              xyz_FactE(i,j,k) = 0.0_DP
            end if
          end do
        end do
      end do
      do k = 1, kmax
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QH2OVapSat(i,j,k) > xyz_QH2OVap(i,j,k) ) then
              xyz_FactE(i,j,k) = xyz_FactE(i,j,k) + xyz_CloudCover(i,j,k)**2 * MixCoef * ( xyz_QH2OVapSat(i,j,k) - xyz_QH2OVap(i,j,k) ) / ( xyz_QCloudWater(i,j,k) + 1.0d-100 )
            end if
          end do
        end do
      end do


      ! Rain used for next k loop
      xy_Rain = 0.0_DP

      do k = kmax, 1, -1

!!$        xy_RainConvFactor = 1.0_DP / ( CloudLifeTime + 1.0d-100 )
        !
        xy_FactCo = 1.0_DP + C1 * sqrt( xy_Rain )
        xy_FactBF = 1.0_DP + C2 * sqrt( max( 268.0_DP - xyz_Temp(:,:,k), 0.0_DP ) )
        xy_QCloudWaterConvThreshold = QCloudWaterConvThreshold / ( xy_FactCo * xy_FactBF )
        xy_RainConvFactor = C0 * xy_FactCo * xy_FactBF * ( 1.0_DP - exp( - ( xyz_QCloudWater(:,:,k) / ( ( xyz_CloudCover(:,:,k) + 1.0d-100 ) * xy_QCloudWaterConvThreshold ) )**2 ) )
        !
        xyz_FactB(:,:,k) = xy_RainConvFactor
        !
        do j = 1, jmax
           do i = 0, imax-1
              if ( xyz_FactB(i,j,k) < 1.0_DP / 1.0d10 ) then
                 xyz_FactB(i,j,k) = 1.0_DP / 1.0d10
              end if
           end do
        end do


        ! Values at next time step are calculated.
        !
        xyz_QCloudWater(:,:,k) = xyz_QCloudWater(:,:,k) * exp( - xyz_FactB(:,:,k) * 2.0_DP * DelTime ) + xyz_FactA(:,:,k) / xyz_FactB(:,:,k) * ( 1.0_DP - exp( - xyz_FactB(:,:,k) * 2.0_DP * DelTime ) )
        !   The value of cloud water amount is checked, and xyz_FactA 
        !   is updated.
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < 0.0_DP ) then
              xyz_FactA2(i,j,k) = - xyz_FactA1(i,j,k) - xyz_FactB(i,j,k) * xyz_QCloudWaterB(i,j,k) * exp( - xyz_FactB(i,j,k) * 2.0_DP * DelTime ) / ( 1.0_DP - exp( - xyz_FactB(i,j,k) * 2.0_DP * DelTime ) )
              xyz_QCloudWater(i,j,k) = 0.0_DP
            end if
          end do
        end do
        xyz_FactA(:,:,k) = xyz_FactA1(:,:,k) + xyz_FactA2(:,:,k)
        !
        xyz_CloudCover(:,:,k) = xyz_CloudCover(:,:,k) * exp( - ( xyz_FactC(:,:,k) + xyz_FactD(:,:,k) ) * 2.0_DP * DelTime ) + ( xyz_FactC(:,:,k) + xyz_FactE(:,:,k) ) / ( xyz_FactC(:,:,k) + xyz_FactD(:,:,k) + 1.0d-100 ) * ( 1.0_DP - exp( - ( xyz_FactC(:,:,k) + xyz_FactD(:,:,k) ) * 2.0_DP * DelTime ) )
        !
        xyz_QH2OVap(:,:,k) = xyz_QH2OVap(:,:,k) - xyz_FactA2(:,:,k) * 2.0_DP * DelTime
        !
        xyz_Temp(:,:,k) = xyz_Temp(:,:,k) + xyz_FactA2(:,:,k) * 2.0_DP * DelTime * LatentHeat / CpDry


        ! Rain
        !
        xyz_Rain(:,:,k) = xyz_FactA(:,:,k) * 2.0_DP * DelTime + xyz_QCloudWaterB(:,:,k) * ( 1.0_DP - exp( - xyz_FactB(:,:,k) * 2.0_DP * DelTime ) ) - xyz_FactA(:,:,k) / xyz_FactB(:,:,k) * ( 1.0_DP - exp( - xyz_FactB(:,:,k) * 2.0_DP * DelTime ) )
        xyz_Rain(:,:,k) = xyz_Rain(:,:,k) / ( 2.0_DP * DelTime )

        ! Rain used for next k loop
        xy_Rain = xy_Rain + xyz_Rain(:,:,k) * ( xyr_Press(:,:,k-1) - xyr_Press(:,:,k) ) / Grav

        ! Evaporation
        !
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < QCloudWaterEvapThreshold ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do

        ! Cloud cover is restricted.
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_CloudCover(i,j,k) > 1.0_DP ) then
              xyz_CloudCover(i,j,k) = 1.0_DP
            else if ( xyz_CloudCover(i,j,k) < 0.0_DP ) then
              xyz_CloudCover(i,j,k) = 0.0_DP
            end if
          end do
        end do


        ! Check values
        do j = 1, jmax
          do i = 0, imax-1
            if ( xyz_QCloudWater(i,j,k) < 0.0_DP ) then
              write( 6, * ) i, j, k, xyz_QCloudWater(i,j,k)
            end if
          end do
        end do

      end do ! k loop


    end if


  end subroutine SetCloudCloudWaterAdjustWithCloudCoverV4
Subroutine :
xyz_CloudCover(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(in )
xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax) :real(DP), intent(inout)

[Source]

  subroutine SetCloudSmearCloudOptDep( xyz_CloudCover, xyz_DelCloudOptDep )

    ! USE statements
    !

    real(DP), intent(in   ) :: xyz_CloudCover    (0:imax-1, 1:jmax, 1:kmax)
    real(DP), intent(inout) :: xyz_DelCloudOptDep(0:imax-1, 1:jmax, 1:kmax)


    ! 実行文 ; Executable statement
    !

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


    ! Cloud optical depth is scaled by the way of Kiehl et al. (1994).

    xyz_DelCloudOptDep = xyz_DelCloudOptDep * xyz_CloudCover**1.5_DP


  end subroutine SetCloudSmearCloudOptDep
TempIceLim
Variable :
TempIceLim :real(DP), save
TempWatLim
Variable :
TempWatLim :real(DP), save
module_name
Constant :
module_name = ‘set_cloud :character(*), parameter
: モジュールの名称. Module name
set_cloud_inited
Variable :
set_cloud_inited = .false. :logical, save
: 初期設定フラグ. Initialization flag
version
Constant :
version = ’$Name: dcpam5-20120813 $’ // ’$Id: set_cloud.f90,v 1.15 2012-08-13 12:33:49 yot Exp $’ :character(*), parameter
: モジュールのバージョン Module version