!---------------------------------------------------------------------
!     Copyright (C) GFD Dennou Club, 2006. All rights reserved.
!---------------------------------------------------------------------
!
!= 2 次元 (xy 方向) 等間隔交互格子 有限差分モデル用 基本モジュール
!
!* 履歴
!  * 2007/07/15 (小高正嗣) : 3D deepconv へ移植, dc_types を Use.
!  * 2006/12/26 (小高正嗣) : 平均関数の名前を元に戻す
!  * 2006/06/17 (小高正嗣) : 積分, 平均関数を追加
!  * 2006/06/02 (小高正嗣) : 新規作成
!
module xy_base_module
  != 2 次元 (xy 方向) 等間隔交互格子 有限差分モデル用 基本モジュール
  !
  !== 概要
  !
  ! xy_base_module は, 2 次元 (xy 方向) 等間隔交互格子を用いた有限差分
  ! 法に基づく数値モデルのための, 基本的な Fortran90 副プログラムおよび
  ! 関数を提供する. 
  !
  ! このモジュールは xy_module の下位モジュールである. 下請けモジュール
  ! として data_type, x_base_module, y_base_module モジュールを用いている. 
  ! 
  !
  !== 備考
  !
  ! * 例外処理がない
  ! * その場合のメッセージ表示がない
  ! * 初期化の値はマシンイプシロン値としておくべき
  !

  use dc_types,      only : DBKIND => DP
  use x_base_module, only : im, imin, imax, xmargin, x_axis_init,   &
    &                       x_X, p_X, x_dx, p_dx, x_avr_p, p_avr_x, &
    &                       IntX_p, IntX_x, AvrX_p, AvrX_x
  use y_base_module, only : jm, jmin, jmax, ymargin, y_axis_init,   &
    &                       y_Y, q_Y, y_dy, q_dy, y_avr_q, q_avr_y, &
    &                       IntY_q, IntY_y, AvrY_q, AvrY_y
  implicit none

  private
  public :: im, imin, imax, xmargin, x_X, p_X, x_dx, p_dx     ! Cascaded
  public :: jm, jmin, jmax, ymargin, y_Y, q_Y, y_dy, q_dy     ! Cascaded
  public :: x_avr_p, p_avr_x, IntX_p, IntX_x, AvrX_p, AvrX_x  ! Cascaded
  public :: y_avr_q, q_avr_y, IntY_q, IntY_y, AvrY_q, AvrY_y  ! Cascaded
  public :: xy_X, xy_Y, xy_dx, py_dx, xy_dy, xq_dy
  public :: xy_axis_init
  public :: xy_avr_py, xq_avr_pq, py_avr_xy, pq_avr_xq
  public :: xy_avr_xq, py_avr_pq, xq_avr_xy, pq_avr_py
  public :: y_IntX_xy, q_IntX_xq, y_IntX_py, x_IntY_xy, p_IntY_py, x_IntY_xq
  public :: y_AvrX_xy, q_AvrX_xq, y_AvrX_py, x_AvrY_xy, p_AvrY_py, x_AvrY_xq
  public :: IntXY_xy, IntXY_xq, IntXY_py, AvrXY_xy, AvrXY_py, AvrXY_xq


  real(DBKIND),allocatable :: xy_X(:,:)   ! x 座標(半整数格子)
  real(DBKIND),allocatable :: xy_Y(:,:)   ! y 座標(半整数格子)
  real(DBKIND),allocatable :: xy_dx(:,:)  ! x 方向格子間隔(半整数格子)
  real(DBKIND),allocatable :: py_dx(:,:)  ! x 方向格子間隔(整数格子)
  real(DBKIND),allocatable :: xy_dy(:,:)  ! y 方向格子間隔(半整数格子)
  real(DBKIND),allocatable :: xq_dy(:,:)  ! y 方向格子間隔(整数格子)

  save xy_X, xy_Y, xy_dx, py_dx, xy_dy, xq_dy

  interface xy_avr_py
    module procedure xa_avr_pa
  end interface 

  interface xq_avr_pq
    module procedure xa_avr_pa
  end interface 

  interface py_avr_xy
    module procedure pa_avr_xa
  end interface 

  interface pq_avr_xq
    module procedure pa_avr_xa
  end interface 

  interface xy_avr_xq
    module procedure ay_avr_aq
  end interface 

  interface py_avr_pq
    module procedure ay_avr_aq
  end interface 

  interface xq_avr_xy
    module procedure aq_avr_ay
  end interface 

  interface pq_avr_py
    module procedure aq_avr_ay
  end interface 

  interface y_IntX_xy
    module procedure a_IntX_xa
  end interface 

  interface q_IntX_xq
    module procedure a_IntX_xa
  end interface 

  interface x_IntY_xy
    module procedure a_IntY_ay
  end interface 

  interface p_IntY_py
    module procedure a_IntY_ay
  end interface 

  interface y_AvrX_xy
    module procedure a_AvrX_xa
  end interface 

  interface q_AvrX_xq
    module procedure a_AvrX_xa
  end interface 

  interface x_AvrY_xy
    module procedure a_AvrY_ay
  end interface 

  interface p_AvrY_py
    module procedure a_AvrY_ay
  end interface 

  contains
