!---------------------------------------------------------------------- ! Copyright (c) 2008 Shin-ichi Takehiro. All rights reserved. !---------------------------------------------------------------------- ! !表題 eee_mpi_module テストプログラム ! !履歴 2008/05/21 竹広真一 ! program eee_mpi_test_intavr use dc_message, only : MessageNotify use eee_mpi_module use mpi implicit none !include 'mpif.h' !---- 空間解像度設定 ---- integer, parameter :: im=32, jm=16, km=64 ! 格子点の設定(X,Y,Z) integer, parameter :: lm=10, mm=5, nm=21 ! 切断波数の設定(X,Y,Z) !---- 変数 ---- real(8), allocatable :: zxv_Data(:,:,:) ! 格子データ real(8), allocatable :: xv_Data(:,:) ! 格子データ real(8), allocatable :: zx_Data(:,:) ! 格子データ real(8), allocatable :: zv_Data(:,:) ! 格子データ real(8), allocatable :: x_Data(:) ! 格子データ real(8), allocatable :: v_Data(:) ! 格子データ real(8), allocatable :: z_Data(:) ! 格子データ !---- 座標変数など ---- real(8), parameter :: pi=3.1415926535897932385D0 real(8), parameter :: eps = 1.0d-10 ! 判定誤差 integer :: np, ip, ierr call MessageNotify('M','eee_mpi_test_intavr', & 'eee_mpi_module integral/average function tests') !---------------- MPI スタート --------------------- call MPI_INIT(IERR) call MPI_COMM_RANK(MPI_COMM_WORLD,IP,IERR) call MPI_COMM_SIZE(MPI_COMM_WORLD,NP,IERR) !---------------- 座標値の設定 --------------------- call eee_mpi_initial(im,jm,km,lm,mm,nm) !---------------- 変数の割付け --------------------- allocate(zxv_Data(0:km-1,0:im-1,js(ip):je(ip))) allocate(xv_Data(0:im-1,js(ip):je(ip))) allocate(zv_Data(0:km-1,js(ip):je(ip))) allocate(zx_Data(0:km-1,0:im-1)) allocate(z_Data(0:km-1)) allocate(x_Data(0:im-1)) allocate(v_Data(js(ip):je(ip))) !-------------------------- 3 次元データ ------------------------------ zxv_Data = sin(zxv_X) * cos(zxv_Y) * cos(zxv_Z) write(6,*) 'f = sin(X)*cos(Y)*cos(Z)' call check0d(IntZXV_zxv(zxv_Data)-0.0, eps, 'IntZXV_zxv') call check0d(AvrZXV_zxv(zxv_Data)-0.0, eps, 'AvrZXV_zxv') call check1d(x_IntZV_zxv(zxv_Data) - 0.0, eps, 'x_IntZV_zxv') call check1d(x_AvrZV_zxv(zxv_Data) - 0.0, eps, 'x_AvrZV_zxv') call check1d(v_IntZX_zxv(zxv_Data) - 0.0, eps, 'v_IntZX_zxv') call check1d(v_AvrZX_zxv(zxv_Data) - 0.0, eps, 'v_AvrZX_zxv') call check1d(z_IntXV_zxv(zxv_Data) - 0.0, eps, 'x_IntXV_zxv') call check1d(z_AvrXV_zxv(zxv_Data) - 0.0, eps, 'x_AvrXV_zxv') call check2d(xv_IntZ_zxv(zxv_Data) - 0.0, eps, 'xv_IntZ_zxv') call check2d(xv_AvrZ_zxv(zxv_Data) - 0.0, eps, 'xv_AvrZ_zxv') call check2d(zx_IntV_zxv(zxv_Data) - 0.0, eps, 'zx_IntV_zxv') call check2d(zx_AvrV_zxv(zxv_Data) - 0.0, eps, 'zx_AvrV_zxv') call check2d(zv_IntX_zxv(zxv_Data) - 0.0, eps, 'zv_IntX_zxv') call check2d(zv_AvrX_zxv(zxv_Data) - 0.0, eps, 'zv_AvrX_zxv') zxv_Data = (1 - sin(zxv_X)) * cos(zxv_Y)**2 * cos(zxv_Z)**2 write(6,*) write(6,*) 'f = (1-sin(X)) * cos(Y)**2 * cos(Z)**2' call check0d(IntZXV_zxv(zxv_Data)-2.0*pi**3, eps, 'IntZXV_zxv') call check0d(AvrZXV_zxv(zxv_Data)-0.25, eps, 'AvrZXV_zxv') call check1d(x_IntZV_zxv(zxv_Data) - (1 - sin(x_X))*pi**2, eps, 'x_IntZV_zxv') call check1d(x_AvrZV_zxv(zxv_Data) - (1 - sin(x_X))/4, eps, 'x_AvrZV_zxv') call check1d(v_IntZX_zxv(zxv_Data) - cos(v_Y)**2 * 2*pi**2, eps, 'v_IntZX_zxv') call check1d(v_AvrZX_zxv(zxv_Data) - cos(v_Y)**2 /2 , eps, 'v_AvrZX_zxv') call check1d(z_IntXV_zxv(zxv_Data) - cos(z_Z)**2 *2*pi**2 , eps, 'x_IntXV_zxv') call check1d(z_AvrXV_zxv(zxv_Data) - cos(z_Z)**2 /2, eps, 'x_AvrXV_zxv') call check2d(xv_IntZ_zxv(zxv_Data) & - (1 - sin(zxv_X(1,:,:)))*cos(zxv_Y(1,:,:))**2*pi, eps, 'xv_IntZ_zxv') call check2d(xv_AvrZ_zxv(zxv_Data) & - (1 - sin(zxv_X(1,:,:)))*cos(zxv_Y(1,:,:))**2/2,eps, 'xv_AvrZ_zxv') call check2d(zv_IntX_zxv(zxv_Data) & - cos(zxv_Y(:,1,:))**2 * cos(zxv_Z(:,1,:))**2*2*pi, eps, 'zv_IntX_zxv') call check2d(zv_AvrX_zxv(zxv_Data) & - cos(zxv_Y(:,1,:))**2 * cos(zxv_Z(:,1,:))**2, eps, 'zv_AvrX_zxv') call check2d(zx_IntV_zxv(zxv_Data) & - (1 - sin(zxv_X(:,:,js(ip))))*cos(zxv_Z(:,:,js(ip)))**2*pi, eps, 'zx_IntV_zxv') call check2d(zx_AvrV_zxv(zxv_Data) & - (1 - sin(zxv_X(:,:,js(ip))))*cos(zxv_Z(:,:,js(ip)))**2/2, eps,'zx_AvrV_zxv') !-------------------------- 2 次元データ ------------------------------ xv_Data = sin(zxv_X(1,:,:)) * cos(zxv_Y(1,:,:)) write(6,*) 'f = sin(X)*cos(Y)' call check0d(IntXV_xv(xv_Data)-0.0, eps, 'IntXV_xv') call check0d(AvrXV_xv(xv_Data)-0.0, eps, 'AvrXV_xv') call check1d(v_IntX_xv(xv_Data) - 0.0, eps, 'v_IntX_xv') call check1d(v_AvrX_xv(xv_Data) - 0.0, eps, 'v_AvrX_xv') call check1d(x_IntV_xv(xv_Data) - 0.0, eps, 'x_IntV_xv') call check1d(x_AvrV_xv(xv_Data) - 0.0, eps, 'x_AvrV_xv') xv_Data = (1 - sin(zxv_X(1,:,:))) * cos(zxv_Y(1,:,:))**2 write(6,*) write(6,*) 'f = (1-sin(X)) * cos(Y)**2' call check0d(IntXV_xv(xv_Data) - 2.0*pi**2, eps, 'IntXV_xv') call check0d(AvrXV_xv(xv_Data) - 0.5, eps, 'AvrXV_xv') call check1d(v_IntX_xv(xv_Data) - 2*cos(v_Y)**2*pi, eps, 'v_IntX_xv') call check1d(v_AvrX_xv(xv_Data) - cos(v_Y)**2, eps, 'v_AvrX_xv') call check1d(x_IntV_xv(xv_Data) - (1 - sin(x_X))*pi, eps, 'x_IntV_xv') call check1d(x_AvrV_xv(xv_Data) - 0.5*(1 - sin(x_X)), eps, 'x_AvrV_xv') zx_Data = sin(zxv_X(:,:,js(ip))) * cos(zxv_Z(:,:,js(ip))) write(6,*) 'f = sin(X)*cos(Y)' call check0d(IntZX_zx(zx_Data)-0.0, eps, 'IntZX_zx') call check0d(AvrZX_zx(zx_Data)-0.0, eps, 'AvrZX_zx') call check1d(z_IntX_zx(zx_Data) - 0.0, eps, 'z_IntX_zx') call check1d(z_AvrX_zx(zx_Data) - 0.0, eps, 'z_AvrX_zx') call check1d(x_IntZ_zx(zx_Data) - 0.0, eps, 'x_IntZ_zx') call check1d(x_AvrZ_zx(zx_Data) - 0.0, eps, 'x_AvrZ_zx') zx_Data = (1 - sin(zxv_X(:,:,js(ip)))) * cos(zxv_Z(:,:,js(ip)))**2 write(6,*) write(6,*) 'f = (1-sin(X)) * cos(Z)**2' call check0d(IntZX_zx(zx_Data) - 2.0*pi**2, eps, 'IntZX_zx') call check0d(AvrZX_zx(zx_Data) - 0.5, eps, 'AvrZX_zx') call check1d(z_IntX_zx(zx_Data) - 2*cos(z_Z)**2*pi, eps, 'z_IntX_zx') call check1d(z_AvrX_zx(zx_Data) - cos(z_Z)**2, eps, 'z_AvrX_zx') call check1d(x_IntZ_zx(zx_Data) - (1 - sin(x_X))*pi, eps, 'x_IntZ_zx') call check1d(x_AvrZ_zx(zx_Data) - 0.5*(1 - sin(x_X)), eps, 'x_AvrZ_zx') zv_Data = sin(zxv_Y(:,1,:)) * cos(zxv_Z(:,1,:)) write(6,*) 'f = sin(Y)*cos(Y)' call check0d(IntZV_zv(zv_Data)-0.0, eps, 'IntZV_zv') call check0d(AvrZV_zv(zv_Data)-0.0, eps, 'AvrZV_zv') call check1d(z_IntV_zv(zv_Data) - 0.0, eps, 'z_IntV_zv') call check1d(z_AvrV_zv(zv_Data) - 0.0, eps, 'z_AvrV_zv') call check1d(v_IntZ_zv(zv_Data) - 0.0, eps, 'v_IntZ_zv') call check1d(v_AvrZ_zv(zv_Data) - 0.0, eps, 'v_AvrZ_zv') zv_Data = (1 - sin(zxv_Y(:,1,:))) * cos(zxv_Z(:,1,:))**2 write(6,*) write(6,*) 'f = (1-sin(Y)) * cos(Z)**2' call check0d(IntZV_zv(zv_Data) - 2.0*pi**2, eps, 'IntZV_zv') call check0d(AvrZV_zv(zv_Data) - 0.5, eps, 'AvrZV_zv') call check1d(z_IntV_zv(zv_Data) - 2*cos(z_Z)**2*pi, eps, 'z_IntV_zv') call check1d(z_AvrV_zv(zv_Data) - cos(z_Z)**2, eps, 'z_AvrV_zv') call check1d(v_IntZ_zv(zv_Data) - (1 - sin(v_Y))*pi, eps, 'v_IntZ_zv') call check1d(v_AvrZ_zv(zv_Data) - 0.5*(1 - sin(v_Y)), eps, 'v_AvrZ_zv') !-------------------------- 1 次元データ ------------------------------ x_Data = 1 - sin(2*x_X) write(6,*) write(6,*) 'f = 1-sin(2*X)' call check0d(IntX_x(x_Data) - 2.0d0*pi, eps, 'IntX_x') call check0d(AvrX_x(x_Data) - 1.0d0, eps, 'AvrX_x') v_Data = cos(v_Y)**2 write(6,*) write(6,*) 'f = cos(Y)**2' call check0d(IntV_v(v_Data) - 1.0*pi, eps, 'IntV_v') call check0d(AvrV_v(v_Data) - 0.5, eps, 'AvrV_v') z_Data = 1 - sin(2*z_Z) write(6,*) write(6,*) 'f = 1-sin(2*Z)' call check0d(IntZ_z(z_Data) - 2.0d0*pi, eps, 'IntZ_z') call check0d(AvrZ_z(z_Data) - 1.0d0, eps, 'AvrZ_z') call MessageNotify('M','eee_mpi_test_intavr', & 'eee_mpi_module integral/average function tests succeeded!') call MPI_FINALIZE(IERR) stop contains subroutine check2d(var,eps,funcname) ! 絶対値が eps 以上の var の要素を出力 real(8) :: var(:,:) ! 判定する配列 real(8) :: eps ! 誤差 character(len=*), optional :: funcname character(len=3) ::cip integer i,j if ( present(funcname) )then write(cip,'(I3)') IP write(6,*) ' Checking ', funcname, ' for IP='//trim(adjustl(cip)) endif do j=1,size(var,2) do i=1,size(var,1) if (abs(var(i,j)) .gt. eps ) then write(6,*) ' Value larger than EPS : i= ', i,' j= ',j, var(i,j) call MessageNotify('E','eee_test_intavr', & 'Integral/average error too large') endif enddo enddo end subroutine check2d subroutine check1d(var,eps,funcname) ! 絶対値が eps 以上の var の要素を出力 real(8) :: var(:) ! 判定する配列 real(8) :: eps ! 誤差 character(len=*), optional :: funcname character(len=3) ::cip integer i if ( present(funcname) )then write(cip,'(I3)') IP write(6,*) ' Checking ', funcname, ' for IP='//trim(adjustl(cip)) endif do i=1,size(var) if (abs(var(i)) .gt. eps ) then write(6,*) ' Value larger than EPS : i= ', i, var(i) call MessageNotify('E','eee_test_intavr', & 'Integral/average error too large') endif enddo end subroutine check1d subroutine check0d(var,eps,funcname) ! var の絶対値が eps 以上だと出力 real(8) :: var ! 判定する配列 real(8) :: eps ! 誤差 real(8) :: vartmp(1) ! character(len=*), optional :: funcname vartmp(1) = var if ( present(funcname) ) then call check1d(vartmp,eps,funcname) else call check1d(vartmp,eps) endif end subroutine check0d end program eee_mpi_test_intavr