aboutsummaryrefslogtreecommitdiff
path: root/src/WaveToy.F90
diff options
context:
space:
mode:
authorgoodale <goodale@daab82bb-f315-4ad1-b6d0-9353ff8b6e27>2000-03-29 07:42:26 +0000
committergoodale <goodale@daab82bb-f315-4ad1-b6d0-9353ff8b6e27>2000-03-29 07:42:26 +0000
commit466765830e3ee8e6b4fb4b9a9c232171f7b55bf3 (patch)
treedfb5eed3b1766f6e4751f408012bbafa85dddc51 /src/WaveToy.F90
parente07ae62012ec2ccb3d283670a90252ea6d4bf6cd (diff)
Changed to use timelevels - thanks for putting the rotation in PUGH Ed.
Note that the no-boundary condition tests assumed that phi_new had zero at its boundaries. This is, or course, no longer true with timelevel rotation. Tom git-svn-id: http://svn.cactuscode.org/arrangements/CactusWave/WaveToyFreeF90/trunk@16 daab82bb-f315-4ad1-b6d0-9353ff8b6e27
Diffstat (limited to 'src/WaveToy.F90')
-rw-r--r--src/WaveToy.F9023
1 files changed, 10 insertions, 13 deletions
diff --git a/src/WaveToy.F90 b/src/WaveToy.F90
index 1a02912..6ba1074 100644
--- a/src/WaveToy.F90
+++ b/src/WaveToy.F90
@@ -46,6 +46,8 @@ subroutine WaveToyFreeF90_Evolution(CCTK_ARGUMENTS)
CCTK_REAL :: dx,dy,dz,dt
CCTK_REAL :: dx2,dy2,dz2,dt2
CCTK_REAL :: dx2i,dy2i,dz2i
+
+ CCTK_REAL :: factor
! Set up shorthands
! -----------------
@@ -70,6 +72,8 @@ subroutine WaveToyFreeF90_Evolution(CCTK_ARGUMENTS)
iend = cctk_lsh(1)-1
jend = cctk_lsh(2)-1
kend = cctk_lsh(3)-1
+
+ factor = 2*(1 - (dt2)*(dx2i + dy2i + dz2i))
! Do the evolution
! ----------------
@@ -77,10 +81,8 @@ subroutine WaveToyFreeF90_Evolution(CCTK_ARGUMENTS)
do j = jstart, jend
do i = istart, iend
- phi_tmp(i,j,k) = &
- 2.0*(1.0 - (dt2)*(dx2i + &
- dy2i + dz2i))*phi(i,j,k) - &
- phi_old(i,j,k) + (dt2) * &
+ phi_n(i,j,k) = factor*phi(i,j,k) - &
+ phi_p(i,j,k) + (dt2) * &
((phi(i+1,j,k)+phi(i-1,j,k))*dx2i &
+(phi(i,j+1,k)+phi(i,j-1,k))*dy2i &
+(phi(i,j,k+1)+phi(i,j,k-1))*dz2i)
@@ -91,17 +93,12 @@ subroutine WaveToyFreeF90_Evolution(CCTK_ARGUMENTS)
! Synchronize
! -----------
- call CCTK_SyncGroup(cctkGH,"wavetoyFreeF90::temps")
+ call CCTK_SyncGroup(cctkGH,"wavetoy::scalarevolve")
! Apply boundary conditions
! -------------------------
call WaveToyFreeF90_Boundaries(CCTK_PASS_FTOF)
- ! Update timeslices
- ! -----------------
- phi_old = phi
- phi = phi_tmp
-
end subroutine WaveToyFreeF90_Evolution
@@ -134,13 +131,13 @@ subroutine WaveToyFreeF90_Boundaries(CCTK_ARGUMENTS)
CCTK_REAL, parameter :: zero = 0.0
CCTK_REAL, parameter :: one = 1.0
- call CartSymBCGroup(ierr,cctkGH,"wavetoyfreef90::temps")
+ call CartSymBCGroup(ierr,cctkGH,"wavetoy::scalarevolve")
if (CCTK_EQUALS(bound,"flat")) then
- call FlatBCVar(ierr,cctkGH,sw,"wavetoyfreef90::phi_tmp")
+ call FlatBCVar(ierr,cctkGH,sw,"wavetoy::phi")
else if (CCTK_EQUALS(bound,"radiation")) then
call RadiativeBCVar(ierr,cctkGH,zero,one,sw, &
- "wavetoyfreef90::phi_tmp","wavetoy::phi")
+ "wavetoy::phi","wavetoy::phi")
end if
if (ierr < 0) then