subroutine PhyCumulusAdjustClose( phy_cumad, err )
    !
    ! PHYCUMAD 型の変数の終了処理を行います.
    ! なお, 与えられた *phy_cumad* が Create によって初期設定
    ! されていない場合, プログラムはエラーを発生させます.
    !
    ! Deconstructor of "PHYCUMAD".
    ! Note that if *phy_cumad* is not initialized by "Create" yet,
    ! error is occurred.
    !
    use dc_trace, only: BeginSub, EndSub
    use dc_string, only: PutLine, Printf
    use dc_types, only: DP, STRING, TOKEN, STDOUT
    use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
!!$    use phy_saturate_nha92, only: PhySatNhaClose
    implicit none
    type(PHYCUMAD), intent(inout):: phy_cumad
    logical, intent(out), optional:: err
                              ! 例外処理用フラグ.
                              ! デフォルトでは, この手続き内でエラーが
                              ! 生じた場合, プログラムは強制終了します.
                              ! 引数 *err* が与えられる場合,
                              ! プログラムは強制終了せず, 代わりに
                              ! *err* に .true. が代入されます.
                              !
                              ! Exception handling flag. 
                              ! By default, when error occur in 
                              ! this procedure, the program aborts. 
                              ! If this *err* argument is given, 
                              ! .true. is substituted to *err* and 
                              ! the program does not abort. 

    !-----------------------------------
    !  作業変数
    !  Work variables
    integer:: stat
    character(STRING):: cause_c
    character(*), parameter:: subname = 'PhyCumulusAdjustClose'
  continue
    call BeginSub( subname )
    stat = DC_NOERR
    cause_c = ''

    !-----------------------------------------------------------------
    !  初期設定のチェック
    !  Check initialization
    !-----------------------------------------------------------------
    if ( .not. phy_cumad % initialized ) then
      stat = DC_ENOTINIT
      cause_c = 'PHYCUMAD'
      goto 999
    end if

    !-----------------------------------------------------------------
    !  "PHYCUMAD" の設定の消去
    !  Clear the settings for "PHYCUMAD"
    !-----------------------------------------------------------------
    deallocate( phy_cumad % TempSatMax )

!!$    !-----------------------------------------------------------------
!!$    !  飽和比湿計算用オブジェクトの終了処理
!!$    !  Terminate object for calculation of saturation specific humidity
!!$    !-----------------------------------------------------------------
!!$    call PhySatNhaClose( &
!!$      & phy_sat_nha = phy_cumad % phy_sat, & ! (inout)
!!$      & err = err )                          ! (out)

    !-----------------------------------------------------------------
    !  終了処理, 例外処理
    !  Termination and Exception handling
    !-----------------------------------------------------------------
    phy_cumad % initialized = .false.
999 continue
    call StoreError( stat, subname, err, cause_c )
    call EndSub( subname )
  end subroutine PhyCumulusAdjustClose

  subroutine PhyCumulusAdjustCreate( phy_cumad, imax, jmax, kmax, Grav, Cp, RAir, EL, RVap, EpsV, ES0, DelTime, nmlfile, err )
    !
    ! PHYCUMAD 型の変数の初期設定を行います.
    ! 他のサブルーチンを使用する前に必ずこのサブルーチンによって
    ! PHYCUMAD 型の変数を初期設定してください.
    !
    ! なお, 与えられた *phy_cumad* が既に初期設定されている場合,
    ! プログラムはエラーを発生させます.
    !
    ! NAMELIST を利用する場合には引数 *nmlfile* に NAMELIST ファイル名
    ! を与えてください. NAMELIST 変数群の詳細に関しては 
    ! NAMELIST#phy_cumulus_adjust_nml を参照してください. 
    !
    ! Constructor of "PHYCUMAD".
    ! Initialize *phy_cumad* by this subroutine, 
    ! before other procedures are used, 
    !
    ! Note that if *phy_cumad* is already initialized 
    ! by this procedure, error is occurred.
    !
    ! In order to use NAMELIST, specify a NAMELIST filename to 
    ! argument *nmlfile*. See "NAMELIST#phy_cumulus_adjust_nml"
    ! for details about a NAMELIST group.
    !
    use dc_trace, only: BeginSub, EndSub
    use dc_string, only: PutLine, Printf
    use dc_types, only: DP, STRING, TOKEN, STDOUT
    use dc_present, only: present_and_not_empty, present_and_true
    use dc_message, only: MessageNotify
    use dc_error, only: StoreError, DC_NOERR, DC_EALREADYINIT, DC_EARGLACK, DC_ENEGATIVE, DC_ENOFILEREAD
