!---------------------------------------------------------------------- ! Copyright (c) 2016 Shin-ichi Takehiro. All rights reserved. !---------------------------------------------------------------------- ! !表題 ut_mpi_module_svpack テストプログラム :: スペクトル関数のテスト ! !履歴 2016/03/30 竹広真一 ! program ut_mpi_module_svpack_spectrum_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_Torvel(:,:,:) real(8), allocatable :: ut_Torvel(:,:) real(8), allocatable :: xyh_Polvel(:,:,:) real(8), allocatable :: ut_Polvel(:,:) real(8), allocatable :: xyh_Vlon(:,:,:) real(8), allocatable :: xyh_Vlat(:,:,:) real(8), allocatable :: xyh_Vrad(:,:,:) real(8), dimension(0:nm,-nm:nm,0:km) :: nmz_EkTor real(8), dimension(0:nm,-nm:nm,0:km) :: nmz_EkPol real(8), dimension(0:nm,0:km) :: nz_EkTor real(8), dimension(0:nm,0:km) :: nz_EkPol real(8), dimension(0:km) :: z_Data ! 判定誤差設定 integer, parameter :: check_digits = 10 integer, parameter :: ignore = -11 real(8) :: PI integer :: n,m integer :: iproc, np, ierr PI = 4.0D0*Atan(1.0D0) !---------------- MPI スタート --------------------- call MPI_INIT(IERR) call MPI_COMM_RANK(MPI_COMM_WORLD,IPROC,IERR) call MPI_COMM_SIZE(MPI_COMM_WORLD,NP,IERR) call ut_mpi_Initial(im,jm,km,nm,lm,ri,ro) call MessageNotify('M','ut_mpi_module_svpack_spectrum_test',& 'ut_mpi_module_svpack spectrum functions tests.') allocate(xyh_Torvel(0:im-1,jm,kc)) allocate(ut_Torvel(nc,0:lm)) allocate(xyh_Polvel(0:im-1,jm,kc)) allocate(ut_Polvel(nc,0:lm)) allocate(xyh_Vlon(0:im-1,jm,kc)) allocate(xyh_Vlat(0:im-1,jm,kc)) allocate(xyh_Vrad(0:im-1,jm,kc)) ! (2Y_1^0+3Y_2^1) * r xyh_Torvel =( 2*sqrt(3.0D0)*sin(xyh_Lat) & + 3*sqrt(5.0D0/6)*3.0*sin(xyh_Lat)*cos(xyh_Lat) * cos(xyh_Lon) )& * xyh_Rad ! (3*Y_1^0 + 2Y_1^1 + 2Y_2^0 +Y_2^-2) * r^2 xyh_Polvel =( 3*sqrt(3.0D0)*sin(xyh_Lat) & + 2* sqrt(3.0D0/2)*cos(xyh_Lat)*cos(xyh_Lon) & + 2*sqrt(5.0D0)*(3.0/2*sin(xyh_Lat)**2-1/2.0) & - sqrt(5.0D0/24)*3.0*cos(xyh_Lat)**2*sin(2*xyh_Lon) ) & * xyh_Rad**2 ut_Torvel = ut_xyh(xyh_Torvel) ut_Polvel = ut_xyh(xyh_Polvel) !--------- Energy spectrum --------- nmz_EkTor=ut_VMiss nmz_EkPol=ut_VMiss do n=0,nm do m=-n,n nmz_EkTor(n,m,:) = 0.0D0 nmz_EkPol(n,m,:) = 0.0D0 enddo enddo nmz_EkTor(1,0,:) = 4.0D0 * z_Rad**2 * (4*pi)*z_Rad**2 nmz_EkTor(2,1,:) = 2*3*(3.0D0/2)**2/2 * z_Rad**2 * (4*pi)*z_Rad**2 nmz_EkTor(2,-1,:)= nmz_EkTor(2,1,:) nmz_EkPol(1,0,:) = 9.0D0 * (9*z_Rad**4 + 2.0D0*z_Rad**4) * (4*pi) nmz_EkPol(1,1,:) = 1*2*(2/2)**2/2 & * (9*z_Rad**4 + 1*2* z_Rad**4) * (4*pi) nmz_EkPol(1,-1,:)= nmz_EkPol(1,1,:) nmz_EkPol(2,0,:) = 12.0D0 *(9*z_Rad**4 + 6.0D0*z_Rad**4) * (4*pi) nmz_EkPol(2,-2,:)= 2*3*(1.0D0/2)**2/2 & *(9*z_Rad**4 + 2*3* z_Rad**4 ) * (4*pi) nmz_EkPol(2,2,:) = nmz_EkPol(2,-2,:) call AssertEqual(& message='nmz_ToroidalEnergySpectrum_ut', & answer = nmz_ToroidalEnergySpectrum_ut(ut_Torvel), & check = nmz_EkTor, & significant_digits = check_digits, ignore_digits = ignore & ) call AssertEqual(& message='nmz_PoloidalEnergySpectrum_ut', & answer = nmz_PoloidalEnergySpectrum_ut(ut_Polvel), & check = nmz_EkPol, & significant_digits = check_digits, ignore_digits = ignore & ) nz_EkTor = 0.0D0 nz_EkTor(1,:) = 4.0D0 * z_Rad**2 * (4*pi)*z_Rad**2 nz_EkTor(2,:) = 2 * 3 *(3.0D0/2)**2 * z_Rad**2 * (4*pi)*z_Rad**2 nz_EkPol = 0.0D0 nz_EkPol(1,:) = (9.0D0+1*2*(2/2)**2/2*2 )*(9*z_Rad**4+2*z_Rad**4)*(4*pi) nz_EkPol(2,:) = (12.0D0+2*3*(1.0D0/2)**2)*(9*z_Rad**4+2*3*z_Rad**4)*(4*pi) call AssertEqual(& message='nz_ToroidalEnergySpectrum_ut', & answer = nz_ToroidalEnergySpectrum_ut(ut_Torvel), & check = nz_EkTor, & significant_digits = check_digits, ignore_digits = ignore & ) call AssertEqual(& message='nz_PoloidalEnergySpectrum_ut', & answer = nz_PoloidalEnergySpectrum_ut(ut_Polvel), & check = nz_EkPol, & significant_digits = check_digits, ignore_digits = ignore & ) !!$ call ut_Potential2VectorMPI(xyh_Vlon,xyh_Vlat,xyh_Vrad,ut_Torvel,ut_Polvel) !!$ !!$ nmz_EkTor = nmz_ToroidalEnergySpectrum_ut(ut_Torvel) !!$ nmz_EkPol = nmz_PoloidalEnergySpectrum_ut(ut_Polvel) !!$ !!$ z_Data = 0.0D0 !!$ do n=0,nm !!$ do m=-n,n !!$ z_Data = z_Data + nmz_EkTor(n,m,:) + nmz_EkPol(n,m,:) !!$ enddo !!$ enddo !!$ !!$ call AssertEqual(& !!$ message='Total energy by nmz-spectrum', & !!$ answer = IntLonLatRad_xyh(xyh_VLon**2+xyh_VLat**2+xyh_VRad**2)/2.0D0, & !!$ check = IntRad_z(z_Data/z_Rad**2), & !!$ significant_digits = check_digits, ignore_digits = ignore & !!$ ) !!$ !!$ nz_EkTor = nz_ToroidalEnergySpectrum_ut(ut_Torvel) !!$ nz_EkPol = nz_PoloidalEnergySpectrum_ut(ut_Polvel) !!$ !!$ call AssertEqual(& !!$ message='Total energy by nz-spectrum', & !!$ answer = IntLonLatRad_xyh(xyh_VLon**2+xyh_VLat**2+xyh_VRad**2)/2.0D0, & !!$ check = IntRad_z(sum(nz_EkTor,1)/z_Rad**2) & !!$ + IntRad_z(sum(nz_EkPol,1)/z_Rad**2), & !!$ significant_digits = check_digits, ignore_digits = ignore & !!$ ) call MessageNotify('M','ut_mpi_module_svpack_spectrum_test',& 'ut_mpi_module_svpack spectrum functions tests suceeded!') CALL MPI_FINALIZE(IERR) end program ut_mpi_module_svpack_spectrum_test