!= 定常状態における \sigma 座標での質量流線関数計算
!
!= Calculation of mass stream function in the \sigma coordinate system
!  under steady-state conditions
!
! Authors::   Shin-ich Takehiro
! Version::   $Id: streamfunc.F90,v 1.0 2025/12/29 00:00:00 takepiro Exp $ 
! Tag Name::  $Name:  $
! Copyright:: Copyright (C) GFD Dennou Club, 2008. All rights reserved.
! License::   See COPYRIGHT[link:../../../COPYRIGHT]
!
module mass_streamfunc
  
  ! 格子点設定
  ! Grid points settings
  !
  use gridset, only: imax, & ! 経度格子点数. 
                             ! Number of grid points in longitude
    &                jmax, & ! 緯度格子点数. 
                             ! Number of grid points in latitude
    &                kmax    ! 鉛直層数. 
                             ! Number of vertical level

  ! 種別型パラメタ
  ! Kind type parameter
  !
  use dc_types, only: DP, &      ! 倍精度実数型. Double precision. 
    &                 STRING     ! 文字列.       Strings. 


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

  ! 宣言文 ; Declaration statements
  !
  implicit none
  private

  ! 公開手続き
  ! Public procedure
  !
  public yr_MassStreamFunc
  
contains
  function yr_MassStreamFunc( xy_Ps, xyz_V ) 

    ! 物理定数設定
    ! Physical constants settings
    !
    use constants, only: &
      & RPlanet, &
                              ! $ a $ [m]. 
                              ! 惑星半径. 
                              ! Radius of planet
      & Grav
                              ! $ g $ [m s-2]. 
                              ! 重力加速度. 
                              ! Gravitational acceleration

    ! 座標データ設定
    ! Axes data settings
    !
    use axesset, only: &
      & y_Lat,  &             ! $ \varphi $ [rad.] .
                              ! 緯度. Latitude
      & z_DelSigma            ! $ \Delta \sigma $ (整数). 
                              ! $ \Delta \sigma $ (full)

    ! 積分と平均の操作
    ! Operation for integral and average
    !
#ifdef AXISYMMETRY
    use wa_zonal_module, only: IntLon_x
#elif LIB_MPI
    use ua_mpi_module, only: IntLon_x => IntLon_p
#else
    use wa_module, only: IntLon_x
#endif
      
    real(DP), intent(in):: xyz_V (0:imax-1, 1:jmax, 1:kmax)
                              ! $ v $ . 南北風速. 
                              ! Northward wind
    real(DP), intent(in):: xy_Ps (0:imax-1, 1:jmax)
                              ! $ p_s $ . 地表面気圧. 
                              ! Surface pressure
    real(DP)            :: yr_MassStreamFunc(1:jmax, 0:kmax)
                              ! $ \Psi $ . 質量流線関数
                              ! Mass stream function

    integer  :: j, k

    yr_MassStreamFunc(:,kmax) = 0.0D0

    do j=1, jmax
       do k=kmax-1,0,-1
          yr_MassStreamFunc(j,k) &
               = yr_MassStreamFunc(j,k+1) &
               + IntLon_x(xy_Ps(:,j) * xyz_V(:,j,k+1)) * z_DelSigma(k+1)
       end do
       yr_MassStreamFunc(j,:) &
            = yr_MassStreamFunc(j,:) * Rplanet * cos(y_Lat(j))/Grav
    end do

  end function yr_MassStreamFunc
  
end module mass_streamfunc
