!---------------------------------------------------------------------
!     Copyright (C) GFD Dennou Club, 2006. All rights reserved.
!---------------------------------------------------------------------

module physics_ps_correction_mod
  !
  != 物理過程 表面気圧調整モジュール
  !
  !== 概要
  !
  ! 凝結・蒸発による大気質量の補正を行う.
  ! 
  !
  use type_mod,    only : REKIND, DBKIND, INTKIND, TOKEN, STRING

  implicit none

  private
  public :: physics_ps_correction

contains

  subroutine physics_ps_correction( &
    & xyz_Qvap   , & !(inout)
    & xy_SurfPress, & !(inout)
    & xyr_QvapFlux , & !(in)
    & xyz_DCumulusQvapDt, & !(in)
    & xyr_Press  , &  !(in)
    & DelTimePhy          ) !(in) 
    ! 凝結・蒸発による大気質量の補正

    use type_mod,      only: REKIND, DBKIND, INTKIND, TOKEN, STRING
    use grid_3d_mod,   only: im, jm, km
    use constants_mod, only: EL    ,& ! 水の凝結の潜熱 
         &                   Grav     ! 重力加速度
    use dc_trace,      only: SetDebug, BeginSub, EndSub, DbgMessage, DataDump

    implicit none

    real(DBKIND), intent(inout) :: xyz_Qvap    ( im*jm, km )  ! 比湿 
    real(DBKIND), intent(inout) :: xy_SurfPress   ( im*jm )   !  表面気圧
    real(DBKIND), intent(in) :: xyr_Press   ( im*jm, km+1 ) ! 気圧(半整数)
    real(DBKIND), intent(in) :: xyr_QvapFlux   ( im*jm, km+1 ) ! q フラックス
    real(DBKIND), intent(in) :: xyz_DCumulusQvapDt   ( im*jm, km ) ! q 凝結量
    real(DBKIND), intent(in) :: DelTimePhy                        ! 2Δt

    real(DBKIND) :: xyz_DelPress  ( im*jm, km )   ! Δｐ
    real(DBKIND) :: xy_DelSurfPress ( im*jm )   ! Ps 変化量
    character(STRING),  parameter:: subname = "physics_ps_correction"
    integer(INTKIND) :: ij ! ループ変数
    integer(INTKIND) :: k ! ループ変数

    continue

    ! 開始処理
    call BeginSub(subname)

    ! 準備
    do k = 1, km
      do ij = 1, im*jm
            xyz_DelPress( IJ,K ) = xyr_Press( IJ,K ) - xyr_Press( IJ,K+1 )
      enddo
    enddo

    ! 大気増加量
    do ij = 1, im*jm
       xy_DelSurfPress(IJ) = xyr_QvapFlux(IJ,1)/EL*GRAV *DelTimePhy * 0.01
       do K = 1, km
         xy_DelSurfPress(IJ) = xy_DelSurfPress(IJ) &
           & + xyz_DCumulusQvapDt(IJ,K) * DelTimePhy &
           &   * xyz_DelPress(IJ,K) * 0.01
       enddo
    enddo

    ! q の値の補正
    do k = 1, km
      do ij = 1, im*jm
         xyz_Qvap ( ij,k ) &
           &  = xy_SurfPress(ij) /( xy_SurfPress(ij) + xy_DelSurfPress(ij) )&
           &    * xyz_Qvap(IJ,K)
      enddo
    enddo

    ! Ps の値の補正
    do ij = 1, im*jm
       xy_SurfPress ( IJ ) = xy_SurfPress(IJ) + xy_DelSurfPress(IJ)
    enddo 

    ! 終了処理
    call EndSub(subname)

  end subroutine physics_ps_correction

end module physics_ps_correction_mod