!----------------------------------------------------------------------
!     Copyright (c) 2012--2013 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!表題  wa_base_mpi_module_supack
!
!  spml/wa_base_mpi_module_supack モジュールは球面上での流体運動を
!  球面調和函数を用いたスペクトル法と MPI 並列化によって
!  数値計算するための モジュール wa_mpi_module_supack の下部モジュールであり, 
!  スペクトル計算の基本的な Fortran90 関数を提供する. 
!
!  球面上の 1 層モデル用 w_base_mpi_module_supack モジュールを多層モデル用に
!  拡張したものであり, 同時に複数個のスペクトルデータ, 格子点データに
!  対する変換が行える.
!
!  内部で ISPACK2 の SUPACK の Fortran77 サブルーチンを呼んでいる. 
!  スペクトルデータおよび格子点データの格納方法や変換の詳しい計算法に
!  ついては ISPACK2/SUPACK のマニュアルを参照されたい.
!
!
!履歴  2016/03/15  竹広真一  wa_base_mpi_module_sjpack.f90 を supack 化
!
module wa_base_mpi_module_supack
  !
  ! wa_base_mpi_module_supack
  !
  !  spml/wa_base_mpi_module_supack モジュールは球面上での流体運動を
  !  球面調和函数を用いたスペクトル法と MPI 並列化によって
  !  数値計算するための モジュール wa_mpi_module_supack の下部モジュールであり, 
  !  スペクトル計算の基本的な Fortran90 関数を提供する. 
  !
  !  球面上の 1 層モデル用 w_base_mpi_module_supack モジュールを多層モデル用に
  !  拡張したものであり, 同時に複数個のスペクトルデータ, 格子点データに
  !  対する変換が行える.
  !
  !  内部で ISPACK2 の SUPACK の Fortran77 サブルーチンを呼んでいる. 
  !  スペクトルデータおよび格子点データの格納方法や変換の詳しい計算法に
  !  ついては ISPACK2/SUPACK のマニュアルを参照されたい.
  !
  !
  use dc_message
  use w_base_mpi_module_supack, only : im, nm=>nn, nc, jc=>jl, w_xv, xv_w, &
                                       w_StreamPotential2VectorMPI, &
                                       w_Vector2VorDivMPI,          &
                                       w_VectorCosLat2VorDivMPI

  implicit none

  private
  private im, nm                               ! Intel Fortran 対策

  public wa_base_mpi_Initial                   ! 初期化サブルーチン
  public wa_base_mpi_Finalize                  ! 終了処理

  public xva_wa, wa_xva                        ! 変換関数
  public wa_StreamPotential2VectorMPI          ! 流線ポテンシャルから速度場計算
  public wa_Vector2VorDivMPI                   ! 速度場から渦度発散を計算
  public wa_VectorCosLat2VorDivMPI             ! 速度場から渦度発散を計算

  contains
  !--------------- 初期化 -----------------
    subroutine wa_base_mpi_initial
      ! 
      ! SNPACK 用 wa_base_mpi_initial の互換のためのダミーサブルーチン
      !
      call MessageNotify('M','wa_base_mpi_initial',&
                         'wa_base_mpi_module_supack is initialized (dummy, 2016/03/15)')

    end subroutine wa_base_mpi_initial

  !--------------- 基本変換 -----------------

    function xva_wa(wa_data,ipow,iflag)    ! 球面調和関数スペクトル -> 格子点
      !
      ! スペクトルデータから格子データへ変換する(多層用).
      !
      real(8), intent(in)   :: wa_data(:,:)
      !(in) スペクトルデータ

      real(8)               :: xva_wa(0:im-1,jc,size(wa_data,2))
      !(out) 格子点データ

      integer, intent(in), optional  :: ipow
      !(in) 作用させる 1/cosφ の次数. 省略時は 0. 

      integer, intent(in), optional  :: iflag
      !(in) 変換の種類
      !     0 : 通常の正変換
      !     1 : 経度微分を作用させた正変換
      !    -1 : 緯度微分を作用させた正変換
      !     2 : sinφを作用させた正変換
      !     省略時は 0.
      !

      integer, parameter  :: ipow_default  = 0
      integer, parameter  :: iflag_default = 0

      integer ipval, ifval
      integer k

      if (present(ipow)) then
         ipval = ipow
      else
         ipval = ipow_default
      endif

      if (present(iflag)) then
         ifval = iflag
      else
         ifval = iflag_default
      endif

      do k=1,size(wa_data,2)
        xva_wa(:,:,k) = xv_w(wa_data(:,k),iflag=ifval,ipow=ipval)
      enddo

    end function xva_wa

    function wa_xva(xva_data,ipow,iflag) ! 格子点 -> 球面調和関数スペクトル
      !
      ! 格子データからスペクトルデータへ(正)変換する(多層用).
      !
      real(8), intent(in)   :: xva_data(:,:,:)
      !(in) 格子点データ(im,jm,*)

      real(8)               :: wa_xva(nc,size(xva_data,3))
      !(out) スペクトルデータ

      integer, intent(in), optional  :: ipow
      !(in) 変換時に同時に作用させる 1/cosφ の次数. 省略時は 0.

      integer, intent(in), optional  :: iflag
      ! 変換の種類
      !     0 : 通常の正変換
      !     1 : 経度微分を作用させた正変換
      !    -1 : 緯度微分を作用させた正変換
      !     2 : sinφを作用させた正変換
      !   省略時は 0.

      integer, parameter  :: ipow_default  = 0      ! スイッチデフォルト値
      integer, parameter  :: iflag_default = 0      ! スイッチデフォルト値

      integer ipval, ifval
      integer k

      if (present(ipow)) then
         ipval = ipow
      else
         ipval = ipow_default
      endif

      if (present(iflag)) then
         ifval = iflag
      else
         ifval = iflag_default
      endif

      do k=1,size(xva_data,3)
         wa_xva(:,k) = w_xv(xva_data(:,:,k),iflag=ifval,ipow=ipval)
      enddo

    end function wa_xva

  !----------- 速度, 渦度・発散, 流線・ポテンシャル計算 -------------

    subroutine wa_StreamPotential2VectorMPI(wa_Psi, wa_Chi, xva_U, xva_V)
      !
      ! 流線・ポテンシャル(スペクトルデータ)から速度場(格子データ)に
      ! (逆)変換する(多層用)
      !
      ! スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ. 
      !
      !   u cosφ =      ∂χ/∂λ - cosφ∂ψ/∂φ, 
      !   v cosφ = cosφ∂χ/∂φ +      ∂ψ/∂λ 
      !
      real(8), intent(in)   :: wa_Psi(:,:)
      !(in) 流線関数((nm+1)*(nm+1),:)

      real(8), intent(in)   :: wa_Chi(:,:)
      !(in) 速度ポテンシャル((nm+1)*(nm+1),:)

      real(8), intent(out)   :: xva_U(0:im-1,1:jc,size(wa_Psi,2))
      !(out) 速度経度成分(0:im-1,1:jc,:)

      real(8), intent(out)   :: xva_V(0:im-1,1:jc,size(wa_Psi,2))
      !(out) 速度緯度成分(0:im-1,1:jc,:)

      integer :: k

      do k=1,size(wa_Psi,2)
         call w_StreamPotential2VectorMPI &
              ( wa_Psi(:,k), wa_Chi(:,k), xva_U(:,:,k), xva_V(:,:,k) )
      enddo

    end subroutine wa_StreamPotential2VectorMPI

    subroutine wa_Vector2VorDivMPI(xva_U, xva_V, wa_Vor, wa_Div)
      !
      ! 速度場(格子データ)から渦度・発散(スペクトルデータ)に
      ! (正)変換する(多層用)
      ! 
      ! スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ. 
      !
      !   ζ = 1/cosφ∂v/∂λ - 1/cosφ ∂(u cosφ)/∂φ 
      !    D = 1/cosφ∂u/∂λ + 1/cosφ ∂(v cosφ)/∂φ
      !
      real(8), intent(in)   :: xva_U(0:,:,:)
      !(in) 速度経度成分(0:im-1,1:jm,:)

      real(8), intent(in)   :: xva_V(0:,:,:)
      !(in) 速度緯度成分(0:im-1,1:jm,:)

      real(8), intent(out)   :: wa_Vor(nc,size(xva_U,3))
      !(out) 流線関数
      real(8), intent(out)   :: wa_Div(nc,size(xva_U,3))
      !(out) 速度ポテンシャル

      integer :: k

      do k=1,size(xva_U,3)
         call w_Vector2VorDivMPI &
               (xva_U(:,:,k), xva_V(:,:,k), wa_Vor(:,k), wa_Div(:,k))
      enddo
    end subroutine wa_Vector2VorDivMPI

    subroutine wa_VectorCosLat2VorDivMPI(xva_UCosLat, xva_VCosLat, wa_Vor, wa_Div)
      !
      ! 速度場(格子データ)から渦度・発散(スペクトルデータ)に
      ! (正)変換する(多層用)
      ! 
      ! スペクトル変換を用いず微分を計算するために, 変換回数が 2 回ですむ. 
      !
      !   ζ = 1/cosφ∂v/∂λ - 1/cosφ ∂(u cosφ)/∂φ 
      !    D = 1/cosφ∂u/∂λ + 1/cosφ ∂(v cosφ)/∂φ
      !
      real(8), intent(in)   :: xva_UCosLat(0:,:,:)
      !(in) 速度経度成分 * cos(lat) (0:im-1,1:jm,:)

      real(8), intent(in)   :: xva_VCosLat(0:,:,:)
      !(in) 速度緯度成分 * cos(lat) (0:im-1,1:jm,:)

      real(8), intent(out)   :: wa_Vor(nc,size(xva_UCosLat,3))
      !(out) 流線関数
      real(8), intent(out)   :: wa_Div(nc,size(xva_UCosLat,3))
      !(out) 速度ポテンシャル

      integer :: k

      do k=1,size(xva_UCosLat,3)
         call w_VectorCosLat2VorDivMPI &
               ( xva_UCosLat(:,:,k), xva_VCosLat(:,:,k), &
                 wa_Vor(:,k), wa_Div(:,k) )
      enddo
    end subroutine wa_VectorCosLat2VorDivMPI

  !--------------- 終了処理 -----------------
    subroutine wa_base_mpi_Finalize
      ! 
      ! モジュールの終了処理(割り付け配列の解放)をおこなう. 
      !
      ! このサブルーチンを単独で用いるのでなく, 
      ! 上位サブルーチン wa_mpi_Finalize を使用すること.
      !
      call MessageNotify('M','wa_base_mpi_finalize',&
           'wa_base_mpi_module_supack is finalized (dummy, 2016/03/15)')

    end subroutine wa_base_mpi_Finalize

end module wa_base_mpi_module_supack
