Included style suggestions from Erik Schnetter (Cactus/1104)

Used a more accurate representation of pi (hence new testsuite files)

Wrote some minimal documentation


git-svn-id: http://svn.cactuscode.org/arrangements/CactusWave/WaveBinarySource/trunk@36 b9286e40-80fe-41ab-903a-d6b447012e1e

5 files changed, 164 insertions, 166 deletions

diff --git a/README b/README
index 6f74b98..d2a2d46 100644
--- a/README
+++ b/README
@@ -3,5 +3,7 @@ Authors    : Gerd Lanfermann
 CVS info   : $Header$
 --------------------------------------------------------------------------
 
-Purpose of the thorn: provides a source term to the scalar field evolution
+Purpose of the thorn: 
+
+Provides a source term to the scalar field evolution
 for two rotating binary charges.
diff --git a/doc/README b/doc/README
index 666abab..7827791 100644
--- a/doc/README
+++ b/doc/README
@@ -1,7 +1,10 @@
+Cactus Code Thorn WaveBinarySource
+Authors    : Cactus Maintainers (cactusmaint@cactuscode.org)
+CVS info   : $Header$
+--------------------------------------------------------------------------
 
+Purpose of the thorn:
+
+This thorn adds a source term to the 3D scalar wave equation. The 
+source term represents two orbiting scalar sources.
 
