!----------------------------------------------------------------------
!     Copyright (c) 2012 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!
!表題  wa_deriv_mpi_module テストプログラム :: 微分関数のテスト
!
!履歴  2012/04/03  竹広真一
!
program wa_deriv_mpi_module_test

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

  integer, parameter :: im=32, jm=16, nm=10, km=2

  real(8), allocatable  ::  xva_data1(:,:,:)              ! 元の関数
  real(8), allocatable  ::  xva_data2(:,:,:)              ! 元の関数
  real(8), allocatable  ::  xva_ddata(:,:,:)              ! 微分の正解
  real(8), allocatable  ::  xv_mu(:,:)                    ! μ=sinφ

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

  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','wa_deriv_mpi_module_test', &
                         'wa_deriv_mpi_module function tests') 

  call wa_mpi_Initial( nm, im, jm, km )

  allocate(xva_data1(0:im-1,jc,km))
  allocate(xva_data2(0:im-1,jc,km))
  allocate(xva_ddata(0:im-1,jc,km))
  allocate(xv_mu(0:im-1,jc))

  !============== 微分計算 (lon,lat座標系用) のテスト ==============

  !---- Y_1^-1 Y_2^1 のテスト ----
  xva_data1(:,:,1) = -cos(xv_Lat)*sin(xv_Lon)             ! Y_1^{-1}
  xva_data1(:,:,2)  = sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon) ! Y_2^1

  xva_ddata(:,:,1) = 2*cos(xv_Lat)*sin(xv_Lon)             ! wa_Lapla_wa
  xva_ddata(:,:,2) = -6*sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon) ! Y_2^1
  call AssertEqual(&
    message='wa_Lapla_wa with Y_1^-1 and Y^2_1',                  &
    answer = xva_wa(wa_Lapla_wa(wa_xva(xva_data1))),              &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = 1.0/2.0*cos(xv_Lat)*sin(xv_Lon)      ! wa_LaplaInv_wa
  xva_ddata(:,:,2) = -1.0D0/6.0*sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon)
  call AssertEqual(&
    message='wa_LaplaInv_wa with Y_1^-1 and Y^2_1',               &
    answer = xva_wa(wa_LaplaInv_wa(wa_xva(xva_data1))),           &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = -cos(xv_Lat)*cos(xv_Lon)         ! wa_DLon_wa
  xva_ddata(:,:,2) = -sin(xv_Lat)*cos(xv_Lat) * sin(xv_Lon)
  call AssertEqual(&
    message='wa_DLon_wa with Y_1^-1 and Y^2_1',                   &
    answer = xva_wa(wa_DLon_wa(wa_xva(xva_data1))),               &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = -cos(xv_Lon)                     ! xva_GradLon_wa
  xva_ddata(:,:,2) = -sin(xv_Lat) * sin(xv_Lon)
  call AssertEqual(&
    message='xva_GradLon_wa with Y_1^-1 and Y^2_1',               &
    answer = xva_GradLon_wa(wa_xva(xva_data1)),                   &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = sin(xv_Lat)*sin(xv_Lon)         ! xva_GradLat_wa
  xva_ddata(:,:,2) = cos(2*xv_Lat) * cos(xv_Lon)
  call AssertEqual(&
    message='xva_GradLat_wa with Y_1^-1 and Y^2_1',               &
    answer = xva_GradLat_wa(wa_xva(xva_data1)),                   &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  !---- Y_1^-1 cosφ, Y_2^1 cosφ のテスト ----
  xva_data1(:,:,1) = -cos(xv_Lat)**2*sin(xv_Lon)              ! Y_1^-1 cosφ
  xva_data1(:,:,2) = sin(xv_Lat)*cos(xv_Lat)**2 * cos(xv_Lon) ! Y_2^1 cosφ

  xva_ddata(:,:,1) = -cos(xv_Lat)*cos(xv_Lon)                ! wa_DivLon_xv
  xva_ddata(:,:,2) = -sin(xv_Lat)*cos(xv_Lat)*sin(xv_Lon)    ! wa_DivLon_xv
  call AssertEqual(&
    message='wa_DivLon_xva with Y_1^-1 cosφ and Y^2_1 cosφ',    &
    answer = xva_wa(wa_DivLon_xva(xva_data1)),                    &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = 3*sin(xv_Lat)*cos(xv_Lat)*sin(xv_Lon)        !wa_DivLat_wa
  xva_ddata(:,:,2) = cos(xv_Lat)*(1-4*sin(xv_Lat)**2)*cos(xv_Lon) !wa_DivLat_wa
  call AssertEqual(&
    message='wa_DivLat_xva with Y_1^-1 cosφ and Y^2_1 cosφ',    &
    answer = xva_wa(wa_DivLat_xva(xva_data1)),                    &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  !---- Jacobian のテスト ----
  xva_data1(:,:,1) = -cos(xv_Lat)*sin(xv_Lon)                ! Y_1^{-1}
  xva_data1(:,:,2) = -cos(xv_Lat)*sin(xv_Lon)                ! Y_1^{-1}
  xva_data2(:,:,1) = -cos(xv_Lat)*sin(xv_Lon)                ! Y_1^{-1}
  xva_data2(:,:,2) = sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon)   ! Y_2^1

  xva_ddata(:,:,1) = 0.0
  xva_ddata(:,:,2) = sin(xv_Lat)**2 - cos(xv_Lat)**2*cos(xv_Lon)**2
   call AssertEqual(&
    message='wa_JacobianMPI_wa_wa Y_1^-1 and Y^2_1',                 &
    answer = xva_wa(wa_JacobianMPI_wa_wa(wa_xva(xva_data1),wa_xva(xva_data2))), &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )


  !============== 微分計算 (λ,μ座標系用) のテスト ==============
  xv_mu = sin(xv_Lat)

  !----- Y_2^0, Y_1^1 のテスト -----
  xva_data1(:,:,1) = 3*xv_mu**2-1                              ! Y_2^0
  xva_data1(:,:,2) = sqrt(1-xv_mu**2)*cos(xv_Lon)                  ! Y_1^1

  xva_ddata(:,:,1) = 0.0
  xva_ddata(:,:,2) = -sqrt(1-xv_mu**2)*sin(xv_Lon)
  call AssertEqual(&
    message='xva_GradLambda_wa with Y_2^0 and Y^1_1',             &
    answer = xva_GradLambda_wa(wa_xva(xva_data1)),                &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = 6*xv_mu*(1-xv_mu**2)
  xva_ddata(:,:,2) = -xv_mu*sqrt(1-xv_mu**2)*cos(xv_Lon)
  call AssertEqual(&
    message='xva_GradMu_wa with Y_2^0 and Y^1_1',                 &
    answer = xva_GradMu_wa(wa_xva(xva_data1)),                    &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  !----- Y_2^0(1-μ^2), Y_1^1 (1-μ^2) のテスト -----
  xva_data1(:,:,1) = (3*xv_mu**2-1)*(1-xv_mu**2)       ! Y_2^0 (1-μ^2)
  xva_data1(:,:,2) = (1-xv_mu**2)**(3.0/2)*cos(xv_Lon) ! Y_1^1 (1-μ^2)

  xva_ddata(:,:,1) = 0.0
  xva_ddata(:,:,2) = -(1-xv_mu**2)**(1.0d0/2)*sin(xv_Lon)
  call AssertEqual(&
    message='wa_DivLambda_xva with Y_2^0 (1-μ^2) and Y^1_1 √(1-μ^2)',&
    answer = xva_wa(wa_DivLambda_xva(xva_data1)),                 &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xva_ddata(:,:,1) = (2-3*xv_mu**2)*4*xv_mu
  xva_ddata(:,:,2) = -3.0D0*xv_mu*(1-xv_mu**2)**(1.0D0/2)*cos(xv_Lon)
  call AssertEqual(&
    message='wa_DivMu_xva with Y_2^0 (1-μ^2) and Y^1_1 √(1-μ^2)',&
    answer = xva_wa(wa_DivMu_xva(xva_data1)),                     &
    check = xva_ddata,                                            &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  call MessageNotify('M','wa_deriv_mpi_module_test', &
                         'wa_deriv_mpi_module function tests succeeded!') 

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

  call MPI_FINALIZE(IERR)

end program wa_deriv_mpi_module_test

