!----------------------------------------------------------------------
!     Copyright (c) 2016 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!
!表題  wtq_mpi_module_supack テストプログラム
!
!      磁場トロイダルポテンシャルの境界値問題
!
!履歴  2016/04/04  竹広真一  wtq_mpi_module_sjpack_tormagbctau.f90 より改変
!
program wtq_mpi_supack_tormagbctau

  use dc_message, only : MessageNotify
  use dc_test, only : AssertEqual
  use mpi
  use wtq_mpi_module_supack

  implicit none

  integer,parameter  :: im=32, jm=16         ! 格子点の設定(経度, 緯度, 動径)
  integer,parameter  :: kmo=16, kmi=8        ! 格子点の設定(球殻動径, 球動径)
  integer,parameter  :: nm=10                ! 切断波数の設定(水平)
  integer,parameter  :: lmo=10, lmi=15       ! 切断波数の設定(球殻動径, 球動径)
  real(8),parameter  :: ri=0.5, ro=1.5       ! 内外半径

  real(8), allocatable  :: xvz_TORMAG(:,:,:)
  real(8), allocatable  :: wt_TORMAG(:,:)
  real(8), allocatable  :: xvr_TORMAG(:,:,:)
  real(8), allocatable  :: wq_TORMAG(:,:)
  real(8), allocatable  :: xvz_DTORDR(:,:,:)
  real(8), allocatable  :: xvr_DTORDR(:,:,:)

  real(8), allocatable  :: xv_Null(:,:)
  
  real(8), parameter  :: Pmo = 1.0           ! 球殻の磁気プランドル数
  real(8), parameter  :: Pmi = 2.0           ! 内球の磁気プランドル数

  real(8), parameter  :: pi=3.1415926535897932385D0
  ! 判定誤差設定
  integer, parameter :: check_digits = 11
  integer, parameter :: ignore = -12

  integer :: iproc, np, ierr

 !---------------- MPI スタート ---------------------
  call MPI_INIT(IERR)
  call MPI_COMM_RANK(MPI_COMM_WORLD,IPROC,IERR)
  call MPI_COMM_SIZE(MPI_COMM_WORLD,NP,IERR)

  call MessageNotify('M','wtq_mpi_module_supack_tormagbctau_test', &
       'wq_mpi_module_supack  wtq_TormagBoundariesTau subroutine test')

  call wtq_mpi_Initial(im,jm,kmi,kmo,nm,lmi,lmo,ri,ro)

  allocate(xvz_TORMAG(0:im-1,jc,0:kmo))
  allocate(wt_TORMAG(nc,0:lmo))
  allocate(xvr_TORMAG(0:im-1,jc,kmi))
  allocate(wq_TORMAG(nc,0:lmi))
  allocate(xvz_DTORDR(0:im-1,jc,0:kmo))
  allocate(xvr_DTORDR(0:im-1,jc,kmi))

  allocate(xv_Null(0:im-1,jc))

  xv_Null=0.0D0
  
  !==================== wtq_TormagBoundariesTau =========================
  ! P_10

  xvz_TORMAG = sin(xvz_lat) * cos( pi*(xvz_rad-ri)/(ro-ri) )
  xvr_TORMAG = sin(xvr_lat) * cos( pi*(xvr_rad-ri)/ri ) * xvr_Rad

  wt_TORMAG = wt_xvz(xvz_TORMAG)
  wq_TORMAG = wq_xvr(xvr_TORMAG)
  call wtq_TormagBoundariesTau(wt_TORMAG,wq_TORMAG,Pmo,Pmi)
  xvz_TORMAG = xvz_wt(wt_TORMAG)
  xvr_TORMAG = xvr_wq(wq_TORMAG)

  xvz_DTORDR = xvz_wt(wt_DRad_wt(wt_TORMAG))
  xvr_DTORDR = xvr_wq(wq_RadDRad_wq(wq_TORMAG))/xvr_Rad

  call checkresults

  ! P_1_1
  xvz_TORMAG = cos(xvz_lat)*cos(xvz_lon)* cos( pi*(xvz_rad-ri)/(ro-ri) )
  xvr_TORMAG = cos(xvr_lat)*cos(xvr_lon)* cos( pi*(xvr_rad-ri)/ri )*xvr_Rad

  wt_TORMAG = wt_xvz(xvz_TORMAG)
  wq_TORMAG = wq_xvr(xvr_TORMAG)
  call wtq_TormagBoundariesTau(wt_TORMAG,wq_TORMAG,Pmo,Pmi)
  xvz_TORMAG = xvz_wt(wt_TORMAG)
  xvr_TORMAG = xvr_wq(wq_TORMAG)

  xvz_DTORDR = xvz_wt(wt_DRad_wt(wt_TORMAG))
  xvr_DTORDR = xvr_wq(wq_RadDRad_wq(wq_TORMAG))/xvr_Rad

  call checkresults

  call MessageNotify('M','wtq_mpi_module_supack_tormagbctau_test', &
       'wtq_mpi_module_supack  wtq_TormagBoundariesTau subroutine test succeded')

 !------ MPIの終了 ------

  call MPI_FINALIZE(IERR)

contains

  subroutine checkresults

    call AssertEqual(&
      message='Top B.C',                                            &
      answer = xvz_TORMAG(:,:,0),                                   &
      check = xv_Null,                                              &
      significant_digits = check_digits, ignore_digits = ignore     &
      )

    call AssertEqual(&
      message='Inner B.C (value conituity) ',                       &
      answer = xvz_TORMAG(:,:,kmo),                                 &
      check = xvr_TORMAG(:,:,kmi),                                  &
      significant_digits = check_digits, ignore_digits = ignore     &
      )

    call AssertEqual(&
      message='Inner B.C (derivative conituity) ',                  &
      answer = Pmo*xvz_DTORDR(:,:,kmo),                             &
      check = Pmi*xvr_DTORDR(:,:,kmi),                              &
      significant_digits = check_digits, ignore_digits = ignore     &
      )
  end subroutine checkresults

end program wtq_mpi_supack_tormagbctau
