!---------------------------------------------------------------------- ! Copyright (c) 2019 Shin-ichi Takehiro. All rights reserved. !---------------------------------------------------------------------- ! !表題 ua_deriv_mpi_module_sypack テストプログラム :: 微分関数のテスト ! !履歴 2019/04/30 竹広真一 ! program ua_deriv_mpi_module_sypack_test use dc_message, only : MessageNotify use dc_test, only : AssertEqual use ua_mpi_module_sypack !!$ use ua_base_mpi_module_sypack, only : ua_base_mpi_Initial, kp_kl, & !!$ xvb_ua, ua_xvb, jc, kc, xv_Lon, xv_Lat !!$ use ua_deriv_mpi_module_sypack use mpi implicit none integer, parameter :: im=64, jm=32, nm=21, km=2, npv=2 real(8), allocatable :: xvb_data1(:,:,:) ! 元の関数 real(8), allocatable :: xvb_data2(:,:,:) ! 元の関数 real(8), allocatable :: xvb_ddata(:,:,:) ! 微分の正解 real(8), allocatable :: xv_mu(:,:) ! μ=sinφ ! 判定誤差設定 integer, parameter :: check_digits = 10 integer, parameter :: ignore = -11 integer :: ip, np, ierr !---------------- MPI スタート --------------------- call MPI_INIT(IERR) call MPI_COMM_RANK(MPI_COMM_WORLD,IP,IERR) call MPI_COMM_SIZE(MPI_COMM_WORLD,NP,IERR) call MessageNotify('M','ua_deriv_mpi_module_sypack_test', & 'ua_deriv_mpi_module_sypack function tests') !!$ call ua_base_mpi_Initial( nm, im, jm, km, npv=npv ) !!$ call ua_deriv_mpi_Initial call ua_mpi_Initial( nm, im, jm, km, npv=npv ) allocate(xvb_data1(0:im-1,jc,kc)) allocate(xvb_data2(0:im-1,jc,kc)) allocate(xvb_ddata(0:im-1,jc,kc)) allocate(xv_mu(0:im-1,jc)) !============== 微分計算 (lon,lat座標系用) のテスト ============== !---- Y_1^-1 Y_2^1 のテスト ---- if ( kp_kl(1) /= 0 ) & xvb_data1(:,:,kp_kl(1)) = -cos(xv_Lat)*sin(xv_Lon) ! Y_1^{-1} if ( kp_kl(2) /= 0 ) & xvb_data1(:,:,kp_kl(2)) = sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon) ! Y_2^1 if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 2*cos(xv_Lat)*sin(xv_Lon) ! ua_Lapla_ua if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -6*sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon) ! Y_2^1 call AssertEqual(& message='ua_Lapla_ua with Y_1^-1 and Y^2_1', & answer = xvb_ua(ua_Lapla_ua(ua_xvb(xvb_data1))), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 1.0/2.0*cos(xv_Lat)*sin(xv_Lon) ! ua_LaplaInv_ua if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -1.0D0/6.0*sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon) call AssertEqual(& message='ua_LaplaInv_ua with Y_1^-1 and Y^2_1', & answer = xvb_ua(ua_LaplaInv_ua(ua_xvb(xvb_data1))), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = -cos(xv_Lat)*cos(xv_Lon) ! ua_DLon_ua if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -sin(xv_Lat)*cos(xv_Lat) * sin(xv_Lon) call AssertEqual(& message='ua_DLon_ua with Y_1^-1 and Y^2_1', & answer = xvb_ua(ua_DLon_ua(ua_xvb(xvb_data1))), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = -cos(xv_Lon) ! xvb_GradLon_ua if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -sin(xv_Lat) * sin(xv_Lon) call AssertEqual(& message='xvb_GradLon_ua with Y_1^-1 and Y^2_1', & answer = xvb_GradLon_ua(ua_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = sin(xv_Lat)*sin(xv_Lon) ! xvb_GradLat_ua if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = cos(2*xv_Lat) * cos(xv_Lon) call AssertEqual(& message='xvb_GradLat_ua with Y_1^-1 and Y^2_1', & answer = xvb_GradLat_ua(ua_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) !---- Y_1^-1 cosφ, Y_2^1 cosφ のテスト ---- if ( kp_kl(1) /= 0 ) & xvb_data1(:,:,kp_kl(1)) = -cos(xv_Lat)**2*sin(xv_Lon) ! Y_1^-1 cosφ if ( kp_kl(2) /= 0 ) & xvb_data1(:,:,kp_kl(2)) = sin(xv_Lat)*cos(xv_Lat)**2 * cos(xv_Lon) ! Y_2^1 cosφ if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = -cos(xv_Lat)*cos(xv_Lon) ! ua_DivLon_xv if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -sin(xv_Lat)*cos(xv_Lat)*sin(xv_Lon) ! ua_DivLon_xv call AssertEqual(& message='ua_DivLon_xvb with Y_1^-1 cosφ and Y^2_1 cosφ', & answer = xvb_ua(ua_DivLon_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 3*sin(xv_Lat)*cos(xv_Lat)*sin(xv_Lon) !ua_DivLat_ua if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = cos(xv_Lat)*(1-4*sin(xv_Lat)**2)*cos(xv_Lon) !ua_DivLat_ua call AssertEqual(& message='ua_DivLat_xvb with Y_1^-1 cosφ and Y^2_1 cosφ', & answer = xvb_ua(ua_DivLat_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) !---- Jacobian のテスト ---- if ( kp_kl(1) /= 0 ) & xvb_data1(:,:,kp_kl(1)) = -cos(xv_Lat)*sin(xv_Lon) ! Y_1^{-1} if ( kp_kl(2) /= 0 ) & xvb_data1(:,:,kp_kl(2)) = -cos(xv_Lat)*sin(xv_Lon) ! Y_1^{-1} if ( kp_kl(1) /= 0 ) & xvb_data2(:,:,kp_kl(1)) = -cos(xv_Lat)*sin(xv_Lon) ! Y_1^{-1} if ( kp_kl(2) /= 0 ) & xvb_data2(:,:,kp_kl(2)) = sin(xv_Lat)*cos(xv_Lat) * cos(xv_Lon) ! Y_2^1 if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 0.0 if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = sin(xv_Lat)**2 - cos(xv_Lat)**2*cos(xv_Lon)**2 call AssertEqual(& message='ua_JacobianMPI_ua_ua Y_1^-1 and Y^2_1', & answer = xvb_ua(ua_JacobianMPI_ua_ua(ua_xvb(xvb_data1),ua_xvb(xvb_data2))), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) !============== 微分計算 (λ,μ座標系用) のテスト ============== xv_mu = sin(xv_Lat) !----- Y_2^0, Y_1^1 のテスト ----- if ( kp_kl(1) /= 0 ) & xvb_data1(:,:,kp_kl(1)) = 3*xv_mu**2-1 ! Y_2^0 if ( kp_kl(2) /= 0 ) & xvb_data1(:,:,kp_kl(2)) = sqrt(1-xv_mu**2)*cos(xv_Lon) ! Y_1^1 if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 0.0 if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -sqrt(1-xv_mu**2)*sin(xv_Lon) call AssertEqual(& message='xvb_GradLambda_ua with Y_2^0 and Y^1_1', & answer = xvb_GradLambda_ua(ua_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 6*xv_mu*(1-xv_mu**2) if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -xv_mu*sqrt(1-xv_mu**2)*cos(xv_Lon) call AssertEqual(& message='xvb_GradMu_ua with Y_2^0 and Y^1_1', & answer = xvb_GradMu_ua(ua_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) !----- Y_2^0(1-μ^2), Y_1^1 (1-μ^2) のテスト ----- if ( kp_kl(1) /= 0 ) & xvb_data1(:,:,kp_kl(1)) = (3*xv_mu**2-1)*(1-xv_mu**2) ! Y_2^0 (1-μ^2) if ( kp_kl(2) /= 0 ) & xvb_data1(:,:,kp_kl(2)) = (1-xv_mu**2)**(3.0/2)*cos(xv_Lon) ! Y_1^1 (1-μ^2) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = 0.0 if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -(1-xv_mu**2)**(1.0d0/2)*sin(xv_Lon) call AssertEqual(& message='ua_DivLambda_xvb with Y_2^0 (1-μ^2) and Y^1_1 √(1-μ^2)',& answer = xvb_ua(ua_DivLambda_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) if ( kp_kl(1) /= 0 ) & xvb_ddata(:,:,kp_kl(1)) = (2-3*xv_mu**2)*4*xv_mu if ( kp_kl(2) /= 0 ) & xvb_ddata(:,:,kp_kl(2)) = -3.0D0*xv_mu*(1-xv_mu**2)**(1.0D0/2)*cos(xv_Lon) call AssertEqual(& message='ua_DivMu_xvb with Y_2^0 (1-μ^2) and Y^1_1 √(1-μ^2)',& answer = xvb_ua(ua_DivMu_xvb(xvb_data1)), & check = xvb_ddata, & significant_digits = check_digits, ignore_digits = ignore & ) call MessageNotify('M','ua_deriv_mpi_module_sypack_test', & 'ua_deriv_mpi_module_sypack function tests succeeded!') !------ MPIの終了 ------ call MPI_FINALIZE(IERR) end program ua_deriv_mpi_module_sypack_test