!!$    use phy_saturate_nha92, only: PhySatNhaCreate
    implicit none
    type(PHYCUMAD), intent(inout):: phy_cumad
    integer, intent(in):: imax ! 経度格子点数. 
                   ! Number of grid points in longitude
    integer, intent(in):: jmax ! 緯度格子点数. 
                   ! Number of grid points in latitude
    integer, intent(in):: kmax ! 鉛直層数. 
                   ! Number of vertical level
    real(DP), intent(in):: Grav      ! $ g $ .      重力加速度.     Gravitational acceleration
    real(DP), intent(in):: Cp        ! $ C_p $ .    大気定圧比熱.   Specific heat of air at constant pressure
    real(DP), intent(in):: RAir      ! $ R $ .      大気気体定数.   Gas constant of air
    real(DP), intent(in):: EL        ! $ L $ .      水の凝結の潜熱. Latent heat of condensation of water vapor
    real(DP), intent(in):: RVap      ! $ R_v $ .    水蒸気気体定数. Gas constant of water vapor
    real(DP), intent(in):: EpsV      ! $ \epsilon_v $ .        水蒸気分子量比. Molecular weight ratio of water vapor
    real(DP), intent(in):: ES0       ! $ e^{*} $ (273K) .      0 ℃での飽和蒸気圧. Saturation vapor pressure at 0 degrees C
    real(DP), intent(in):: DelTime    ! $ \Delta t $ . タイムステップ. Time step
    character(*), intent(in), optional:: nmlfile
                              ! NAMELIST ファイルの名称. 
                              ! この引数に空文字以外を与えた場合, 
                              ! 指定されたファイルから 
                              ! NAMELIST 変数群を読み込みます. 
                              ! ファイルを読み込めない場合にはエラーを
                              ! 生じます.
                              !
                              ! NAMELIST 変数群の詳細に関しては 
                              ! NAMELIST#phy_cumulus_adjust_nml 
                              ! を参照してください. 
                              !
                              ! NAMELIST file name. 
                              ! If nonnull character is specified to
                              ! this argument, 
                              ! NAMELIST group name is loaded from the 
                              ! file. 
                              ! If the file can not be read, 
                              ! an error occurs.
                              ! 
                              ! See "NAMELIST#phy_cumulus_adjust_nml" 
                              ! for details about a NAMELIST group.
                              ! 
    logical, intent(out), optional:: err
                              ! 例外処理用フラグ.
                              ! デフォルトでは, この手続き内でエラーが
                              ! 生じた場合, プログラムは強制終了します.
                              ! 引数 *err* が与えられる場合,
                              ! プログラムは強制終了せず, 代わりに
                              ! *err* に .true. が代入されます.
                              !
                              ! Exception handling flag. 
                              ! By default, when error occur in 
                              ! this procedure, the program aborts. 
                              ! If this *err* argument is given, 
                              ! .true. is substituted to *err* and 
                              ! the program does not abort. 

    !-----------------------------------
    !  作業変数
    !  Work variables
    integer:: stat
    character(STRING):: cause_c
    character(*), parameter:: subname = 'PhyCumulusAdjustCreate'
  continue
    call BeginSub( subname, version )
    stat = DC_NOERR
    cause_c = ''

    !-----------------------------------------------------------------
    !  初期設定のチェック
    !  Check initialization
    !-----------------------------------------------------------------
    if ( phy_cumad % initialized ) then
      stat = DC_EALREADYINIT
      cause_c = 'PHYCUMAD'
      goto 999
    end if

    !-----------------------------------------------------------------
    !  引数の正当性のチェック
    !  Validate arguments
    !-----------------------------------------------------------------
    if (imax < 1) then
      stat = DC_ENEGATIVE
      cause_c = 'imax'
      goto 999
    end if
    if (jmax < 1) then
      stat = DC_ENEGATIVE
      cause_c = 'jmax'
      goto 999
    end if
    if (kmax < 1) then
      stat = DC_ENEGATIVE
      cause_c = 'kmax'
      goto 999
    end if
    if (DelTime < 0.0_DP) then
      stat = DC_ENEGATIVE
      cause_c = 'DelTime'
      goto 999
    end if

    !-----------------------------------------------------------------
    !  "PHYCUMAD" の設定
    !  Configure the settings for "PHYCUMAD"
    !-----------------------------------------------------------------
    phy_cumad % imax    = imax   
    phy_cumad % jmax    = jmax   
    phy_cumad % kmax    = kmax   
    phy_cumad % Grav    = Grav   
    phy_cumad % Cp      = Cp     
    phy_cumad % RAir    = RAir   
    phy_cumad % EL      = EL     
    phy_cumad % RVap    = RVap   
    phy_cumad % EpsV    = EpsV   
    phy_cumad % ES0     = ES0    
    phy_cumad % DelTime = DelTime

    !-------------------------
    !  デフォルト値
    !  Default values
!!$    phy_cumad % param_r = 0.0_DP
!!$    phy_cumad % param_c = 'hogehoge'

    !-------------------------
    !  オプショナル引数からの値
    !  Values from optional arguments
!!$    phy_cumad % param_i = param_i
!!$    if ( present(param_r) )  phy_cumad % param_r = param_r
!!$    if ( present(param_c) )  phy_cumad % param_c = param_c

    !-------------------------
    !  NAMELIST からの値
    !  Values from NAMELIST

