!----------------------------------------------------------------------
!     Copyright (c) 2011--2012 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!
!表題  wtq_module テストプログラム
!
!      磁場トロイダルポテンシャルの境界値問題
!
!履歴  2011/09/09  竹広真一   wq_module_tormagbc_test.f90 を MPI 化
!      2011/09/14  竹広真一   wtq_mpi_module 用に改造
!      2012/04/02  竹広真一   wtq_mpi_module_sjpack 用に改造
!      2012/04/03  竹広真一   wtq_module_sjpack 用に改造
!      2012/04/03  竹広真一   wtq_module 用に改造
!      2012/06/06  竹広真一   メッセージ修正
!
program wtq_module_tormagbc_wq_test

  use dc_message, only : MessageNotify
  use dc_test, only : AssertEqual
  use wtq_module

  implicit none

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

  real(8), dimension(0:im-1,1:jm,kmi)       :: xyr_TORMAG
  real(8), dimension(0:im-1,1:jm,kmi)       :: xyr_TORMAG_orig
  real(8), dimension((nm+1)*(nm+1),0:lmi)   :: wq_TORMAG
  real(8), dimension((nm+1)*(nm+1),0:lmi)   :: wq_TORMAG_orig
  real(8), dimension((nm+1)*(nm+1),kmi)     :: wr_TORMAG

  real(8), dimension(0:im-1,1:jm,kmi)       :: xyr_Zero = 0.0D0

  ! 判定誤差設定
  integer, parameter :: check_digits = 8
  integer, parameter :: ignore = -9

  real(8), parameter  :: pi=3.1415926535897932385D0

  call MessageNotify('M','wtq_module_tormagbc_wq_test', &
       'wtq_module wq_TormagBoundary subroutine test')

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

 !==================== wq_TormagBoundary =========================
  ! P_10
  xyr_TORMAG = sin(xyr_lat) * cos( pi*(xyr_rad-ri)/ri )*xyr_Rad
  ! P_1_1
  !xyr_TORMAG = cos(xyr_lat)*cos(xyr_lon)* cos( pi*(xyr_rad-ri)/(ro-ri) )
  !xyr_TORMAG = 2*sin(xyr_lat)**2 * cos( pi*(xyr_rad-ri)/(ro-ri) )

  wq_TORMAG = wq_xyr(xyr_TORMAG)
  wq_TORMAG_orig = wq_TORMAG

  call wq_TormagBoundary(wq_TORMAG)
  xyr_TORMAG = xyr_wq(wq_TORMAG)

  call AssertEqual(&
       message='wq_TormagBoundary (Top B.C.)',                   &
       answer = xyr_Zero(:,:,kmi),                                   &
       check = xyr_TORMAG(:,:,kmi),                                  &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  call AssertEqual(&
       message='wq_TormagBoundary (internal value)',                 &
       answer = wq_Tormag(:,0:lmi-2),                                &
       check = wq_Tormag_orig(:,0:lmi-2),                            &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

 !==================== wq_TormagBoundaryGrid =========================
  ! P_10
  !xyr_TORMAG = sin(xyr_lat) * cos( pi*(xyr_rad-ri)/ri )
  ! P_1_1
  xyr_TORMAG = cos(xyr_lat)*cos(xyr_lon)* cos( pi*(xyr_rad-ri)/ri )*xyr_Rad
  !xyr_TORMAG = 2*sin(xyr_lat)**2 * cos( pi*(xyr_rad-ri)/ri )

  xyr_TORMAG_orig = xyr_TORMAG
  wr_TORMAG = wr_xyr(xyr_TORMAG)
  call wr_TormagBoundaryGrid(wr_TORMAG)
  xyr_TORMAG = xyr_wr(wr_TORMAG)

  call AssertEqual(&
       message='wr_TormagBoundaryGrid (Top B.C.)',                   &
       answer = xyr_Zero(:,:,kmi),                                   &
       check = xyr_TORMAG(:,:,kmi),                                  &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  call AssertEqual(&
       message='wr_TormagBoundaryGrid (internal value)',             &
       answer = xyr_TORMAG_orig(:,:,1:kmi-1),                        &
       check = xyr_TORMAG(:,:,1:kmi-1),                              &
       significant_digits = check_digits, ignore_digits = ignore     &
       )


  call MessageNotify('M','wtq_module_tormagbc_wq_test', &
       'wtq_module wq_TormagBoundary subroutine test succeeded!')

end program wtq_module_tormagbc_wq_test