!--------------------------------------------------------------------
    subroutine xy_axis_init(i, j, xmg, ymg, xmin, xmax, ymin, ymax)
      ! xy 方向の座標値と格子点間隔を設定する

      integer,intent(in) :: i     ! x 方向格子点数
      integer,intent(in) :: j     ! y 方向格子点数
      integer,intent(in) :: xmg   ! x 方向糊代格子点数
      integer,intent(in) :: ymg   ! y 方向糊代格子点数
      real(DBKIND),intent(in)     :: xmin  ! x 座標最小値     
      real(DBKIND),intent(in)     :: xmax  ! x 座標最大値  
      real(DBKIND),intent(in)     :: ymin  ! y 座標最小値     
      real(DBKIND),intent(in)     :: ymax  ! y 座標最大値  

      ! 配列の上下限の値, 座標値と格子点間隔を設定
      ! * 1 次元用のサブルーチンを用いる
      !
      call x_axis_init(i, xmg, xmin, xmax)
      call y_axis_init(j, ymg, ymin, ymax)

      ! 2 次元座標配列の設定
      ! * 組み込み関数 spread を用いて 1 次元配列を 2 次元に拡張する.
      ! 
      allocate(xy_X(imin:imax,jmin:jmax))
      allocate(xy_Y(imin:imax,jmin:jmax))
      allocate(xy_dx(imin:imax,jmin:jmax))
      allocate(py_dx(imin:imax,jmin:jmax))
      allocate(xy_dy(imin:imax,jmin:jmax))
      allocate(xq_dy(imin:imax,jmin:jmax))

      xy_X  = spread(x_X,2,size(y_Y))
      xy_Y  = spread(y_Y,1,size(x_X))
      xy_dx = spread(x_dx,2,size(y_Y))
      py_dx = spread(p_dx,2,size(y_Y))
      xy_dy = spread(y_dy,1,size(x_X))
      xq_dy = spread(q_dy,1,size(x_X))

    end subroutine xy_axis_init  
!--------------------------------------------------------------------
    function xa_avr_pa(pa_Var)
      ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す
  
      real(DBKIND),intent(in) :: pa_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: xa_avr_pa(imin:imax,jmin:jmax) ! 出力
      integer                 :: jy                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      xa_avr_pa = 0.0d0

      ! 平均操作
      ! * 関数 x_avr_p を用いて計算
      !
      do jy = jmin, jmax
        xa_avr_pa(:,jy) = x_avr_p(pa_Var(:,jy))
      end do

    end function xa_avr_pa
!--------------------------------------------------------------------
    function pa_avr_xa(xa_Var)
      ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す
  
      real(DBKIND),intent(in) :: xa_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: pa_avr_xa(imin:imax,jmin:jmax) ! 出力
      integer                 :: jy                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      pa_avr_xa = 0.0d0

      ! 平均操作
      ! * 関数 p_x を用いて計算
      !
      do jy = jmin, jmax
        pa_avr_xa(:,jy) = p_avr_x(xa_Var(:,jy))
      end do

    end function pa_avr_xa
!--------------------------------------------------------------------
    function ay_avr_aq(aq_Var)
      ! 平均操作を行い整数格子点の配列値を半整数格子点上へ返す
  
      real(DBKIND),intent(in) :: aq_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: ay_avr_aq(imin:imax,jmin:jmax) ! 出力
      integer                 :: ix                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      ay_avr_aq = 0.0d0

      ! 平均操作
      ! * 関数 y_q を用いて計算
      !
      do ix = imin, imax
        ay_avr_aq(ix,:) = y_avr_q(aq_Var(ix,:))
      end do

    end function ay_avr_aq
!--------------------------------------------------------------------
    function aq_avr_ay(ay_Var)
      ! 平均操作を行い半整数格子点の配列値を整数格子点上へ返す
  
      real(DBKIND),intent(in) :: ay_Var(imin:imax,jmin:jmax)   ! 入力
      real(DBKIND)            :: aq_avr_ay(imin:imax,jmin:jmax) ! 出力
      integer                 :: ix                            ! ループ添字

      ! 初期化
      ! * 0 割りを防ぐためにはマシンイプシロン値を用いるべき
      !
      aq_avr_ay = 0.0d0

      ! 平均操作
      ! * 関数 q_y を用いて計算.
      ! 
      do ix = imin, imax
        aq_avr_ay(ix,:) = q_avr_y(ay_Var(ix,:))
      end do

    end function aq_avr_ay
!--------------------------------------------------------------------
    function a_IntX_xa(xa_Var)
      ! xa 格子上の配列に対し x 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: xa_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_IntX_xa(jmin:jmax)        ! 出力
      integer                  :: jy                          ! ループ添字

      ! 初期化
      a_IntX_xa = 0.0d0

      ! 積分
      do jy = jmin, jmax
        a_IntX_xa(jy) = IntX_x(xa_Var(:,jy))
      end do
      
    end function a_IntX_xa