!!$    if ( present_and_not_empty(nmlfile) ) then
!!$      call MessageNotify( 'M', subname, &
!!$        & 'Loading NAMELIST file "%c" ...', &
!!$        & c1 = trim(nmlfile) )
!!$      call NmlRead ( nmlfile = nmlfile, &      ! (in)
!!$        & param_i = phy_cumad % param_i, &   ! (inout)
!!$        & param_r = phy_cumad % param_r, &   ! (inout)
!!$        & param_c_ = phy_cumad % param_c, &  ! (inout)
!!$        & err = err )                          ! (out)
!!$      if ( present_and_true(err) ) then
!!$        call MessageNotify( 'W', subname, &
!!$          & '"%c" can not be read.', &
!!$          & c1 = trim(nmlfile) )
!!$        stat = DC_ENOFILEREAD
!!$        cause_c = nmlfile
!!$        goto 999
!!$      end if
!!$    end if

    !-----------------------------------------------------------------
    !  臨界相対湿度や不安定の許容誤差の設定
    !  Configure critical relative humidity and admissible error of unstability
    !-----------------------------------------------------------------
    phy_cumad % CrtlRH = 0.990_DP
    phy_cumad % IterationMax = 10

    allocate( phy_cumad % TempSatMax (1:phy_cumad % IterationMax) )

    phy_cumad % TempSatMax(1:phy_cumad % IterationMax) = (/ 0.01_DP, 0.02_DP, 0.02_DP, 0.05_DP, 0.05_DP, 0.10_DP, 0.10_DP, 0.20_DP, 0.20_DP, 0.40_DP  /)

!!$    !-----------------------------------------------------------------
!!$    !  飽和比湿計算用オブジェクトの初期化
!!$    !  Initialize object for calculation of saturation specific humidity
!!$    !-----------------------------------------------------------------
!!$    call PhySatNhaCreate( &
!!$      & phy_sat_nha = phy_cumad % phy_sat, &     ! (out)
!!$      & imax = imax, jmax = jmax, kmax = kmax, & ! (in)
!!$      & EpsV = EpsV, &                           ! (in)
!!$      & err = err )                              ! (out)

    !-----------------------------------------------------------------
    !  設定値の正当性のチェック
    !  Validate setting values
    !-----------------------------------------------------------------
!!$    if ( phy_cumad % param_i < 0 ) then
!!$      stat = DC_ENEGATIVE
!!$      cause_c = 'param_i'
!!$      goto 999
!!$    end if


    !-----------------------------------------------------------------
    !  終了処理, 例外処理
    !  Termination and Exception handling
    !-----------------------------------------------------------------
    phy_cumad % initialized = .true.
999 continue
    call StoreError( stat, subname, err, cause_c )
    call EndSub( subname )
  end subroutine PhyCumulusAdjustCreate

  subroutine PhyCumulusAdjustCumulus( phy_cumad, xyz_Temp, xyz_QVap, xyz_Press, xyr_Press, xy_Rain, xyz_DTempDt, xyz_DQVapDt, err )
    !
    ! 対流調節スキームにより, 温度と比湿を調節します. 
    !
    ! なお, 与えられた *phy_cumad* が Create によって初期設定
    ! されていない場合, プログラムはエラーを発生させます.
    !
    ! Adjust temperature and specific humidity by 
    ! convective adjustment scheme.
    !
    ! If *phy_cumad* is not initialized by "Create" yet,
    ! error is occurred.
    !
    use dc_trace, only: BeginSub, EndSub
    use dc_string, only: PutLine, Printf
    use dc_types, only: DP, STRING, TOKEN, STDOUT
    use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