-This thorn adds a source term to the evolution equation. It is
-scheduled to run after wavetoy evolution and adds data corresponding
-to two rotating binary charges. It simply overwrites the scalarfield
-grid variable where the binary sources are located. The binary sources 
-are not evolved. they are rotated based on the evolution time.
diff --git a/doc/documentation.tex b/doc/documentation.tex
index b097bbe..db988fc 100644
--- a/doc/documentation.tex
+++ b/doc/documentation.tex
@@ -1,65 +1,3 @@
-% *======================================================================*
-%  Cactus Thorn template for ThornGuide documentation
-%  Author: Ian Kelley
-%  Date: Sun Jun 02, 2002
-%  $Header$                                                             
-%
-%  Thorn documentation in the latex file doc/documentation.tex 
-%  will be included in ThornGuides built with the Cactus make system.
-%  The scripts employed by the make system automatically include 
-%  pages about variables, parameters and scheduling parsed from the 
-%  relevent thorn CCL files.
-%  
-%  This template contains guidelines which help to assure that your     
-%  documentation will be correctly added to ThornGuides. More 
-%  information is available in the Cactus UsersGuide.
-%                                                    
-%  Guidelines:
-%   - Do not change anything before the line
-%       % BEGIN CACTUS THORNGUIDE",
-%     except for filling in the title, author, date etc. fields.
-%   - You can define your own macros are OK, but they must appear after
-%     the BEGIN CACTUS THORNGUIDE line, and do not redefine standard 
-%     latex commands.
-%   - To avoid name clashes with other thorns, 'labels', 'citations', 
-%     'references', and 'image' names should conform to the following 
-%     convention:          
-%       ARRANGEMENT_THORN_LABEL
-%     For example, an image wave.eps in the arrangement CactusWave and 
-%     thorn WaveToyC should be renamed to CactusWave_WaveToyC_wave.eps
-%   - Graphics should only be included using the graphix package. 
-%     More specifically, with the "includegraphics" command. Do
-%     not specify any graphic file extensions in your .tex file. This 
-%     will allow us (later) to create a PDF version of the ThornGuide
-%     via pdflatex. |
-%   - References should be included with the latex "bibitem" command.  
-%   - For the benefit of our Perl scripts, and for future extensions, 
-%     please use simple latex.     
-%
-% *======================================================================* 
-% 
-% Example of including a graphic image:
-%    \begin{figure}[ht]
-% 	\begin{center}
-%    	   \includegraphics[width=6cm]{MyArrangement_MyThorn_MyFigure}
-% 	\end{center}
-% 	\caption{Illustration of this and that}
-% 	\label{MyArrangement_MyThorn_MyLabel}
-%    \end{figure}
-%
-% Example of using a label:
-%   \label{MyArrangement_MyThorn_MyLabel}
-%
-% Example of a citation:
-%    \cite{MyArrangement_MyThorn_Author99}
-%
-% Example of including a reference
-%   \bibitem{MyArrangement_MyThorn_Author99}
-%   {J. Author, {\em The Title of the Book, Journal, or periodical}, 1 (1999), 
-%   1--16. {\tt http://www.nowhere.com/}}
-%
-% *======================================================================* 
-
 % If you are using CVS use this line to give version information
 % $Header$
 
@@ -74,15 +12,15 @@
 \begin{document}
 
 % The author of the documentation
-\author{} 
+\author{Gabrielle Allen} 
 
 % The title of the document (not necessarily the name of the Thorn)
-\title{}
+\title{Orbiting Scalar Charges}
 
 % the date your document was last changed, if your document is in CVS, 
 % please us:
-%    \date{$ $Date$ $}
-\date{}
+\date{$ $Date$ $}
+
 
 \maketitle
 
@@ -102,34 +40,32 @@
 
 \section{Introduction}
 
-\section{Physical System}
-
-\section{Numerical Implementation}
-
-\section{Using This Thorn}
-
-\subsection{Obtaining This Thorn}
+This thorn provides a source term to the scalar field evolution
+for two rotating binary charges.
 
-\subsection{Basic Usage}
 
-\subsection{Special Behaviour}
-
-\subsection{Interaction With Other Thorns}
-
-\subsection{Support and Feedback}
-
-\section{History}
-
-\subsection{Thorn Source Code}
+\section{Physical System}
 
-\subsection{Thorn Documentation}
+The 3D scalar wave equation with a source term $\rho(t,x,y,z)$ is written
+$$
+\nabla \phi = 4 \pi \rho
+$$
 
-\subsection{Acknowledgements}
+Each scalar source with charge $Q$ and radius $R$ contributes 
+$$
+\rho = \frac{3Q}{4\pi R^3}
+$$
 
+\section{Numerical Implementation}
 
-\begin{thebibliography}{9}
+The only involved part of this thorn arise in working out where the
+sources are located (if at all) on each local grid for a multiprocessor
+run. The source terms are not numerically evolved, but are calculated 
+exactly, based on the physical time and their orbital velocity.
 
-\end{thebibliography}
+A routine is scheduled to run {\it after} the homogeneous equation
+for the scalar field has been evolved, and simply updates the value
+of the scalar field by adding on the source contribution.
 
 % Do not delete next line
 % END CACTUS THORNGUIDE
diff --git a/src/CoordinateStuff.c b/src/CoordinateStuff.c
index e9a2b9e..025fc01 100644
--- a/src/CoordinateStuff.c
+++ b/src/CoordinateStuff.c
@@ -1,3 +1,13 @@
+ /*@@
+   @file      CoordinateStuff.c
+   @date      
+   @author    Cactus Maintainers
+   @desc 
+   Extra functions for coordinate information
+   @enddesc 
+   @version $Header$
+ @@*/
+
 #include <stdio.h>
 #include <stdlib.h>
 #include <string.h>
@@ -9,12 +19,36 @@ static const char *rcsid = "$Header$";
 
 CCTK_FILEVERSION(CactusWave_WaveBinarySource_CoordinateStuff_c)
 
-int IndexFloorC(cGH *GH, CCTK_REAL coord_value, int d);
+/********************************************************************
+ *********************     Local Data Types   ***********************
+ ********************************************************************/
+
+
+/********************************************************************
+ ********************* Local Routine Prototypes *********************
+ ********************************************************************/
+
+
+/********************************************************************
+ *********************     Local Data   *****************************
+ ********************************************************************/
+
+
+/********************************************************************
+ *****************   External Routines Prototype ********************
+ ********************************************************************/
+
 int IndexCeilC(cGH *GH, CCTK_REAL coord_value, int d);
+int IndexFloorC(cGH *GH, CCTK_REAL coord_value, int d);
+
+
+/********************************************************************
+ **********************   External Routines  ************************
+ ********************************************************************/
 
 
  /*@@
-   @routine    IndexFloor
+   @routine    IndexFloorC
    @date       Fri Jan  7 10:34:29 2000
    @author     Gerd Lanfermann
    @desc 
@@ -24,6 +58,8 @@ int IndexCeilC(cGH *GH, CCTK_REAL coord_value, int d);
 
      Other return values: 
      -1 : gridindex not on this gridpatch
+     -2 : illegal dimension
+     -3 : coordinate value outside grid
    @enddesc 
    @calls     
    @calledby   
@@ -35,51 +71,43 @@ int IndexCeilC(cGH *GH, CCTK_REAL coord_value, int d);
 
 int IndexFloorC(cGH *GH, CCTK_REAL coord_value, int d)
 { 
-  int index_low,idummy;
+  int index_low;
   CCTK_REAL cmin,cmax;
 
-  if (d==0)
-  {
-    CCTK_CoordRange(GH,&cmin,&cmax, 1, "x", "cart3d");
-  }
-  else if (d==1)
-  {
-    CCTK_CoordRange(GH,&cmin,&cmax, 2, "y", "cart3d");
-  }
-  else if (d==2)
+  if (d>=0 || d<3)
   {
-    CCTK_CoordRange(GH,&cmin,&cmax, 3, "z", "cart3d");
+    CCTK_CoordRange(GH,&cmin,&cmax, d+1, NULL, "cart3d");
+
+    if ((coord_value<cmin)||(coord_value>cmax)) 
+    {
+      CCTK_VWarn(2,__LINE__,__FILE__,CCTK_THORNSTRING, 
+		 "IndexFloor: coordinate value outside grid [%f,%f]:  %f \n", 
+		 cmin, cmax, coord_value); 
+      index_low = -3;
+    }
+    else
+    {
+      index_low    = floor((coord_value-cmin)/ GH->cctk_delta_space[d])
+	  - GH->cctk_lbnd[d];
+
+      if (index_low<0 || index_low>=GH->cctk_lsh[d])
+      {
+	index_low = -1;
+      }
+    }
   }
   else 
   {
-    CCTK_WARN(0,"IndexFloor: dimension is not valid. Fortran:1..3  C:0..2");
+    CCTK_WARN(1,"IndexFloorC: dimension is not valid");
+    index_low = -2;
   }
  
-  if ((coord_value<cmin)||(coord_value>cmax)) 
-  {
-    CCTK_VWarn(2,__LINE__,__FILE__,CCTK_THORNSTRING, 
-	       "IndexFloor: coordinate value outside grid [%f,%f]:  %f \n", 
-	       cmin, cmax, coord_value); 
-    return(-3);
-  }
-  
-  idummy    = floor((coord_value-cmin)/ GH->cctk_delta_space[d]);
-  index_low = idummy-GH->cctk_lbnd[d];
-
-  if (index_low<0)
-  {
-    index_low = -1;
-  }
-  if (index_low>=GH->cctk_lsh[d])
-  {
-    index_low = -1;
-  } 
   return(index_low);
 }
 
 
  /*@@
-   @routine    IndexCeil
+   @routine    IndexCeilC
    @date       Fri Jan  7 10:34:29 2000
    @author     Gerd Lanfermann
    @desc 
@@ -88,6 +116,8 @@ int IndexFloorC(cGH *GH, CCTK_REAL coord_value, int d)
      including ghostzones.
      Other return values: 
      -1 : grid index not on this gridpatch
+     -2 : illegal dimension
+     -3 : coordinate value outside grid
    @enddesc 
    @calls     
    @calledby   
@@ -100,44 +130,36 @@ int IndexFloorC(cGH *GH, CCTK_REAL coord_value, int d)
 
 int IndexCeilC(cGH *GH, CCTK_REAL coord_value, int d)
 {  
-  int index_up,idummy;
+  int index_up;
   CCTK_REAL cmin,cmax;
 
-  if (d==0)
-  {
-    CCTK_CoordRange(GH,&cmin,&cmax, -1, "x", "cart3d");
-  }
-  else if (d==1)
+  if (d>=0 || d<3)
   {
-    CCTK_CoordRange(GH,&cmin,&cmax, -1, "y", "cart3d");
-  }
-  else if (d==2)
-  {
-    CCTK_CoordRange(GH,&cmin,&cmax, -1, "z", "cart3d");
+    CCTK_CoordRange(GH,&cmin,&cmax, d+1, NULL, "cart3d");
+
+    if ((coord_value<cmin)||(coord_value>cmax)) 
+    {
+      CCTK_VWarn(2,__LINE__,__FILE__,CCTK_THORNSTRING, 
+		 "IndexCeil: coordinate value outside grid [%f,%f]:  %f \n", 
+		 cmin, cmax, coord_value); 
+      index_up = -3;
+    }
+    else
+    {
+      index_up   = ceil((coord_value-cmin)/ GH->cctk_delta_space[d])
+	  -GH->cctk_lbnd[d];
+  
+      if (index_up<0 || index_up>=GH->cctk_lsh[d])
+      {
+	index_up = -1;
+      }
+    }
   }
   else 
   {
-    CCTK_WARN(0,"IndexCeilC: dimension is not valid. Fortran:1..3  C:0..2");
+    CCTK_WARN(1,"IndexCeilC: dimension is not valid");
+    index_up = -2;
   }
 
-  if ((coord_value<cmin)||(coord_value>cmax)) 
-  { 
-    CCTK_VWarn(2,__LINE__,__FILE__,CCTK_THORNSTRING,
-	       "IndexCeil: coordinate value outside grid [%f,%f]:  %f \n", 
-	       cmin,cmax, coord_value);
-    return(-3);
-  }
-  
-  idummy   = ceil((coord_value-cmin)/ GH->cctk_delta_space[d]);
-  index_up = idummy-GH->cctk_lbnd[d];
-  
-  if (index_up<0)
-  {
-    index_up = -1;
-  }
-  if (index_up>=GH->cctk_lsh[d])
-  {
-    index_up = -1;
-  } 
   return(index_up);
 }
diff --git a/src/WaveBinary.c b/src/WaveBinary.c
index 269a024..ab8f3ec 100644
--- a/src/WaveBinary.c
+++ b/src/WaveBinary.c
@@ -1,8 +1,17 @@
+ /*@@
+   @file      WaveBinary.c
+   @date      
+   @author    Cactus Maintainers
+   @desc 
+   Add a binary source term to the 3D scalar wave equation
+   @enddesc 
+   @version $Header$
+ @@*/
+
 #include <stdio.h>
 #include <stdlib.h>
 #include <math.h>
 
-/* Using Cactus infrastructure */
 #include "cctk.h" 
 #include "cctk_Parameters.h"
 #include "cctk_Arguments.h"
@@ -11,11 +20,35 @@ static const char *rcsid = "$Header$";
 
 CCTK_FILEVERSION(CactusWave_WaveBinarySource_WaveBinary_c)
 
+/********************************************************************
+ *********************     Local Data Types   ***********************
+ ********************************************************************/
+
+
+/********************************************************************
+ ********************* Local Routine Prototypes *********************
+ ********************************************************************/
 
-void WaveBinaryC(CCTK_ARGUMENTS);
 int IndexCeilC(cGH *GH, CCTK_REAL coord_value, int d);
 int IndexFloorC(cGH *GH, CCTK_REAL coord_value, int d);
 
+
+/********************************************************************
+ *********************     Local Data   *****************************
+ ********************************************************************/
+
+
+/********************************************************************
+ *****************   External Routines Prototype ********************
+ ********************************************************************/
+
+void WaveBinaryC(CCTK_ARGUMENTS);
+
+
+/********************************************************************
+ **********************   External Routines  ************************
+ ********************************************************************/
+
 void WaveBinaryC(CCTK_ARGUMENTS)
 {
   
@@ -38,16 +71,18 @@ void WaveBinaryC(CCTK_ARGUMENTS)
   /* some flags */
   int onthisproc, xproc, sign;
   static int firstcall;
-  
+
   CCTK_REAL xs,ys,zs,rad;
   CCTK_REAL zmin,zmax,charge_factor;
+  const CCTK_REAL binary_pi=3.141592653589793;
 
   /* Because binary_charge and binary_size are a steerable parameters now,
      charge_factor needs to be recomputed at every iteration. */
   charge_factor = 3.0*binary_charge/
-    (4.0*3.1415*binary_size*binary_size*binary_size);
+    (4.0*binary_pi*binary_size*binary_size*binary_size);
 
   /* Initialize the range arrays */
+
   for (d=0;d<3;d++)
   {
     mingp[d] = 0;
@@ -68,7 +103,7 @@ void WaveBinaryC(CCTK_ARGUMENTS)
         
 
   /*  we have two charges opposite, origin is the center
-      sign mulitplies the charge xy-positions by +1 and -1 */
+      sign multiplies the charge xy-positions by +1 and -1 */
   for (sign=1; sign>=-1; sign-=2)
   {
 
@@ -96,8 +131,8 @@ void WaveBinaryC(CCTK_ARGUMENTS)
     for (d=0;d<3;d++)
     {   
      
-      /* check where our binaries are:
-         - the can be on our grid chunk completly
+      /* check where our binary sources are:
+         - they can be on our grid chunk completly
          - we have a lower part of the source on our grid chunk
          - we have a upper part
          - the source covers a grid chunk completetly 
@@ -157,7 +192,7 @@ void WaveBinaryC(CCTK_ARGUMENTS)
              maxc[0],maxc[1],maxc[2]);
     }
       
-    /* we have the source, no loop over the vindex boudary */
+    /* we have the source, now loop over the vindex boundary */
     if (onthisproc==3)
     {
       for (k=lowerloc[2];k<=upperloc[2];k++)