Developing

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

4 files changed, 41 insertions, 18 deletions

diff --git a/interface.ccl b/interface.ccl
index 1067c00..e937f0f 100644
--- a/interface.ccl
+++ b/interface.ccl
@@ -12,7 +12,10 @@ cctk_real ellcoeffs type = GF
   Mcoeff,Ncoeff
 } "Coefficients of the elliptic equation"
 
- 
+cctk_real outvars type = GF
+{
+  temp
+} "Temporary vars" 
  
 
 
diff --git a/param.ccl b/param.ccl
index 0b40989..5b5f597 100644
--- a/param.ccl
+++ b/param.ccl
@@ -21,4 +21,7 @@ REAL charge "Charge of uniformly charged sphere"
  0:* ::
 } 1.0
 
+BOOLEAN output_tmp "Output intermediate variables"
+{
+} "no"
 
diff --git a/schedule.ccl b/schedule.ccl
index 1205472..e410728 100644
--- a/schedule.ccl
+++ b/schedule.ccl
@@ -3,7 +3,7 @@
 
 schedule UniformCharge at CCTK_INITIAL
 {
-  STORAGE:       wavetoy::scalarevolve,ellcoeffs
+  STORAGE:       wavetoy::scalarevolve,ellcoeffs,outvars
   COMMUNICATION: wavetoy::scalarevolve
   LANG:          Fortran
 } "Initial data (Uniform Charge) for 3D wave equation"
diff --git a/src/SourceData.F b/src/SourceData.F
index d16c0de..2913efb 100644
--- a/src/SourceData.F
+++ b/src/SourceData.F
@@ -16,19 +16,17 @@ c     where Q is the total charge and R is the sphere radius
       DECLARE_CCTK_FARGUMENTS
       DECLARE_CCTK_PARAMETERS
 
-      CCTK_REAL tolerance,pi
+      CCTK_REAL pi
       CCTK_REAL AbsTol(3), RelTol(3)
 
       integer iphi,iMcoeff,iNcoeff
-      integer i,j,k
+      integer i,j,k,ierr
       CCTK_REAL charge_factor
 
       pi = 4.0*atan(1.0)
 
       charge_factor = 4.0d0*pi*charge*3.0d0/(4.0d0*pi*radius**3)
 
-      tolerance = 1.0d-5
-
       do k=1, cctk_lsh(3)
          do j=1, cctk_lsh(2)
             do i=1, cctk_lsh(1)
@@ -51,8 +49,6 @@ c     where Q is the total charge and R is the sphere radius
       call CCTK_VarIndex (iNcoeff,  "idscalarwaveelliptic::Ncoeff")
       call CCTK_VarIndex (iphi,"wavetoy::phi")
       
-      write (*,*) iMcoeff, iNcoeff, iphi
-
       AbsTol(1)=1.0d-5
       AbsTol(2)=1.0d-5
       AbsTol(3)=1.0d-5
@@ -62,30 +58,51 @@ c     where Q is the total charge and R is the sphere radius
       RelTol(3)=-1
 
 c     Call to elliptic solver will go here
-c      call Ell_LinConfMetricSolver(cctkGH, 
-c     $     metpsi_index, afield_index, Mlin_index, Nsrc_index, AbsTol, RelTol,
-c     $     "petsc")
+      call Ell_LinFlatSolver(cctkGH, 
+     &     iphi, 
+     &     iMcoeff, iNcoeff, 
+     &     AbsTol, RelTol,
+     &     "sor")
+
+c     Set up last timestep
+      do k=1, cctk_lsh(3)
+         do j=1, cctk_lsh(2)
+            do i=1, cctk_lsh(1)
+                     
+               phi_old(i,j,k) = phi(i,j,k)
+
+            end do
+         end do
+      end do
+
 
-c     Cheat, and use exact solution
+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
-                  phi(i,j,k) = charge/r(i,j,k)
+                  temp(i,j,k) = charge/r(i,j,k)
                else
-                  phi(i,j,k) = charge/(2.0d0*radius**3)*
+                  temp(i,j,k) = charge/(2.0d0*radius**3)*
      &                 (3.0d0*radius**2-r(i,j,k)**2)
                end if
 
-c              Should put in an option for time symmetric first step
-               phi_old(i,j,k) = phi(i,j,k)
-
             end do
          end do
       end do
 
-      
+      if (output_tmp==1) then
+         call CCTK_OutputVarAsByMethod(ierr,cctkGH,
+     &        "idscalarwaveelliptic::temp",
+     &        "IOASCII_1D","phi_exact")
+      end if      
 
       return
       end
+
+
+
+
+
+