Initial revision

darcs-hash:20010301114010-f6438-12fb8a9ffcc80e86c0a97e37b5b0dae0dbc59b79.gz

3 files changed, 201 insertions, 0 deletions

diff --git a/CarpetExtra/IDSpaceTimeToy/src/InitialData.F77 b/CarpetExtra/IDSpaceTimeToy/src/InitialData.F77
new file mode 100644
index 000000000..118f0d075
--- /dev/null
+++ b/CarpetExtra/IDSpaceTimeToy/src/InitialData.F77
@@ -0,0 +1,177 @@
+c     -*-Fortran-*-
+c     $Header: /home/eschnett/C/carpet/Carpet/CarpetExtra/IDSpaceTimeToy/src/InitialData.F77,v 1.8 2003/11/05 16:18:40 schnetter Exp $
+      
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Functions.h"
+#include "cctk_Parameters.h"
+      
+      subroutine IDSpaceTimeToy_InitialData (CCTK_ARGUMENTS)
+      
+      implicit none
+      
+      DECLARE_CCTK_ARGUMENTS
+      DECLARE_CCTK_PARAMETERS
+      DECLARE_CCTK_FUNCTIONS
+      
+      CCTK_REAL pi
+      CCTK_REAL omega
+      CCTK_REAL dt
+      CCTK_REAL x,y,z, r
+      integer i,j,k
+      
+      pi = 4*atan(1.d0)
+      
+      omega = sqrt(kx**2+ky**2+kz**2)
+      
+      dt = CCTK_DELTA_TIME
+      
+      if (CCTK_EQUALS(initial_data,"plane")) then 
+         
+         do k=1,cctk_lsh(3)
+            do j=1,cctk_lsh(2)
+               do i=1,cctk_lsh(1)
+                  
+                  x = cart3d_x(i,j,k)
+                  y = cart3d_y(i,j,k)
+                  z = cart3d_z(i,j,k)
+                  
+                  phi(i,j,k) = amplitude
+     $                 * cos((kx*x + ky*y + kz*z + omega*cctk_time) * pi)
+                  psi(i,j,k) = - amplitude
+     $                 * sin((kx*x + ky*y + kz*z + omega*cctk_time) * pi)
+     $                 * pi * omega
+                  
+                  phi_p(i,j,k) = amplitude
+     $                 * cos((kx*x + ky*y + kz*z + omega*(cctk_time - dt)) * pi)
+                  psi_p(i,j,k) = - amplitude
+     $                 * sin((kx*x + ky*y + kz*z + omega*(cctk_time - dt)) * pi)
+     $                 * pi * omega
+                  
+                  phi_p_p(i,j,k) = amplitude
+     $                 * cos((kx*x + ky*y + kz*z + omega*(cctk_time - 2*dt)) * pi)
+                  psi_p_p(i,j,k) = - amplitude
+     $                 * sin((kx*x + ky*y + kz*z + omega*(cctk_time - 2*dt)) * pi)
+     $                 * pi * omega
+                  
+               end do
+            end do
+         end do
+         
+      else if (CCTK_EQUALS(initial_data,"gaussian")) then 
+         
+         do k=1,cctk_lsh(3)
+            do j=1,cctk_lsh(2)
+               do i=1,cctk_lsh(1)
+                  
+                  r = spher3d_r(i,j,k)
+                  
+                  phi(i,j,k) = amplitude / r
+     $                 * exp(- (r - radius - cctk_time)**2 / sigma**2)
+                  psi(i,j,k) = phi(i,j,k)
+     $                 * 2 * (r - radius - cctk_time) / sigma**2
+                  
+                  phi_p(i,j,k) = amplitude / r
+     $                 * exp(- (r - radius - (cctk_time - dt))**2 / sigma**2)
+                  psi_p(i,j,k) = phi(i,j,k)
+     $                 * 2 * (r - radius - (cctk_time - dt)) / sigma**2
+                  
+                  phi_p_p(i,j,k) = amplitude / r
+     $                 * exp(- (r - radius - (cctk_time - 2*dt))**2 / sigma**2)
+                  psi_p_p(i,j,k) = phi(i,j,k)
+     $                 * 2 * (r - radius - (cctk_time - 2*dt)) / sigma**2
+                  
+               end do
+            end do
+         end do
+         
+      else if (CCTK_EQUALS(initial_data, "box")) then
+         
+c     Use kx,ky,kz as number of modes in each direction.
+         
+         do k=1,cctk_lsh(3)
+            do j=1,cctk_lsh(2)
+               do i=1,cctk_lsh(1)
+                  
+                  x = cart3d_x(i,j,k)
+                  y = cart3d_y(i,j,k)
+                  z = cart3d_z(i,j,k)
+                  
+                  phi(i,j,k) = amplitude
+     $                 * sin(kx * (x - 0.5d0) * pi)
+     $                 * sin(ky * (y - 0.5d0) * pi)
+     $                 * sin(kz * (z - 0.5d0) * pi)
+     $                 * cos(omega * cctk_time * pi)
+                  psi(i,j,k) = - amplitude
+     $                 * sin(kx * (x - 0.5d0) * pi)
+     $                 * sin(ky * (y - 0.5d0) * pi)
+     $                 * sin(kz * (z - 0.5d0) * pi)
+     $                 * sin(omega * cctk_time * pi)
+     $                 * omega * pi
+                  
+                  phi_p(i,j,k) = amplitude
+     $                 * sin(kx * (x - 0.5d0) * pi)
+     $                 * sin(ky * (y - 0.5d0) * pi)
+     $                 * sin(kz * (z - 0.5d0) * pi)
+     $                 * cos(omega * (cctk_time - dt) * pi)
+                  psi_p(i,j,k) = - amplitude
+     $                 * sin(kx * (x - 0.5d0) * pi)
+     $                 * sin(ky * (y - 0.5d0) * pi)
+     $                 * sin(kz * (z - 0.5d0) * pi)
+     $                 * sin(omega * (cctk_time - dt) * pi)
+     $                 * omega * pi
+                  
+                  phi_p_p(i,j,k) = amplitude
+     $                 * sin(kx * (x - 0.5d0) * pi)
+     $                 * sin(ky * (y - 0.5d0) * pi)
+     $                 * sin(kz * (z - 0.5d0) * pi)
+     $                 * cos(omega * (cctk_time - 2*dt) * pi)
+                  psi_p_p(i,j,k) = - amplitude
+     $                 * sin(kx * (x - 0.5d0) * pi)
+     $                 * sin(ky * (y - 0.5d0) * pi)
+     $                 * sin(kz * (z - 0.5d0) * pi)
+     $                 * sin(omega * (cctk_time - 2*dt) * pi)
+     $                 * omega * pi
+                  
+               end do
+            end do
+         end do
+         
+      else
+         
+         do k=1,cctk_lsh(3)
+            do j=1,cctk_lsh(2)
+               do i=1,cctk_lsh(1)
+                  
+                  phi(i,j,k) = 0
+                  psi(i,j,k) = 0
+                  
+                  phi_p(i,j,k) = 0
+                  psi_p(i,j,k) = 0
+                  
+                  phi_p_p(i,j,k) = 0
+                  psi_p_p(i,j,k) = 0
+                  
+               end do
+            end do
+         end do 
+         
+      end if
+      
+      if (hydrotoy_active.eq.1) then
+         
+         do k=1,cctk_lsh(3)
+            do j=1,cctk_lsh(2)
+               do i=1,cctk_lsh(1)
+                  
+                  psi(i,j,k) = psi(i,j,k) - u(i,j,k)
+                  psi_p(i,j,k) = psi_p(i,j,k) - u_p(i,j,k)
+                  psi_p_p(i,j,k) = psi_p_p(i,j,k) - u_p_p(i,j,k)
+                  
+               end do
+            end do
+         end do
+         
+      end if
+      
+      end 
diff --git a/CarpetExtra/IDSpaceTimeToy/src/Startup.F77 b/CarpetExtra/IDSpaceTimeToy/src/Startup.F77
new file mode 100644
index 000000000..0baf8e999
--- /dev/null
+++ b/CarpetExtra/IDSpaceTimeToy/src/Startup.F77
@@ -0,0 +1,15 @@
+c     -*-Fortran-*-
+c     $Header: /home/eschnett/C/carpet/Carpet/CarpetExtra/IDSpaceTimeToy/src/Startup.F77,v 1.1 2001/03/17 16:05:56 eschnett Exp $
+      
+#include "cctk.h"
+      
+      subroutine IDSpaceTimeToy_Startup
+      
+      implicit none
+      
+      integer ierr
+      
+      call CCTK_RegisterBanner
+     $     (ierr, "IDSpaceTimeToy: Initial data for two Scalar Fields")
+      
+      end 
diff --git a/CarpetExtra/IDSpaceTimeToy/src/make.code.defn b/CarpetExtra/IDSpaceTimeToy/src/make.code.defn
new file mode 100644
index 000000000..ca77c58cd
--- /dev/null
+++ b/CarpetExtra/IDSpaceTimeToy/src/make.code.defn
@@ -0,0 +1,9 @@
+# Main make.code.defn file for thorn IDSpaceTimeToy
+# $Header: /home/eschnett/C/carpet/Carpet/CarpetExtra/IDSpaceTimeToy/src/make.code.defn,v 1.1 2001/03/17 16:05:56 eschnett Exp $
+
+# Source files in this directory
+SRCS = InitialData.F77 Startup.F77
+
+# Subdirectories containing source files
+SUBDIRS = 
+