!!$    use phy_saturate_nha92, only: CalcQVapSat, CalcDQVapSatDTemp
    implicit none
    type(PHYCUMAD), intent(inout):: phy_cumad
    real(DP), intent(inout):: xyz_Temp (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! $ T $ .     温度. Temperature
    real(DP), intent(inout):: xyz_QVap (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! $ q $ .     比湿. Specific humidity
    real(DP), intent(in):: xyz_Press (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! $ p $ . 気圧 (整数レベル). 
                              ! Surface pressure (full level)
    real(DP), intent(in):: xyr_Press (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax)
                              ! $ p $ . 気圧 (半整数レベル). 
                              ! Surface pressure (half level)
    real(DP), intent(out):: xy_Rain (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1)
                              ! 降水量. 
                              ! Precipitation
    real(DP), intent(out):: xyz_DTempDt (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! 温度変化率. 
                              ! Temperature tendency
    real(DP), intent(out):: xyz_DQVapDt (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! 比湿変化率. 
                              ! Specific humidity tendency
    logical, intent(out), optional:: err
                              ! 例外処理用フラグ.
                              ! デフォルトでは, この手続き内でエラーが
                              ! 生じた場合, プログラムは強制終了します.
                              ! 引数 *err* が与えられる場合,
                              ! プログラムは強制終了せず, 代わりに
                              ! *err* に .true. が代入されます.
                              !
                              ! Exception handling flag. 
                              ! By default, when error occur in 
                              ! this procedure, the program aborts. 
                              ! If this *err* argument is given, 
                              ! .true. is substituted to *err* and 
                              ! the program does not abort. 

    !-----------------------------------
    !  作業変数
    !  Work variables
    integer:: imax ! 経度格子点数. 
                   ! Number of grid points in longitude
    integer:: jmax ! 緯度格子点数. 
                   ! Number of grid points in latitude
    integer:: kmax ! 鉛直層数. 
                   ! Number of vertical level
    real(DP):: Grav      ! $ g $ .      重力加速度.     Gravitational acceleration
    real(DP):: Cp        ! $ C_p $ .    大気定圧比熱.   Specific heat of air at constant pressure
    real(DP):: DelTime    ! $ \Delta t $ . タイムステップ. Time step

    real(DP):: CrtlRH
                              ! 臨界相対湿度. 
                              ! Critical relative humidity
    integer:: IterationMax
                              ! イテレーション回数. 
                              ! Number of iteration

    real(DP), allocatable:: TempSatMax(:)
                              ! 不安定の許容誤差. 
                              ! Admissible error of unstability

    real(DP):: RAir      ! $ R $ .      大気気体定数.   Gas constant of air
    real(DP):: EL        ! $ L $ .      水の凝結の潜熱. Latent heat of condensation of water vapor
    real(DP):: RVap      ! $ R_v $ .    水蒸気気体定数. Gas constant of water vapor
    real(DP):: EpsV      ! $ \epsilon_v $ .        水蒸気分子量比. Molecular weight ratio of water vapor
    real(DP):: ES0       ! $ e^{*} $ (273K) .      0 ℃での飽和蒸気圧. Saturation vapor pressure at 0 degrees C

    real(DP):: RUniv
                              ! $ R $ [J mol-1 K-1].
                              ! 普遍気体定数. Universal gas constant
    real(DP):: LatHeat
                              ! $ l $ [J mol-1].
                              ! 水の凝結の潜熱. Latent heat of condensation of water vapor
    real(DP):: P0
                              ! $ p_0^{*} $ [Pa].
                              ! 水蒸気飽和曲線の定数. constant for water vapor saturation curve



    real(DP):: xyz_QVapB (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! 調節前の比湿. 
                              ! Specific humidity before adjust. 
    real(DP):: xyz_TempB (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! 調節前の温度. 
                              ! Temperature before adjust. 
    logical:: xy_Adjust (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1)
                              ! 今回調節されたか否か?. 
                              ! Whether it was adjusted this time or not?
    logical:: xy_AdjustB (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1)
                              ! 前回調節されたか否か?. 
                              ! Whether it was adjusted last time or not?
    real(DP):: xyz_DPressDz (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! $ \DD{p}{z} $
                              !
    real(DP):: xyz_QVapSat (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! 飽和比湿. 
                              ! Saturation specific humidity. 
    real(DP):: xyz_DDPressDDPress (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! $ \DD{p_{k}}{p_{k-1}} $
                              !
    real(DP):: xyz_DPFact (0:phy_cumad%imax-1, 0:phy_cumad%jmax-1, 0:phy_cumad%kmax-1)
                              ! $ (R / C_p)
                              ! \frac{p_{k-1} - p_{k}}{2 p_{k-1/2}} $ . 
                              !
                              ! ファクター. 
                              ! Factor
    real(DP):: TempSat        ! $ S_t $ .
                              ! 飽和温度. 
                              ! Saturation temperature
    real(DP):: DelTempSat
                              ! 調節による飽和温度の変化量. 
                              ! Saturation temperature variation by adjustment
    real(DP):: DelQVap
                              ! 調節による比湿の変化量. 
                              ! Specific humidity variation by adjustment
    real(DP):: DelTempUpper
                              ! 調節による温度 (k) の変化量. 
                              ! Temperature (k) variation by adjustment
    real(DP):: DelTempLower
                              ! 調節による温度 (k-1) の変化量. 
                              ! Temperature (k-1) variation by adjustment
    real(DP):: DQVapSatDTempUpper
                              ! $ \DD{q^{*}} (k)}{T} $
    real(DP):: DQVapSatDTempLower
                              ! $ \DD{q^{*}} (k-1)}{T} $
    real(DP):: DHDTempUpper
                              ! $ 1 + \gamma_{k} =
                              ! 1 + \frac{L}{C_p} \DP{q^{*}}{T}_{k} $
    real(DP):: DHDTempLower
                              ! $ 1 + \gamma_{k-1} =
                              ! 1 + \frac{L}{C_p} \DP{q^{*}}{T}_{k-1} $
    logical:: Adjust
                              ! 今回全領域において一度でも調節されたか否か?. 
                              ! Whether it was adjusted even once in global 
                              ! this time or not?

    integer:: Iteration
    integer:: i, j, k         ! DO ループ用作業変数
                              ! Work variables for DO loop

    integer:: stat
    character(STRING):: cause_c
    character(*), parameter:: subname = 'PhyCumulusAdjustCumulus'
  continue
    call BeginSub( subname )
    stat = DC_NOERR
    cause_c = ''

    !-----------------------------------------------------------------
    !  初期設定のチェック
    !  Check initialization
    !-----------------------------------------------------------------
    if ( .not. phy_cumad % initialized ) then
      stat = DC_ENOTINIT
      cause_c = 'PHYCUMAD'
      goto 999
    end if

    !-----------------------------------------------------------------
    !  *phy_cumad* に格納されている設定値の取り出し
    !  Fetch setting values stored in *phy_cumad*
    !-----------------------------------------------------------------
    imax    = phy_cumad % imax   
    jmax    = phy_cumad % jmax   
    kmax    = phy_cumad % kmax   
    Grav    = phy_cumad % Grav   
    Cp      = phy_cumad % Cp     
    DelTime = phy_cumad % DelTime

    CrtlRH       = phy_cumad % CrtlRH
    IterationMax = phy_cumad % IterationMax

    RAir    = phy_cumad % RAir   
    EL      = phy_cumad % EL     
    RVap    = phy_cumad % RVap   
    EpsV    = phy_cumad % EpsV   
    ES0     = phy_cumad % ES0    
    RUniv   = phy_cumad % RUniv  
    LatHeat = phy_cumad % LatHeat
    P0      = phy_cumad % P0     

    if ( allocated(TempSatMax) ) deallocate( TempSatMax )
    allocate( TempSatMax(1:IterationMax) )
    TempSatMax = phy_cumad % TempSatMax

    !-----------------------------------------------------------------
    !  調節前 "QVap", "Temp" の保存
    !  Store "QVap", "Temp" before adjustment
    !-----------------------------------------------------------------
    xyz_QVapB  = xyz_QVap
    xyz_TempB  = xyz_Temp
    
    !-----------------------------------------------------------------
    !  ファクターの計算
    !  Calculate factor
    !-----------------------------------------------------------------
    do k = 0, kmax - 1
      xyz_DPressDz(:,:,k) = xyr_Press(:,:,k) - xyr_Press(:,:,k+1)
    end do

    !-----------------------------------------------------------------
    !  飽和比湿計算
    !  Calculate saturation specific humidity
    !-----------------------------------------------------------------

    !--------------------------------
    !  Nakajima et al. (1992) を用いた飽和比湿の計算
    !  Calculate saturation specific humidity with Nakajima et al. (1992) 
    xyz_QVapSat = EpsV * ( P0 / xyz_Press ) * exp ( - LatHeat / ( RUniv * xyz_Temp ) )

!!$    call CalcQVapSat( phy_cumad % phy_sat, & ! (in)
!!$      &    xyz_Temp = xyz_Temp, &            ! (in)
!!$      &   xyz_Press = xyz_Press, &           ! (in)
!!$      & xyz_QVapSat = xyz_QVapSat, &         ! (out)
!!$      &         err = err )                  ! (out)

!!$    !--------------------------------
!!$    !  クラウジウスクラペイロンの式より飽和比湿計算
!!$    !  * 273 K の時の飽和水蒸気圧を ES0 (= 611 Pa) としている
!!$    !    * 浅井 冨男, 武田 喬男, 木村 竜治, 1981: 
!!$    !      大気科学講座 2 雲や降水を伴う大気, 東京大学出版会, 249pp. 
!!$    !  Calculate saturation specific humidity with Clausius Clapeyron equation
!!$    !  * Saturation vapor pressure is ES0 (= 611 Pa) at 273 K
!!$    xyz_QVapSat = EpsV * ES0  &
!!$      & * exp( EL / RVap * ( 1.0_DP / 273.0_DP - 1.0_DP / xyz_Temp ) ) &
!!$      &    / xyz_Press

    do k = 1, kmax - 1
      xyz_DDPressDDPress(:,:,k) = xyz_DPressDz(:,:,k) / xyz_DPressDz(:,:,k-1)
      xyz_DPFact(:,:,k) = RAir / Cp * ( xyz_Press(:,:,k-1) - xyz_Press(:,:,k) ) / xyr_Press(:,:,k) / 2.0_DP
    end do

    !-----------------------------------------------------------------
    !  調節
    !  Adjustment
    !-----------------------------------------------------------------
    xy_AdjustB = .true.
    
    !-------------------------
    !  イテレーション
    !  iteration
    do Iteration = 1, IterationMax
      xy_Adjust = .false.
      
      do k = 1, kmax-1
        do j = 0, jmax-1
          do i = 0, imax-1
            if ( xy_AdjustB(i,j) ) then
              
              TempSat = xyz_Temp(i,j,k-1) - xyz_Temp(i,j,k) + ( xyz_QVapSat(i,j,k-1) - xyz_QVapSat(i,j,k) ) * EL / Cp - xyz_DPFact(i,j,k) * ( xyz_Temp(i,j,k-1) + xyz_Temp(i,j,k) )
              
              !-------------------------
              !  不安定であるならば
              !  If it is unstable
              if ( TempSat > TempSatMax(Iteration) ) then

                !-------------------------
                !  .. かつ, 飽和しているならば
                !  .. and, if it is saturated
                if ( ( xyz_QVap(i,j,k) / xyz_QVapSat(i,j,k) >= CrtlRH ) .and. ( xyz_QVap(i,j,k-1) / xyz_QVapSat(i,j,k-1) >= CrtlRH ) ) then
                  
                  DelQVap = xyz_DPressDz(i,j,k-1) * (xyz_QVap(i,j,k-1) - xyz_QVapSat(i,j,k-1) ) + xyz_DPressDz(i,j,k) * (xyz_QVap(i,j,k) - xyz_QVapSat(i,j,k) ) 

                  DQVapSatDTempUpper = xyz_QVapSat(i,j,k) * ( LatHeat / ( RUniv * xyz_Temp(i,j,k)**2 ) )
                  
                  DQVapSatDTempLower = xyz_QVapSat(i,j,k-1) * ( LatHeat / ( RUniv * xyz_Temp(i,j,k-1)**2 ) )

!!$                  call CalcDQVapSatDTemp( phy_cumad % phy_sat, &  ! (in)
!!$                    &          Temp = xyz_Temp(i,j,k), &          ! (in)
!!$                    &         Press = xyz_Press(i,j,k), &         ! (in)
!!$                    & DQVapSatDTemp = DQVapSatDTempUpper )        ! (out)
!!$
!!$                  call CalcDQVapSatDTemp( phy_cumad % phy_sat, &  ! (in)
!!$                    &          Temp = xyz_Temp(i,j,k-1), &        ! (in)
!!$                    &         Press = xyz_Press(i,j,k-1), &       ! (in)
!!$                    & DQVapSatDTemp = DQVapSatDTempLower )        ! (out)
!!$
!!$                  DQvapSatDTempUpper = EL * xyz_QvapSat(i,j,k) &
!!$                    & / ( RVap * xyz_Temp(i,j,k) * xyz_Temp(i,j,k) )
!!$                  
!!$                  DQvapSatDTempLower = EL * xyz_QvapSat(i,j,k-1) &
!!$                    & / ( RVap * xyz_Temp(i,j,k-1) * xyz_Temp(i,j,k-1) )
                  
                  DHDTempUpper = 1.0_DP + EL/Cp * DQVapSatDTempUpper
                  DHDTempLower = 1.0_DP + EL/Cp * DQVapSatDTempLower
                  
                  DelTempSat = TempSat + ( 1.0_DP - xyz_DPFact(i,j,k) / DHDTempLower ) * EL/Cp * DelQVap / xyz_DPressDz(i,j,k-1)
                  
                  !-------------------------
                  !  温度の調節
                  !  Adjust temperature
                  DelTempUpper = DelTempSat / ( ( 1.0_DP + xyz_DDPressDDPress(i,j,k) ) * DHDTempUpper + xyz_DPFact(i,j,k) * ( 1.0_DP - xyz_DDPressDDPress(i,j,k) * DHDTempUpper / DHDTempLower ) )
                  
                  DelTempLower = - DHDTempUpper / DHDTempLower * xyz_DDPressDDPress(i,j,k) * DelTempUpper + EL / Cp * DelQVap / ( xyz_DPressDz(i,j,k-1) * DHDTempLower )
                  
                  xyz_Temp(i,j,k)   = xyz_Temp(i,j,k) + DelTempUpper
                  xyz_Temp(i,j,k-1) = xyz_Temp(i,j,k-1) + DelTempLower
                  
                  !-------------------------
                  !  比湿の調節
                  !  Adjust specific humidity
                  xyz_QVap(i,j,k) = xyz_QVapSat(i,j,k) + DQVapSatDTempUpper * DelTempUpper

                  xyz_QVap(i,j,k-1) = xyz_QVapSat(i,j,k-1) + DQVapSatDTempLower * DelTempLower

                  xyz_QVapSat(i,j,k)   = xyz_QVap(i,j,k)
                  xyz_QVapSat(i,j,k-1) = xyz_QVap(i,j,k-1)
                  
                  !-------------------------
                  !  調節したか否か?
                  !  Whether it was adjusted or not?
                  xy_Adjust(i,j) = .true.
                end if
                
              end if
              
            end if
          end do
        end do
      end do
      
      Adjust = .false.
      do i = 0, imax-1
        do j = 0, jmax-1
          xy_AdjustB(i,j) = xy_Adjust(i,j)
          Adjust          = Adjust .or. xy_Adjust(i,j)
        end do
      end do
      
      if ( .not. Adjust ) exit
      
    end do
    
    !-----------------------------------------------------------------
    !  比湿変化率, 温度変化率, 降水量の算出
    !  Calculate specific humidity tendency, temperature tendency, 
    !  precipitation
    !-----------------------------------------------------------------
    xy_Rain     = 0.0_DP
    xyz_DTempDt = 0.0_DP
    xyz_DQVapDt = 0.0_DP

    xyz_DQVapDt = xyz_DQVapDt + ( xyz_QVap - xyz_QVapB ) / ( 2.0_DP * DelTime )

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

    do k = 0, kmax-1
      xy_Rain = xy_Rain + ( xyz_QVapB(:,:,k) - xyz_QVap(:,:,k) ) * xyz_DPressDz(:,:,k) / Grav * EL / ( 2.0_DP * DelTime )
    end do

    !-----------------------------------------------------------------
    !  終了処理, 例外処理
    !  Termination and Exception handling
    !-----------------------------------------------------------------
999 continue
    if ( allocated(TempSatMax) ) deallocate( TempSatMax )
    call StoreError( stat, subname, err, cause_c )
    call EndSub( subname )
  end subroutine PhyCumulusAdjustCumulus

  logical function PhyCumulusAdjustInitialized( phy_cumad ) result(result)
    !
    ! *phy_cumad* が初期設定されている場合には .true. が,
    ! 初期設定されていない場合には .false. が返ります.
    !
    ! If *phy_cumad* is initialized, .true. is returned.
    ! If *phy_cumad* is not initialized, .false. is returned.
    !
    implicit none
    type(PHYCUMAD), intent(in):: phy_cumad
  continue
    result = phy_cumad % initialized
  end function PhyCumulusAdjustInitialized

  subroutine PhyCumulusAdjustNmlRead( nmlfile, err )
    !
    ! NAMELIST ファイル *nmlfile* から値を入力するための
    ! 内部サブルーチンです. Create 内で呼び出されることを
    ! 想定しています.
    !
    ! 値が NAMELIST ファイル内で指定されていない場合には,
    ! 入力された値がそのまま返ります.
    !
    ! なお, *nmlfile* に空文字が与えられた場合, または
    ! 与えられた *nmlfile* を読み込むことができない場合, 
    ! プログラムはエラーを発生させます.
    !
    ! This is an internal subroutine to input values from 
    ! NAMELIST file *nmlfile*. This subroutine is expected to be
    ! called by "Create".
    !
    ! A value not specified in NAMELIST file is returned
    ! without change.
    !
    ! If *nmlfile* is empty, or *nmlfile* can not be read, 
    ! error is occurred.
    !
    use dc_trace, only: BeginSub, EndSub
    use dc_string, only: PutLine, Printf
    use dc_types, only: DP, STRING, TOKEN, STDOUT
    use dc_iounit, only: FileOpen
    use dc_message, only: MessageNotify
    use dc_present, only: present_and_true
    use dc_error, only: StoreError, DC_NOERR, DC_ENOFILEREAD
    implicit none
    character(*), intent(in):: nmlfile
                              ! NAMELIST ファイルの名称. 
                              ! NAMELIST file name
!!$    integer, intent(inout):: param_i
!!$    real(DP), intent(inout):: param_r
!!$    character(*), intent(inout):: param_c_
!!$    character(TOKEN):: param_c
    logical, intent(out), optional:: err
                              ! 例外処理用フラグ.
                              ! デフォルトでは, この手続き内でエラーが
                              ! 生じた場合, プログラムは強制終了します.
                              ! 引数 *err* が与えられる場合,
                              ! プログラムは強制終了せず, 代わりに
                              ! *err* に .true. が代入されます.
                              !
                              ! Exception handling flag. 
                              ! By default, when error occur in 
                              ! this procedure, the program aborts. 
                              ! If this *err* argument is given, 
                              ! .true. is substituted to *err* and 
                              ! the program does not abort. 

!!$    namelist /phy_cumulus_adjust_nml/ &
!!$      & param_i, param_r, param_c
                              ! phy_cumulus_adjust モジュール用
                              ! NAMELIST 変数群名.
                              !
                              ! phy_cumulus_adjust#Create を使用する際に, 
                              ! オプショナル引数 *nmlfile* へ NAMELIST 
                              ! ファイル名を指定することで, そのファイルから
                              ! この NAMELIST 変数群を読み込みます.
                              !
                              ! NAMELIST group name for 
                              ! "phy_cumulus_adjust" module.
                              ! 
                              ! If a NAMELIST filename is specified to 
                              ! an optional argument *nmlfile* 
                              ! when "phy_cumulus_adjust#Create" is used, 
                              ! this NAMELIST group is loaded from 
                              ! the file.

    !-----------------------------------
    !  作業変数
    !  Work variables
    integer:: stat
    character(STRING):: cause_c
    integer:: unit_nml        ! NAMELIST ファイルオープン用装置番号. 
                              ! Unit number for NAMELIST file open
!!$    integer:: iostat_nml      ! NAMELIST 読み込み時の IOSTAT. 
!!$                              ! IOSTAT of NAMELIST read
    character(*), parameter:: subname = 'PhyCumulusAdjustNmlRead'
  continue
    call BeginSub( subname )
    stat = DC_NOERR
    cause_c = ''



    !-----------------------------------------------------------------
    !  文字型引数を NAMELIST 変数群へ代入
    !  Substitute character arguments to NAMELIST group
    !-----------------------------------------------------------------
!!$    param_c = param_c_

    !----------------------------------------------------------------
    !  NAMELIST ファイルのオープン
    !  Open NAMELIST file
    !----------------------------------------------------------------
    call FileOpen( unit = unit_nml, file = nmlfile, mode = 'r', err = err )                   ! (out)
    if ( present_and_true(err) ) then
      stat = DC_ENOFILEREAD
      cause_c = nmlfile
      goto 999
    end if


    !-----------------------------------------------------------------
    !  NAMELIST 変数群の取得
    !  Get NAMELIST group
    !-----------------------------------------------------------------
!!$    read( unit = unit_nml, & ! (in)
!!$      & nml = phy_cumulus_adjust_nml, iostat = iostat_nml ) ! (out)
!!$    if ( iostat_nml == 0 ) then
!!$      call MessageNotify( 'M', subname, &
!!$        & 'NAMELIST group "%c" is loaded from "%c".', &
!!$        & c1 = 'phy_cumulus_adjust_nml', c2 = trim(nmlfile) )
!!$      write(STDOUT, nml = phy_cumulus_adjust_nml)
!!$    else
!!$      call MessageNotify( 'W', subname, &
!!$        & 'NAMELIST group "%c" is not found in "%c" (iostat=%d).', &
!!$        & c1 = 'phy_cumulus_adjust_nml', c2 = trim(nmlfile), &
!!$        & i = (/iostat_nml/) )
!!$    end if

    close( unit_nml )

    !-----------------------------------------------------------------
    !  NAMELIST 変数群を文字型引数へ代入
    !  Substitute NAMELIST group to character arguments
    !-----------------------------------------------------------------
!!$    param_c_ = param_c

    !-----------------------------------------------------------------
    !  終了処理, 例外処理
    !  Termination and Exception handling
    !-----------------------------------------------------------------
999 continue
    call StoreError( stat, subname, err, cause_c )
    call EndSub( subname )
  end subroutine PhyCumulusAdjustNmlRead

  subroutine PhyCumulusAdjustPutLine( phy_cumad, unit, indent, err )
    !
    ! 引数 *phy_cumad* に設定されている情報を印字します.
    ! デフォルトではメッセージは標準出力に出力されます. 
    ! *unit* に装置番号を指定することで, 出力先を変更することが可能です.
    !
    ! Print information of *phy_cumad*.
    ! By default messages are output to standard output.
    ! Unit number for output can be changed by *unit* argument.
    !
    use dc_trace, only: BeginSub, EndSub
    use dc_string, only: PutLine, Printf
    use dc_types, only: DP, STRING, TOKEN, STDOUT
    use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
!!$    use phy_saturate_nha92, only: PhySatNhaPutLine
    implicit none
    type(PHYCUMAD), intent(in):: phy_cumad
    integer, intent(in), optional:: unit
                              ! 出力先の装置番号.
                              ! デフォルトの出力先は標準出力.
                              !
                              ! Unit number for output.
                              ! Default value is standard output.
    character(*), intent(in), optional:: indent
                              ! 表示されるメッセージの字下げ.
                              !
                              ! Indent of displayed messages.
    logical, intent(out), optional:: err
                              ! 例外処理用フラグ.
                              ! デフォルトでは, この手続き内でエラーが
                              ! 生じた場合, プログラムは強制終了します.
                              ! 引数 *err* が与えられる場合,
                              ! プログラムは強制終了せず, 代わりに
                              ! *err* に .true. が代入されます.
                              !
                              ! Exception handling flag. 
                              ! By default, when error occur in 
                              ! this procedure, the program aborts. 
                              ! If this *err* argument is given, 
                              ! .true. is substituted to *err* and 
                              ! the program does not abort. 

    !-----------------------------------
    !  作業変数
    !  Work variables
    integer:: stat
    character(STRING):: cause_c
    integer:: out_unit
    integer:: indent_len
    character(STRING):: indent_str
    character(*), parameter:: subname = 'PhyCumulusAdjustPutLine'
  continue
    call BeginSub( subname )
    stat = DC_NOERR
    cause_c = ''

    !-----------------------------------------------------------------
    !  初期設定のチェック
    !  Check initialization
    !-----------------------------------------------------------------
    if ( present(unit) ) then
      out_unit = unit
    else
      out_unit = STDOUT
    end if

    indent_len = 0
    indent_str = ''
    if ( present(indent) ) then
      if ( len(indent) /= 0 ) then
        indent_len = len(indent)
        indent_str(1:indent_len) = indent
      end if
    end if


    !-----------------------------------------------------------------
    !  "PHYCUMAD" の設定の印字
    !  Print the settings for "PHYCUMAD"
    !-----------------------------------------------------------------
    if ( phy_cumad % initialized ) then
      call Printf( out_unit, indent_str(1:indent_len) // '#<PHYCUMAD:: @initialized=%y', l = (/phy_cumad % initialized/) )

      call Printf( out_unit, indent_str(1:indent_len) // ' @imax=%d @jmax=%d @kmax=%d', i = (/phy_cumad % imax, phy_cumad % jmax, phy_cumad % kmax/) )

      call Printf( out_unit, indent_str(1:indent_len) // ' @Grav=%f @Cp=%f @RAir=%f', d = (/ phy_cumad % Grav, phy_cumad % Cp, phy_cumad % RAir/) )

      call Printf( out_unit, indent_str(1:indent_len) // ' @EL=%f @RVap=%f @EpsV=%f @ES0=%f', d = (/ phy_cumad % EL, phy_cumad % RVap, phy_cumad % EpsV, phy_cumad % ES0/) )

      call Printf( out_unit, indent_str(1:indent_len) // ' @DelTime=%f', d = (/phy_cumad % DelTime/) )

      call Printf( out_unit, indent_str(1:indent_len) // ' @CrtlRH=%f @IterationMax=%d', d = (/ phy_cumad % CrtlRH/), i = (/ phy_cumad % IterationMax/) )

      call PutLine( phy_cumad % TempSatMax, unit = out_unit, lbounds = lbound(phy_cumad % TempSatMax), ubounds = ubound(phy_cumad % TempSatMax), indent = indent_str(1:indent_len) // ' @TempSatMax=' )

      call Printf( unit = out_unit, fmt = indent_str(1:indent_len) // ' @RUniv=%f @LatHeat=%f @P0=%f', d = (/phy_cumad % RUniv, phy_cumad % LatHeat, phy_cumad % P0/) )               ! (in)

!!$      call Printf( out_unit, &
!!$        & indent_str(1:indent_len) // ' @phy_sat=' )
!!$      call PhySatNhaPutLine( &
!!$        & phy_sat_nha = phy_cumad % phy_sat, &         ! (in)
!!$        & unit = out_unit, &                           ! (in)
!!$        & indent = indent_str(1:indent_len) // '  ', & ! (in)
!!$        & err = err )                                  ! (out)

      call Printf( out_unit, indent_str(1:indent_len) // '>' )
    else
      call Printf( out_unit, indent_str(1:indent_len) // '#<PHYCUMAD:: @initialized=%y>', l = (/phy_cumad % initialized/) )
    end if

    !-----------------------------------------------------------------
    !  終了処理, 例外処理
    !  Termination and Exception handling
    !-----------------------------------------------------------------
999 continue
    call StoreError( stat, subname, err, cause_c )
    call EndSub( subname )
  end subroutine PhyCumulusAdjustPutLine