Fixing up but it still doesn't work

git-svn-id: http://svn.cactuscode.org/arrangements/CactusWave/IDScalarWaveElliptic/trunk@8 41e88fdd-2190-4c69-9c84-4659c8cf322e

3 files changed, 34 insertions, 21 deletions

diff --git a/par/charge.par b/par/charge.par
index 948dbcb..0e15d40 100644
--- a/par/charge.par
+++ b/par/charge.par
@@ -10,13 +10,13 @@
 #   @enddesc
 # @@*/          
 
-ActiveThorns = "idscalarwaveelliptic idscalarwave time wavetoyf77 pugh cartgrid3d ioutil ioascii iobasic"
+ActiveThorns = "ellsor ellbase idscalarwaveelliptic idscalarwave time wavetoyf77 pugh cartgrid3d ioutil ioascii iobasic"
  
 time::dtfac = 0.5
 
 idscalarwave::initial_data = "charge"
 idscalarwaveelliptic::output_tmp = "yes"
-idscalarwaveelliptic::radius = 0.5
+idscalarwaveelliptic::radius = 5.5
 idscalarwaveelliptic::charge = 1
 
 wavetoyf77::bound = "radiation"
@@ -29,6 +29,7 @@ driver::global_nx = 60
 driver::global_ny = 60
 driver::global_nz = 60
 
+ellsor::maxit = 100
 cactus::cctk_itlast = 120
 
 IOBasic::outScalar_every = 2
diff --git a/src/SourceData.F b/src/SourceData.F
index 8bd4222..43bfe34 100644
--- a/src/SourceData.F
+++ b/src/SourceData.F
@@ -1,3 +1,11 @@
+ /*@@
+   @file      SourceData.F77
+   @date      
+   @author    Gabrielle Allen
+   @desc 
+              Elliptic initial data for wave equation
+   @enddesc 
+ @@*/
 
 #include "cctk.h"
 #include "cctk_arguments.h"
@@ -26,7 +34,6 @@ c     where Q is the total charge and R is the sphere radius
       pi = 4.0*atan(1.0)
 
       charge_factor = 4.0d0*pi*charge*3.0d0/(4.0d0*pi*radius**3)
-                  print *,charge_factor
 
       do k=1, cctk_lsh(3)
          do j=1, cctk_lsh(2)
@@ -34,8 +41,10 @@ c     where Q is the total charge and R is the sphere radius
                      
                Mcoeff(i,j,k) = 0.0d0
 
+               print *,r(i,j,k),radius
                if (r(i,j,k) <= radius) then
                   Ncoeff(i,j,k) = charge_factor
+                  print *,"Here"
                else
                   Ncoeff(i,j,k) = 0d0
                end if
@@ -58,14 +67,16 @@ c     where Q is the total charge and R is the sphere radius
       RelTol(2)=-1
       RelTol(3)=-1
 
-c     Call to elliptic solver will go here
+c     Call elliptic solver to fill out phi
+
       call Ell_LinFlatSolver(ierr,cctkGH, 
      &     iphi, 
      &     iMcoeff, iNcoeff, 
      &     AbsTol, RelTol,
      &     "sor")
 
-c     Set up last timestep
+c     Set up last timestep ... assume time symmetry
+
       do k=1, cctk_lsh(3)
          do j=1, cctk_lsh(2)
             do i=1, cctk_lsh(1)
@@ -78,23 +89,24 @@ c     Set up last timestep
 
 
 c     Output exact solution if required
-      do k=1, cctk_lsh(3)
-         do j=1, cctk_lsh(2)
-            do i=1, cctk_lsh(1)
-                     
-               if (r(i,j,k) >= radius) then
-                  temp(i,j,k) = charge/r(i,j,k)
-               else
-                  temp(i,j,k) = charge/(2.0d0*radius**3)*
-     &                 (3.0d0*radius**2-r(i,j,k)**2)
-               end if
-
-            end do
-         end do
-      end do
 
       if (output_tmp==1) then
-         call CCTK_OutputVarAsByMethod(ierr,cctkGH,
+        do k=1, cctk_lsh(3)
+           do j=1, cctk_lsh(2)
+              do i=1, cctk_lsh(1)
+                     
+                 if (r(i,j,k) >= radius) then
+                    temp(i,j,k) = charge/r(i,j,k)
+                 else
+                    temp(i,j,k) = charge/(2.0d0*radius**3)*
+     &                   (3.0d0*radius**2-r(i,j,k)**2)
+                 end if
+
+              end do
+           end do
+        end do
+
+        call CCTK_OutputVarAsByMethod(ierr,cctkGH,
      &        "idscalarwaveelliptic::temp",
      &        "IOASCII_1D","phi_exact")
       end if      
diff --git a/src/make.code.defn b/src/make.code.defn
index 0f9d963..0c6b937 100644
--- a/src/make.code.defn
+++ b/src/make.code.defn
@@ -2,7 +2,7 @@
 # $Header$
 
 # Source files in this directory
-SRCS = SourceData.F
+SRCS = SourceData.F77
 
 # Subdirectories containing source files
 SUBDIRS =