!----------------------------------------------------------------------
!     Copyright (c) 2016 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!
!表題  ut_mpi_module_supack テストプログラム
!
!      微分のテスト
!  
!履歴  2016/08/04  竹広真一   wt_module_supack_test_derivavit1.f90 より改造
!
program ut_mpi_module_svpack_derivative1_test

  use dc_message, only : MessageNotify
  use dc_test, only : AssertEqual
  use ut_mpi_module_svpack
  use mpi
  implicit none

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

  real(8), allocatable              :: xyh_Data(:,:,:)
  real(8), allocatable              :: xyh_DData(:,:,:)
                                                      
  real(8), allocatable              :: xyh_VLon(:,:,:)
  real(8), allocatable              :: xyh_VLat(:,:,:)
  real(8), allocatable              :: xyh_VRad(:,:,:)
  real(8), allocatable              :: xyh_Div(:,:,:)
  real(8), allocatable              :: xyh_Psi(:,:,:)

  ! 判定誤差設定
  integer, parameter :: check_digits = 10
  integer, parameter :: ignore = -11

  integer, parameter :: nmin=1, nmax=9

  integer :: n

  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','ut_mpi_module_svpack_derivative_test1', &
       'ut_mpi_module_svpack derivative function test #1')

  call ut_mpi_Initial(im,jm,km,nm,lm,ri,ro)

  allocate(xyh_Data(0:im-1,1:jm,kc))
  allocate(xyh_DData(0:im-1,1:jm,kc))

  allocate(xyh_VLon(0:im-1,1:jm,kc))
  allocate(xyh_VLat(0:im-1,1:jm,kc))
  allocate(xyh_VRad(0:im-1,1:jm,kc))
  allocate(xyh_Div(0:im-1,1:jm,kc))
  allocate(xyh_Psi(0:im-1,1:jm,kc))

  !==================== ut_DivRad_ut ====================
  do n=nmin,nmax
     write(6,*) 'n=',n

     xyh_Data = xyh_Rad**n
     xyh_DData = (n+2)*xyh_Rad**(n-1)

     call AssertEqual(&
       message='ut_Divrat_ut with r^n',                              &
       answer = xyh_DData,                                           &
       check = xyh_ut(ut_DivRad_ut(ut_xyh(xyh_Data))),               &
       significant_digits = check_digits, ignore_digits = ignore     &
       )
  enddo

  !========== xyh_GradLon_ut, xyh_GradLat_ut, ut_Div_xyh_xyh_xyh ==========

  ! ------ r^2 Y_1^1 Toridal field -----
  n = 2
  write(6,*) 'n=',n
  xyh_Psi = xyh_Rad**n * cos(xyh_Lat)*sin(xyh_Lon)   ! r^2 Y_1^1
  xyh_VRad = 0

  xyh_VLon =  xyh_Rad**(n-1)*cos(xyh_Lon)
  call AssertEqual(&
       message='xyh_GradLon_ut with r^2 Y_1^-1',                     &
       answer = xyh_VLon,                                            &
       check = xyh_GradLon_ut(ut_xyh(xyh_Psi)),                      &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  xyh_VLat = -xyh_Rad**(n-1)*sin(xyh_Lat)*sin(xyh_Lon)
  call AssertEqual(&
       message='xyh_GradLat_ut with r^2 Y_1^-1',                     &
       answer = xyh_VLat,                                            &
       check = xyh_GradLat_ut(ut_xyh(xyh_Psi)),                     &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  xyh_Div = - 2 * xyh_Psi/xyh_Rad**2
  call AssertEqual(&
       message='ut_Div_xyh_xyh_xyh_ut with Toroiral r^2 Y_1^-1',     &
       answer = xyh_Div,                                             &
       check = xyh_ut(ut_Div_xyh_xyh_xyh(xyh_VLon,xyh_VLat,xyh_VRad)),&
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  ! ------ r^2Y_2^-1 Toridal field -----
  n = 2
  xyh_VRad = 0
  xyh_Psi = xyh_Rad**n * cos(xyh_Lat)*sin(xyh_Lat) * sin(xyh_Lon) ! r^2Y_2^-1

  xyh_VLon =  xyh_Rad**(n-1)*sin(xyh_Lat)*cos(xyh_Lon)
  call AssertEqual(&
       message='xyh_GradLon_ut with r^2 Y_2^-1',                     &
       answer = xyh_VLon,                                            &
       check =  xyh_GradLon_ut(ut_xyh(xyh_Psi)),                     &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  xyh_VLat = xyh_Rad**(n-1)*cos(2*xyh_Lat)*sin(xyh_Lon)
  call AssertEqual(&
       message='xyh_GradLat_ut with r^2 Y_2^-1',                     &
       answer = xyh_VLat,                                            &
       check = xyh_GradLat_ut(ut_xyh(xyh_Psi)),                      &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  xyh_Div = - 6* xyh_Psi/xyh_Rad**2
  call AssertEqual(&
       message='ut_Div_xyh_xyh_xyh_ut with Toroiral r2 Y_1^-1',      &
       answer = xyh_Div,                                             &
       check = xyh_ut(ut_Div_xyh_xyh_xyh(xyh_VLon,xyh_VLat,xyh_VRad)),&
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  !==================== ut_DivLon_xyz, ut_DivLat_xyz ====================
  n=2
  write(6,*)'n =',n

  xyh_VLon  = xyh_rad**n * cos(xyh_Lat)**2*sin(xyh_Lon)
  xyh_DData = xyh_rad**(n-1)*cos(xyh_Lat)*cos(xyh_Lon)
  call AssertEqual(&
       message='ut_DivLon_xyh with r^n cos^2 φ sinλ',              &
       answer = xyh_DData,                                           &
       check = xyh_ut(ut_DivLon_xyh(xyh_VLon)),                      &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  xyh_VLat   = xyh_rad**n * cos(xyh_Lat)**2*sin(xyh_Lon)
  xyh_DData = -3*xyh_rad**(n-1)*cos(xyh_Lat)*sin(xyh_Lat)*sin(xyh_Lon)
  call AssertEqual(&
       message='ut_DivLat_xyh with r^n cos^2 φ sinλ',              &
       answer = xyh_DData,                                           &
       check = xyh_ut(ut_DivLat_xyh(xyh_VLat)),                      &
       significant_digits = check_digits, ignore_digits = ignore     &
       )

  call MessageNotify('M','ut_mpi_module_supack_derivative1_test', &
       'ut_mpi_module_supack derivative function test #1 succeeded!')

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

  call MPI_FINALIZE(IERR)

end program ut_mpi_module_svpack_derivative1_test