!--------------------------------------------------------------------
    function y_IntX_py(py_Var)
      ! py 格子上の配列に対し x 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: y_IntX_py(jmin:jmax)        ! 出力
      integer                  :: jy                          ! ループ添字

      ! 初期化
      y_IntX_py = 0.0d0

      ! 積分
      do jy = jmin, jmax
        y_IntX_py(jy) = IntX_p(py_Var(:,jy))
      end do
      
    end function y_IntX_py
!--------------------------------------------------------------------
    function a_IntY_ay(ay_Var)
      ! ay 格子上の配列に対し y 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: ay_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_IntY_ay(imin:imax)        ! 出力
      integer                  :: ix                          ! ループ添字

      ! 初期化
      a_IntY_ay = 0.0d0

      ! 積分
      do ix = imin, imax
        a_IntY_ay(ix) = IntY_y(ay_Var(ix,:))
      end do
      
    end function a_IntY_ay
!--------------------------------------------------------------------
    function x_IntY_xq(xq_Var)
      ! ay 格子上の配列に対し y 方向に重み付きの積分を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: x_IntY_xq(imin:imax)        ! 出力
      integer                  :: ix                          ! ループ添字

      ! 初期化
      x_IntY_xq = 0.0d0

      ! 積分
      do ix = imin, imax
        x_IntY_xq(ix) = IntY_q(xq_Var(ix,:))
      end do
      
    end function x_IntY_xq
!--------------------------------------------------------------------
    function a_AvrX_xa(xa_Var)
      ! xa 格子上の配列に対し x 方向に平均を行う

      real(DBKIND), intent(in) :: xa_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_AvrX_xa(jmin:jmax)        ! 出力


      a_AvrX_xa = a_IntX_xa(xa_Var)/sum(x_dx(1:im))
      
    end function a_AvrX_xa
!--------------------------------------------------------------------
    function y_AvrX_py(py_Var)
      ! py 格子上の配列に対し x 方向に平均を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: y_AvrX_py(jmin:jmax)        ! 出力


      y_AvrX_py = y_IntX_py(py_Var)/sum(p_dx(1:im))
      
    end function y_AvrX_py
!--------------------------------------------------------------------
    function a_AvrY_ay(ay_Var)
      ! ay 格子上の配列に対し y 方向に平均を行う

      real(DBKIND), intent(in) :: ay_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: a_AvrY_ay(imin:imax)        ! 出力


      a_AvrY_ay = a_IntY_ay(ay_Var)/sum(y_dy(1:jm))
      
    end function a_AvrY_ay
!--------------------------------------------------------------------
    function x_AvrY_xq(xq_Var)
      ! xq 格子上の配列に対し y 方向に平均を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: x_AvrY_xq(imin:imax)        ! 出力


      x_AvrY_xq = x_IntY_xq(xq_Var)/sum(q_dy(1:jm))
      
    end function x_AvrY_xq
!--------------------------------------------------------------------
    function IntXY_xy(xy_Var)
      ! xy 格子上の配列に対し領域積分を行う

      real(DBKIND), intent(in) :: xy_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: IntXY_xy                    ! 出力


      IntXY_xy = IntY_y(a_IntX_xa(xy_Var))

    end function IntXY_xy
!--------------------------------------------------------------------
    function IntXY_xq(xq_Var)
      ! xq 格子上の配列に対し領域積分を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: IntXY_xq                    ! 出力


      IntXY_xq = IntY_q(a_IntX_xa(xq_Var))

    end function IntXY_xq
!--------------------------------------------------------------------
    function IntXY_py(py_Var)
      ! py 格子上の配列に対し領域積分を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: IntXY_py                    ! 出力


      IntXY_py = IntY_y(y_IntX_py(py_Var))

    end function IntXY_py
!--------------------------------------------------------------------
    function AvrXY_xy(xy_Var)
      ! xy 格子上の配列に対し領域平均を行う

      real(DBKIND), intent(in) :: xy_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: AvrXY_xy                    ! 出力


      AvrXY_xy = IntXY_xy(xy_Var)/(sum(x_dx(1:im))*sum(y_dy(1:jm)))

    end function AvrXY_xy
!--------------------------------------------------------------------
    function AvrXY_xq(xq_Var)
      ! xq 格子上の配列に対し領域平均を行う

      real(DBKIND), intent(in) :: xq_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: AvrXY_xq                    ! 出力


      AvrXY_xq = IntXY_xq(xq_Var)/(sum(x_dx(1:im))*sum(q_dy(1:jm)))

    end function AvrXY_xq
!--------------------------------------------------------------------
    function AvrXY_py(py_Var)
      ! py 格子上の配列に対し領域平均を行う

      real(DBKIND), intent(in) :: py_Var(imin:imax,jmin:jmax) ! 入力
      real(DBKIND)             :: AvrXY_py                    ! 出力


      AvrXY_py = IntXY_py(py_Var)/(sum(p_dx(1:im))*sum(y_dy(1:jm)))

    end function AvrXY_py
!--------------------------------------------------------------------
end module xy_base_module
