!----------------------------------------------------------------------
!     Copyright (c) 2019-2020 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!
!表題  w_base_module テストプログラム :: ベクトル用サブルーチンのテスト
!
!履歴  2019/03/13  竹広真一
!      2020/11/06  竹広真一 セクター計算テスト
!
program w_module_base_vector_mint_test

  use dc_message, only : MessageNotify
  use dc_test, only : AssertEqual
  use w_module_mint
  implicit none

!  integer, parameter :: im=64, jm=32, nm=21, mint=1
  integer, parameter :: im=32, jm=32, nm=21, mint=2  
  integer, parameter :: nc=(nm/mint+1)*(2*(nm+1)-nm/mint*mint)

  real(8), dimension(0:im-1,1:jm)  ::  xy_U              ! 速度経度成分
  real(8), dimension(0:im-1,1:jm)  ::  xy_V              ! 速度緯度成分

  real(8), dimension(0:im-1,1:jm)  ::  xy_UCosLat        ! 速度経度成分
  real(8), dimension(0:im-1,1:jm)  ::  xy_VCosLat        ! 速度緯度成分

  real(8), dimension(0:im-1,1:jm)  ::  xy_Uans           ! 速度経度成分
  real(8), dimension(0:im-1,1:jm)  ::  xy_Vans           ! 速度緯度成分

  real(8), dimension(0:im-1,1:jm)  ::  xy_Psi            ! 流線関数
  real(8), dimension(0:im-1,1:jm)  ::  xy_Chi            ! 速度ポテンシャル

  real(8), dimension(nc)           ::  w_Vor             ! 渦度
  real(8), dimension(nc)           ::  w_Div             ! 発散

  real(8), dimension(0:im-1,1:jm)  ::  xy_Vorans         ! 渦度
  real(8), dimension(0:im-1,1:jm)  ::  xy_Divans         ! 発散


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

  call MessageNotify('M','w_module_base_vector_test', &
                         'w_module_base vector subroutine tests')

  call w_Initial( nm, im, jm, mint=mint )

  !---- Y_1^0 のテスト ----
  xy_Psi = sin(xy_Lat)       ! Y_1^0
  xy_Chi = 0.0D0

  xy_Uans = -cos(xy_Lat)
  xy_Vans = 0.0D0

  xy_Vorans = -2*sin(xy_Lat)
  xy_Divans = 0.0D0

  call w_StreamPotential2Vector( w_xy(xy_Psi), w_xy(xy_Chi), xy_U, xy_V )

  call AssertEqual(&
    message='Psi = Y_1^0 test of w_StreamPotential2Vector(U)',   &
    answer = xy_Uans,                                             &
    check  = xy_U,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_1^0 test of w_StreamPotential2Vector(V)',   &
    answer = xy_Vans,                                             &
    check  = xy_V,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  call w_Vector2VorDiv( xy_U, xy_V, w_Vor, w_Div )

  call AssertEqual(&
    message='Psi = Y_1^0 test of w_Vector2VorDiv(Vor)',          &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_1^0 test of w_Vector2VorDiv(Div)',          &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_UCosLat = xy_U * cos(xy_Lat) ; xy_VCosLat = xy_V * cos(xy_Lat)

  call w_VectorCosLat2VorDiv( xy_UCosLat, xy_VCosLat, w_Vor, w_Div )

  call AssertEqual(&
    message='Psi = Y_1^0 test of w_VectorCosLat2VorDiv(Vor)',    &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_1^0 test of w_VectorCosLat2VorDiv(Div)',    &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_Psi = 0.0D0
  xy_Chi = sin(xy_Lat)       ! Y_1^0

  xy_Uans = 0.0D0
  xy_Vans = cos(xy_Lat)

  xy_Vorans = 0.0D0
  xy_Divans = -2*sin(xy_Lat)

  call w_StreamPotential2Vector( w_xy(xy_Psi), w_xy(xy_Chi), xy_U, xy_V )

  call AssertEqual(&
    message='Chi = Y_1^0 test of w_StreamPotential2Vector(U)',   &
    answer = xy_Uans,                                             &
    check  = xy_U,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Chi = Y_1^0 test of w_StreamPotential2Vector(V)',   &
    answer = xy_Vans,                                             &
    check  = xy_V,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  call w_Vector2VorDiv( xy_U, xy_V, w_Vor, w_Div )

  call AssertEqual(&
    message='Chi = Y_1^0 test of w_Vector2VorDiv(Vor)',          &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Chi = Y_1^0 test of w_Vector2VorDiv(Div)',          &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_UCosLat = xy_U * cos(xy_Lat) ; xy_VCosLat = xy_V * cos(xy_Lat)

  call w_VectorCosLat2VorDiv( xy_UCosLat, xy_VCosLat, w_Vor, w_Div )

  call AssertEqual(&
    message='Chi = Y_1^0 test of w_VectorCosLat2VorDiv(Vor)',    &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Chi = Y_1^0 test of w_VectorCosLat2VorDiv(Div)',    &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  !---- Y_2^2 のテスト ----
  xy_Psi = cos(xy_Lat)**2*cos(2*xy_Lon)       ! Y_2^2
  xy_Chi = 0.0D0

  xy_Uans = 2*sin(xy_Lat)*cos(xy_Lat)*cos(2*xy_Lon)
  xy_Vans = -2*cos(xy_Lat)*sin(2*xy_Lon)

  xy_Vorans = -6*cos(xy_Lat)**2*cos(2*xy_Lon)      ! Y_2^2
  xy_Divans = 0.0D0

  call w_StreamPotential2Vector( w_xy(xy_Psi), w_xy(xy_Chi), xy_U, xy_V )

  call AssertEqual(&
    message='Psi = Y_2^2 test of w_StreamPotential2Vector(U)',    &
    answer = xy_Uans,                                             &
    check  = xy_U,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_2^2 test of w_StreamPotential2Vector(V)',    &
    answer = xy_Vans,                                             &
    check  = xy_V,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  call w_Vector2VorDiv( xy_U, xy_V, w_Vor, w_Div )

  call AssertEqual(&
    message='Psi = Y_2^2 test of w_Vector2VorDiv(Vor)',           &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_2^2 test of w_Vector2VorDiv(Div)',           &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_UCosLat = xy_U * cos(xy_Lat) ; xy_VCosLat = xy_V * cos(xy_Lat)

  call w_VectorCosLat2VorDiv( xy_UCosLat, xy_VCosLat, w_Vor, w_Div )

  call AssertEqual(&
    message='Psi = Y_2^2 test of w_VectorCosLat2VorDiv(Vor)',     &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_2^2 test of w_VectorCosLat2VorDiv(Div)',     &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_Psi = 0.0D0
  xy_Chi = cos(xy_Lat)**2*cos(2*xy_Lon)       ! Y_2^2

  xy_Uans = - 2*cos(xy_Lat)*sin(2*xy_Lon)
  xy_Vans = - 2*sin(xy_Lat)*cos(xy_Lat)*cos(2*xy_Lon)

  xy_Vorans = 0.0D0
  xy_Divans = -6*cos(xy_Lat)**2*cos(2*xy_Lon)        ! Y_2^2

  call w_StreamPotential2Vector( w_xy(xy_Psi), w_xy(xy_Chi), xy_U, xy_V )

  call AssertEqual(&
    message='Chi = Y_2^2 test of w_StreamPotential2Vector(U)',    &
    answer = xy_Uans,                                             &
    check  = xy_U,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Chi = Y_2^2 test of w_StreamPotential2Vector(V)',    &
    answer = xy_Vans,                                             &
    check  = xy_V,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  call w_Vector2VorDiv( xy_U, xy_V, w_Vor, w_Div )

  call AssertEqual(&
    message='Chi = Y_2^2 test of w_Vector2VorDiv(Vor)',           &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Chi = Y_2^2 test of w_Vector2VorDiv(Div)',           &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_UCosLat = xy_U * cos(xy_Lat) ; xy_VCosLat = xy_V * cos(xy_Lat)

  call w_VectorCosLat2VorDiv( xy_UCosLat, xy_VCosLat, w_Vor, w_Div )

  call AssertEqual(&
    message='Chi = Y_2^2 test of w_VectorCosLat2VorDiv(Vor)',     &
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Chi = Y_2^2 test of w_VectorCosLat2VorDiv(Div)',     &
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  !---- Y_4^2, Y_4^-4 のテスト ----
  xy_Psi = cos(xy_Lat)**2*(7*sin(xy_Lat)**2-1)*cos(2*xy_Lon)  ! Y_4^2
  xy_Chi = cos(xy_Lat)**4*sin(4*xy_Lon)                       ! Y_4^-4

  xy_Uans = -( -2*cos(xy_Lat)*sin(xy_Lat)*(7*sin(xy_Lat)**2-1) &
               + cos(xy_Lat)**2*14*sin(xy_Lat)*cos(xy_Lat) )*cos(2*xy_Lon) &
            + cos(xy_Lat)**3*4*cos(4*xy_Lon)
  xy_Vans =   cos(xy_Lat)*(7*sin(xy_Lat)**2-1)*(-2)*sin(2*xy_Lon) &
            + (-4)*cos(xy_Lat)**3*sin(xy_Lat)*sin(4*xy_Lon) 
  xy_Vorans = -20*cos(xy_Lat)**2*(7*sin(xy_Lat)**2-1)*cos(2*xy_Lon)
  xy_Divans = -20*cos(xy_Lat)**4*sin(4*xy_Lon)

  call w_StreamPotential2Vector( w_xy(xy_Psi), w_xy(xy_Chi), xy_U, xy_V )

  call AssertEqual(&
    message='Psi = Y_4^-4, Chi = Y_4^2 test of w_StreamPotential2Vector(U)',     &
    answer = xy_Uans,                                             &
    check  = xy_U,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_4^-4, Chi = Y_4^2 test of w_StreamPotential2Vector(V)',     &
    answer = xy_Vans,                                             &
    check  = xy_V,                                                &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  call w_Vector2VorDiv( xy_U, xy_V, w_Vor, w_Div )

  call AssertEqual(&
    message='Psi = Y_4^-4, Chi=Y_4^2 test of w_Vector2VorDiv(Vor)',&
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_4^-4, Chi=Y_4^2 test of w_Vector2VorDiv(Div)',&
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )

  xy_UCosLat = xy_U * cos(xy_Lat) ; xy_VCosLat = xy_V * cos(xy_Lat)

  call w_VectorCosLat2VorDiv( xy_UCosLat, xy_VCosLat, w_Vor, w_Div )

  call AssertEqual(&
    message='Psi = Y_4^-4, Chi=Y_4^2 test of w_VectorCosLat2VorDiv(Vor)',&
    answer = xy_Vorans,                                           &
    check  = xy_w(w_Vor),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )
  call AssertEqual(&
    message='Psi = Y_4^-4, Chi=Y_4^2 test of w_VectorCosLat2VorDiv(Div)',&
    answer = xy_Divans,                                           &
    check  = xy_w(w_Div),                                         &
    significant_digits = check_digits, ignore_digits = ignore     &
    )


  call MessageNotify('M','w_module_base_vector_mint_test', &
                         'w_module_base vector subroutine tests succeeded!')

end program w_module_base_vector_mint_test
