mirror of
https://github.com/JHUAPL/kaiju.git
synced 2026-01-10 07:38:00 -05:00
244 lines
7.8 KiB
Plaintext
244 lines
7.8 KiB
Plaintext
module testCases
|
|
use testHelper
|
|
use gamapp
|
|
use voltapp
|
|
use uservoltic
|
|
|
|
implicit none
|
|
|
|
|
|
contains
|
|
|
|
@before
|
|
subroutine firstSerial()
|
|
end subroutine firstSerial
|
|
|
|
@after
|
|
subroutine lastSerial()
|
|
end subroutine lastSerial
|
|
|
|
subroutine caseTestInit(gameraApp, caseFile)
|
|
class(gamApp_T), intent(inout) :: gameraApp
|
|
character(len=*), intent(in) :: caseFile
|
|
|
|
type(XML_Input_T) :: xmlInp
|
|
|
|
gameraApp%gOptions%userInitFunc => initUser
|
|
xmlInp = New_XML_Input(trim(caseFile),'Kaiju',.true.)
|
|
call gameraApp%InitModel(xmlInp)
|
|
|
|
end subroutine
|
|
|
|
@test
|
|
subroutine testLoop2D()
|
|
! a magnetic field should be moved along by the fluid
|
|
type(gamApp_T) :: gameraApp
|
|
real(rp), dimension(:,:,:), allocatable :: initialMagneticField
|
|
real(rp) :: totalError, initialTotalMagneticField, finalTotalMagneticField
|
|
|
|
call caseTestInit(gameraApp,'loop2d.xml')
|
|
|
|
! save initial magnetic field
|
|
allocate(initialMagneticField(gameraApp%Grid%isg:gameraApp%Grid%ieg+1,gameraApp%Grid%jsg:gameraApp%Grid%jeg+1,gameraApp%Grid%ksg:gameraApp%Grid%keg+1))
|
|
initialMagneticField = NORM2(gameraApp%State%magFlux,4)
|
|
|
|
do while ((gameraApp%Model%tFin - gameraApp%Model%t) > 1e-15)
|
|
call stepGamera(gameraApp)
|
|
end do
|
|
write(*,*) 'End time = ', gameraApp%Model%t
|
|
|
|
! compare final magnetic field to initial magnetic field
|
|
initialTotalMagneticField = SUM(initialMagneticField)
|
|
finalTotalMagneticField = SUM(NORM2(gameraApp%State%magFlux,4))
|
|
totalError = SUM(ABS(initialMagneticField - NORM2(gameraApp%State%magFlux,4)))
|
|
|
|
@assertLessThanOrEqual(abs(finalTotalMagneticField - initialTotalMagneticField)/initialTotalMagneticField, 0.005_rp, 'Lost more than 0.5% of total magnetic field')
|
|
@assertLessThanOrEqual(totalError/initialTotalMagneticField, 0.05_rp, 'Dissipation was more than 5% of total magnetic field')
|
|
|
|
end subroutine testLoop2D
|
|
|
|
@test
|
|
subroutine testAdv2D()
|
|
! density should be moved along by the fluid
|
|
type(gamApp_T) :: gameraApp
|
|
real(rp), dimension(:,:,:), allocatable :: initialGasDensity
|
|
real(rp) :: totalError, initialTotalDensity, finalTotalDensity
|
|
|
|
call caseTestInit(gameraApp, 'adv2d.xml')
|
|
|
|
! save initial density
|
|
allocate(initialGasDensity(gameraApp%Grid%isg:gameraApp%Grid%ieg,gameraApp%Grid%jsg:gameraApp%Grid%jeg,gameraApp%Grid%ksg:gameraApp%Grid%keg))
|
|
initialGasDensity = gameraApp%State%Gas(:,:,:,DEN,BLK)
|
|
|
|
do while ((gameraApp%Model%tFin - gameraApp%Model%t) > 1e-15)
|
|
call stepGamera(gameraApp)
|
|
end do
|
|
write(*,*) 'End time = ', gameraApp%Model%t
|
|
|
|
! compare final density to initial density
|
|
initialTotalDensity = SUM(initialGasDensity)
|
|
finalTotalDensity = SUM(gameraApp%State%Gas(:,:,:,DEN,BLK))
|
|
totalError = SUM(ABS(initialGasDensity - gameraApp%State%Gas(:,:,:,DEN,BLK)))
|
|
|
|
@assertLessThanOrEqual(abs(finalTotalDensity-initialTotalDensity)/initialTotalDensity, 0.0001_rp, 'Lost more than 0.01% of total density')
|
|
@assertLessThanOrEqual(totalError/initialTotalDensity, 0.05_rp, 'Dissipation was more than 5% of total density')
|
|
|
|
end subroutine testAdv2D
|
|
|
|
@test
|
|
subroutine testAlfven()
|
|
type(gamApp_T) :: gameraApp
|
|
real(rp), dimension(:,:,:), allocatable :: initialMagneticField
|
|
real(rp) :: totalError, initialTotalMagneticField, finalTotalMagneticField
|
|
|
|
call caseTestInit(gameraApp, 'alfven.xml')
|
|
|
|
! save initial magnetic field
|
|
allocate(initialMagneticField(gameraApp%Grid%isg:gameraApp%Grid%ieg+1,gameraApp%Grid%jsg:gameraApp%Grid%jeg+1,gameraApp%Grid%ksg:gameraApp%Grid%keg+1))
|
|
initialMagneticField = NORM2(gameraApp%State%magFlux,4)
|
|
|
|
do while ((gameraApp%Model%tFin - gameraApp%Model%t) > 1e-15)
|
|
call stepGamera(gameraApp)
|
|
end do
|
|
write(*,*) 'End time = ', gameraApp%Model%t
|
|
|
|
! compare final magnetic field to initial magnetic field
|
|
initialTotalMagneticField = SUM(initialMagneticField)
|
|
finalTotalMagneticField = SUM(NORM2(gameraApp%State%magFlux,4))
|
|
totalError = SUM(ABS(initialMagneticField - NORM2(gameraApp%State%magFlux,4)))
|
|
|
|
@assertLessThanOrEqual(abs(finalTotalMagneticField - initialTotalMagneticField)/initialTotalMagneticField, 0.001_rp, 'Lost more than 0.1% of total magnetic field')
|
|
@assertLessThanOrEqual(totalError/initialTotalMagneticField, 0.01_rp, 'Dissipation was more than 5% of total magnetic field')
|
|
end subroutine testAlfven
|
|
|
|
@test
|
|
subroutine testAdv1D()
|
|
type(gamApp_T) :: gameraApp
|
|
real(rp), dimension(:,:,:), allocatable :: sD,eD !Initial/Final densities
|
|
real(rp) :: L2err,L2Tol
|
|
|
|
L2Tol = 7.0e-3
|
|
|
|
call caseTestInit(gameraApp, 'adv1d.xml')
|
|
|
|
associate(Gr=>gameraApp%Grid,Gas=>gameraApp%State%Gas)
|
|
!Allocate arrays
|
|
allocate(sD(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke))
|
|
allocate(eD(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke))
|
|
!Save initial density
|
|
sD = Gas(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke,DEN,BLK)
|
|
|
|
!Now advance for one full advection cycle
|
|
do while ( (gameraApp%Model%t < gameraApp%Model%tFin) )
|
|
call stepGamera(gameraApp)
|
|
enddo
|
|
|
|
!Save final density
|
|
eD = Gas(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke,DEN,BLK)
|
|
L2err = norm2(sD-eD)/(1.0*Gr%Nip)
|
|
|
|
write(*,*) '1D Advection test'
|
|
write(*,*) 'Cells = ', Gr%Nip
|
|
write(*,*) 'L2 error = ', L2err
|
|
|
|
end associate
|
|
@assertLessThanOrEqual(L2err,L2Tol,"1D Advection test failed!")
|
|
end subroutine testAdv1D
|
|
|
|
@test
|
|
subroutine testAdv1D_8CENT()
|
|
type(gamApp_T) :: gameraApp
|
|
real(rp), dimension(:,:,:), allocatable :: sD,eD !Initial/Final densities
|
|
real(rp) :: L2err,L2Tol
|
|
|
|
L2Tol = 7.0e-3
|
|
|
|
call caseTestInit(gameraApp, 'adv1d_8cent.xml')
|
|
|
|
associate(Gr=>gameraApp%Grid,Gas=>gameraApp%State%Gas)
|
|
!Allocate arrays
|
|
allocate(sD(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke))
|
|
allocate(eD(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke))
|
|
!Save initial density
|
|
sD = Gas(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke,DEN,BLK)
|
|
|
|
!Now advance for one full advection cycle
|
|
do while ( (gameraApp%Model%t < gameraApp%Model%tFin) )
|
|
call stepGamera(gameraApp)
|
|
enddo
|
|
|
|
!Save final density
|
|
eD = Gas(Gr%is:Gr%ie,Gr%js:Gr%je,Gr%ks:Gr%ke,DEN,BLK)
|
|
L2err = norm2(sD-eD)/(1.0*Gr%Nip)
|
|
|
|
write(*,*) '1D Advection test 8CENT'
|
|
write(*,*) 'Cells = ', Gr%Nip
|
|
write(*,*) 'L2 error = ', L2err
|
|
|
|
end associate
|
|
@assertLessThanOrEqual(L2err,L2Tol,"1D Advection test 8CENT failed!")
|
|
end subroutine testAdv1D_8CENT
|
|
|
|
!Use Orszag-Tang to test for maintaining symmetry
|
|
@test
|
|
subroutine testOTSym()
|
|
type(gamApp_T) :: gameraApp
|
|
integer :: i,j,ip,jp
|
|
|
|
real(rp) :: dG
|
|
real(rp), dimension(NVAR) :: Gij,Gijp
|
|
|
|
call caseTestInit(gameraApp, 'ot2d.xml')
|
|
|
|
!Now advance for period of time specified in XML
|
|
do while ( (gameraApp%Model%t < gameraApp%Model%tFin) )
|
|
call stepGamera(gameraApp)
|
|
enddo
|
|
|
|
dG = 0.0
|
|
|
|
associate(Gr=>gameraApp%Grid,Gas=>gameraApp%State%Gas)
|
|
|
|
do j=Gr%js,Gr%je
|
|
do i=Gr%is,Gr%ie
|
|
!Find conjugate points
|
|
ip = Gr%ie-i+Gr%is
|
|
jp = Gr%je-j+Gr%js
|
|
|
|
Gij = Gas(i ,j ,Gr%ks,:,BLK)
|
|
Gijp = Gas(ip,jp,Gr%ks,:,BLK)
|
|
!Flip velocities
|
|
Gijp(MOMX:MOMY) = -Gijp(MOMX:MOMY)
|
|
dG = dG + norm2( Gij - Gijp )
|
|
|
|
enddo
|
|
enddo
|
|
dG = dG/(1.0*Gr%Nip*Gr%Njp)
|
|
|
|
write(*,*) '2D Advection test'
|
|
write(*,*) 'Cells = ', Gr%Nip,Gr%Njp
|
|
write(*,*) 'Symmetry violation = ', dG
|
|
@assertLessThanOrEqual(dG,1.0e-11_rp,"OTSym Test Failed")
|
|
|
|
end associate
|
|
|
|
end subroutine testOTSym
|
|
|
|
@test
|
|
subroutine testRemix()
|
|
type(voltApp_T) :: voltronApp
|
|
character(len=strLen) :: caseInput = 'cmiD.xml'
|
|
|
|
voltronApp%vOptions%gamUserInitFunc => initUser
|
|
call initVoltron(voltronApp, caseInput)
|
|
|
|
! run for one coupling interval
|
|
call stepVoltron(voltronApp, voltronApp%DeepDT)
|
|
|
|
write(*,*) 'End time = ', voltronApp%gApp%Model%t
|
|
|
|
end subroutine testRemix
|
|
|
|
end module testCases
|
|
|