mirror of
https://github.com/JHUAPL/kaiju.git
synced 2026-01-09 15:17:56 -05:00
305 lines
11 KiB
Plaintext
305 lines
11 KiB
Plaintext
module testVMpiDeep
|
|
use testHelperMpi
|
|
use voltapp_mpi
|
|
use gamCouple_mpi_G2V
|
|
use uservoltic
|
|
use ioH5
|
|
|
|
implicit none
|
|
|
|
type(gamCouplerMpi_gam_T), allocatable :: gamCplMpi
|
|
type(voltAppMpi_T), allocatable :: voltAppMpi
|
|
|
|
contains
|
|
|
|
@before
|
|
subroutine setup(this)
|
|
class (MpiTestMethod), intent(inout) :: this
|
|
character(len=strLen) :: caseFile
|
|
|
|
integer :: ierror
|
|
type(MPI_Comm) :: gamComm, voltComm
|
|
type(XML_Input_T) :: xmlInp
|
|
|
|
call setMpiReal()
|
|
|
|
write(caseFile,'(A,I0,A)') 'cmiD_deep_', this%getNumProcesses()-1, '.xml'
|
|
|
|
if(this%getProcessRank() < (this%getNumProcesses()-1)) then
|
|
allocate(gamCplMpi)
|
|
|
|
! make gamera-only mpi communicator
|
|
call MPI_Comm_Split(getMpiF08Communicator(this), gamId, this%getProcessRank(), gamCplMpi%gOptionsMpi%gamComm, ierror)
|
|
|
|
gamCplMpi%gOptionsCplMpiG%couplingPoolComm = getMpiF08Communicator(this)
|
|
gamCplMpi%gOptions%userInitFunc => initUser
|
|
gamCplMpi%gOptionsMpi%doIO = .false.
|
|
xmlInp = New_XML_Input(trim(caseFile),'Kaiju',.true.)
|
|
call gamCplMpi%InitModel(xmlInp)
|
|
else
|
|
allocate(voltAppMpi)
|
|
|
|
! make gamera-only mpi communicator
|
|
call MPI_Comm_Split(getMpiF08Communicator(this), voltId, this%getProcessRank(), voltComm, ierror)
|
|
|
|
voltAppMpi%vOptions%gamUserInitFunc => initUser
|
|
allocate(voltAppMpi%vOptionsMpi%couplingPoolComm)
|
|
voltAppMpi%vOptionsMpi%couplingPoolComm = getMpiF08Communicator(this)
|
|
call initVoltron_mpi(voltAppMpi, trim(caseFile))
|
|
endif
|
|
|
|
end subroutine setup
|
|
|
|
@after
|
|
subroutine teardown(this)
|
|
class (MpiTestMethod), intent(inout) :: this
|
|
|
|
if(allocated(voltAppMpi)) then
|
|
call endVoltronWaits(voltAppMpi)
|
|
deallocate(voltAppMpi)
|
|
endif
|
|
if(allocated(gamCplMpi)) deallocate(gamCplMpi)
|
|
|
|
end subroutine teardown
|
|
|
|
@test(npes=[2,5,9])
|
|
subroutine testMhd2VoltDeepCopy(this)
|
|
class (MpiTestMethod), intent(inout) :: this
|
|
|
|
integer :: i,j,k,l,m
|
|
real(rp) :: checkVal
|
|
character(len=strLen) :: checkMessage
|
|
real(rp) :: iVal, jVal, kVal, lVal, mVal
|
|
|
|
! coefficients chosen to allow easy interpretation of final values
|
|
iVal = 1
|
|
jVal = 20
|
|
kVal = 400
|
|
lVal = 8000
|
|
mVal = 160000
|
|
|
|
if(allocated(gamCplMpi)) then
|
|
! clear all cells
|
|
gamCplMpi%State%Gas(:,:,:,:,:) = 0
|
|
gamCplMpi%State%Bxyz(:,:,:,:) = 0
|
|
|
|
! set values for gas
|
|
do i=gamCplMpi%Grid%is,gamCplMpi%Grid%ie
|
|
do j=gamCplMpi%Grid%js,gamCplMpi%Grid%je
|
|
do k=gamCplMpi%Grid%ks,gamCplMpi%Grid%ke
|
|
do l=1,NVAR
|
|
do m=0,gamCplMpi%Model%nSpc
|
|
gamCplMpi%State%Gas(i,j,k,l,m) = &
|
|
mVal*m + &
|
|
lVal*l + &
|
|
kVal*(k+gamCplMpi%Grid%ijkShift(KDIR)) + &
|
|
jVal*(j+gamCplMpi%Grid%ijkShift(JDIR)) + &
|
|
iVal*(i+gamCplMpi%Grid%ijkShift(IDIR))
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
|
|
! set values for bxyz
|
|
do i=gamCplMpi%Grid%is,gamCplMpi%Grid%ie
|
|
do j=gamCplMpi%Grid%js,gamCplMpi%Grid%je
|
|
do k=gamCplMpi%Grid%ks,gamCplMpi%Grid%ke
|
|
do l=1,NDIM
|
|
gamCplMpi%State%Bxyz(i,j,k,l) = &
|
|
lVal*l + &
|
|
kVal*(k+gamCplMpi%Grid%ijkShift(KDIR)) + &
|
|
jVal*(j+gamCplMpi%Grid%ijkShift(JDIR)) + &
|
|
iVal*(i+gamCplMpi%Grid%ijkShift(IDIR))
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
|
|
! send the data to voltron
|
|
call sendVoltronCplDataMpi(gamCplMpi)
|
|
else
|
|
! clear all cells
|
|
voltAppMpi%gApp%State%Gas(:,:,:,:,:) = 0
|
|
voltAppMpi%gApp%State%Bxyz(:,:,:,:) = 0
|
|
|
|
SELECT type(cpl=>voltAppMpi%gApp)
|
|
TYPE IS (gamCouplerMpi_volt_T)
|
|
! receive the data from gamera
|
|
call recvGameraCplDataMpi(cpl)
|
|
CLASS DEFAULT
|
|
@assertTrue(.false., "Voltron Allocated non-mpi Gamera coupler for MPI Voltron Coupling Test. Failure")
|
|
ENDSELECT
|
|
|
|
! check gas values
|
|
do i=voltAppMpi%gApp%Grid%is,voltAppMpi%gApp%Grid%ie
|
|
do j=voltAppMpi%gApp%Grid%js,voltAppMpi%gApp%Grid%je
|
|
do k=voltAppMpi%gApp%Grid%ks,voltAppMpi%gApp%Grid%ke
|
|
do l=1,NVAR
|
|
do m=0,voltAppMpi%gApp%Model%nSpc
|
|
checkVal = mVal*m + &
|
|
lVal*l + &
|
|
iVal*i + &
|
|
jVal*j + &
|
|
kVal*k
|
|
write (checkMessage,'(A,I0,A,I0,A,I0,A,I0,A,I0,A)') 'voltron gas wrong at (',i,',',j,',',k,',',l,',',m,')'
|
|
@assertEqual(checkVal,voltAppMpi%gApp%State%Gas(i,j,k,l,m), trim(checkMessage))
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
|
|
! check bxyz values
|
|
do i=voltAppMpi%gApp%Grid%is,voltAppMpi%gApp%Grid%ie
|
|
do j=voltAppMpi%gApp%Grid%js,voltAppMpi%gApp%Grid%je
|
|
do k=voltAppMpi%gApp%Grid%ks,voltAppMpi%gApp%Grid%ke
|
|
do l=1,NDIM
|
|
checkVal = lVal*l + &
|
|
iVal*i + &
|
|
jVal*j + &
|
|
kVal*k
|
|
write (checkMessage,'(A,I0,A,I0,A,I0,A,I0,A)') 'voltron bxyz wrong at (',i,',',j,',',k,',',l,')'
|
|
@assertEqual(checkVal, voltAppMpi%gApp%State%Bxyz(i,j,k,l), trim(checkMessage))
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
|
|
endif
|
|
|
|
end subroutine testMhd2VoltDeepCopy
|
|
|
|
@test(npes=[2,5,9])
|
|
subroutine testVolt2MhdDeepCopy(this)
|
|
class (MpiTestMethod), intent(inout) :: this
|
|
|
|
integer :: i,j,k,l,m
|
|
real(rp) :: checkVal
|
|
character(len=strLen) :: checkMessage
|
|
real(rp) :: iVal, jVal, kVal, lVal
|
|
|
|
! coefficients chosen to allow easy interpretation of final values
|
|
iVal = 1
|
|
jVal = 20
|
|
kVal = 400
|
|
lVal = 8000
|
|
|
|
if(allocated(gamCplMpi)) then
|
|
! clear all cells
|
|
gamCplMpi%Grid%Gas0(:,:,:,:) = 0
|
|
|
|
! receive the data from voltron
|
|
call recvDeepCplDataMpi(gamCplMpi)
|
|
|
|
! check gas0 values
|
|
do i=gamCplMpi%Grid%isg,gamCplMpi%Grid%ieg
|
|
do j=gamCplMpi%Grid%jsg,gamCplMpi%Grid%jeg
|
|
do k=gamCplMpi%Grid%ksg,gamCplMpi%Grid%keg
|
|
do l=1,gamCplMpi%Model%nvSrc
|
|
checkVal = lVal*l + &
|
|
kVal*(k+gamCplMpi%Grid%ijkShift(KDIR)) + &
|
|
jVal*(j+gamCplMpi%Grid%ijkShift(JDIR)) + &
|
|
iVal*(i+gamCplMpi%Grid%ijkShift(IDIR))
|
|
write (checkMessage,'(A,I0,A,I0,A,I0,A,I0,A)') 'gamera gas0 wrong at (',i,',',j,',',k,',',l,')'
|
|
@assertEqual(checkVal,gamCplMpi%Grid%Gas0(i,j,k,l), trim(checkMessage))
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
else
|
|
! clear all cells
|
|
voltAppMpi%gApp%Grid%Gas0(:,:,:,:) = 0
|
|
|
|
! set gas0 values
|
|
do i=voltAppMpi%gApp%Grid%isg,voltAppMpi%gApp%Grid%ieg
|
|
do j=voltAppMpi%gApp%Grid%jsg,voltAppMpi%gApp%Grid%jeg
|
|
do k=voltAppMpi%gApp%Grid%ksg,voltAppMpi%gApp%Grid%keg
|
|
do l=1,voltAppMpi%gApp%Model%nvSrc
|
|
voltAppMpi%gApp%Grid%Gas0(i,j,k,l) = lVal*l + &
|
|
iVal*i + &
|
|
jVal*j + &
|
|
kVal*k
|
|
enddo
|
|
enddo
|
|
enddo
|
|
enddo
|
|
|
|
SELECT type(cpl=>voltAppMpi%gApp)
|
|
TYPE IS (gamCouplerMpi_volt_T)
|
|
! send the data to gamera
|
|
call sendDeepCplDataMpi(cpl)
|
|
CLASS DEFAULT
|
|
@assertTrue(.false., "Voltron Allocated non-mpi Gamera coupler for MPI Voltron Coupling Test. Failure")
|
|
ENDSELECT
|
|
|
|
endif
|
|
end subroutine testVolt2MhdDeepCopy
|
|
|
|
@test(npes=[2,5,9])
|
|
subroutine VMPIdipoleTest(this)
|
|
class (MpiTestMethod), intent(inout) :: this
|
|
|
|
! testing that a dipole field from gamera results in 0 FAC in remix
|
|
procedure(VectorField_T), pointer :: Axyz
|
|
real(rp) testValue
|
|
|
|
if(allocated(gamCplMpi)) then
|
|
! create a dipole field in Gamera. This code is copied from prob.F90
|
|
Axyz => VP_Dipole
|
|
call VectorPot2Flux(gamCplMpi%Model,gamCplMpi%Grid,gamCplMpi%State,Axyz)
|
|
|
|
call bFlux2Fld(gamCplMpi%Model, gamCplMpi%Grid, gamCplMpi%State%magFlux, gamCplMpi%State%Bxyz)
|
|
|
|
! send the data to voltron
|
|
call sendVoltronCplDataMpi(gamCplMpi)
|
|
|
|
else
|
|
SELECT type(cpl=>voltAppMpi%gApp)
|
|
TYPE IS (gamCouplerMpi_volt_T)
|
|
! receive the data from gamera
|
|
call recvGameraCplDataMpi(cpl)
|
|
CLASS DEFAULT
|
|
@assertTrue(.false., "Voltron Allocated non-mpi Gamera coupler for MPI Voltron Coupling Test. Failure")
|
|
ENDSELECT
|
|
|
|
call convertGameraToRemix(voltAppMpi%mhd2mix, voltAppMpi%gApp, voltAppMpi%remixApp)
|
|
|
|
call mapGameraToRemix(voltAppMpi%mhd2mix, voltAppMpi%remixApp)
|
|
|
|
! verify that remix input FAC is all zeroes
|
|
testValue = sum(abs(voltAppMpi%remixApp%ion(NORTH)%St%Vars(:,:,FAC))) + sum(abs(voltAppMpi%remixApp%ion(SOUTH)%St%Vars(:,:,FAC)))
|
|
|
|
@assertLessThanOrEqual(testValue, 1e-1_rp, 'Remix did not get all 0 FAC from an input dipole field')
|
|
endif
|
|
|
|
contains
|
|
|
|
! modified to hard-code M2 constant to 1.0 (therefore removed)
|
|
subroutine VP_Dipole(x,y,z,Ax,Ay,Az)
|
|
real(rp), intent(in) :: x,y,z
|
|
real(rp), intent(out) :: Ax,Ay,Az
|
|
|
|
real(rp), dimension(NDIM) :: A,m,r,rhat
|
|
m = [0,0,1]
|
|
r = [x,y,z]
|
|
rhat = r/norm2(r)
|
|
|
|
A = cross(m,rhat)/(dot_product(r,r))
|
|
Ax = A(XDIR)
|
|
Ay = A(YDIR)
|
|
Az = A(ZDIR)
|
|
|
|
end subroutine VP_Dipole
|
|
|
|
end subroutine VMPIdipoleTest
|
|
|
|
!@test(npes=[2,5,9])
|
|
subroutine testDeepUpdate(this)
|
|
class (MpiTestMethod), intent(inout) :: this
|
|
|
|
end subroutine testDeepUpdate
|
|
|
|
end module testVMpiDeep
|
|
|