Initial import of new Einstein stuff. This has the new thorns, but not all

are fully functional yet.  

When I have completed this stage I'll send an updated spec out with some 
questions which have arisen during this process.

Please don't import anything new without checking with me first, as I want 
to play games on the server copying cvs files around to preserve histories 
on files which are only minimally touched.

Tom



git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/ADMAnalysis/trunk@2 80bd93c7-81bc-453a-9e3f-619c0b4f6fe4

11 files changed, 922 insertions, 0 deletions

diff --git a/README b/README
new file mode 100644
index 0000000..4d9f5eb
--- /dev/null
+++ b/README
@@ -0,0 +1,31 @@
+Cactus Code Thorn ADMAnalysis
+Authors    : Tom Goodale
+CVS info   : $Header$
+--------------------------------------------------------------------------
+
+Purpose of the thorn:
+
+This thorn does basic analysis of the metric and extrinsic curvature 
+tensors.
+
+It calculates
+
+  the trace of the extrinsic curvature ({\bf trK})
+
+  the determinant of the metric ({\bf detg})
+
+  the components of the metric in spherical coordinates 
+  (grr,grq,grp,gqq,gqp,gpp)
+
+  the components of the extrinsic curvature in spherical coordinates 
+  (Krr,Krq,Krp,Kqq,Kqp,Kpp)
+
+if output is requested for them.
+
+(q refers to the theta coordinate and p to the phi coordinate.)
+
+If the parameter 'rsquared_in_sphm' is true, the thorn squares the
+radial coordinate when it uses it for the transformation.
+
+This thorn knows how to handle 'physical' and 'static conformal'
+metric types.
diff --git a/doc/documentation.tex b/doc/documentation.tex
new file mode 100644
index 0000000..e3730fc
--- /dev/null
+++ b/doc/documentation.tex
@@ -0,0 +1,46 @@
+\documentclass{article}
+\begin{document}
+
+\title{ADMAnalysis}
+\author{Tom Goodale et al}
+\date{April 2002}
+\maketitle
+
+\abstract{Basic analysis of the metric and extrinsic curvature tensors}
+
+\section{Purpose}
+
+This thorn calculates
+
+\begin{itemize}
+\item
+The trace of the extrinsic curvature ({\bf trK}).
+\item
+The determinant of the metric ({\bf detg}).
+\item
+The components of the metric in spherical coordinates 
+({\bf grr,grq,grp,gqq,gqp,gpp}).
+\item
+The components of the extrinsic curvature in spherical coordinates 
+({\bf Krr,Krq,Krp,Kqq,Kqp,Kpp}).
+\end{itemize}
+
+if output is requested for them.
+
+\section{Comments}
+
+If the parameter {\bf rsquared\_in\_sphm} is set, it squares the
+radial coordinate before applying the tranformation.
+
+In the spherical transormation, the $\theta$ coordinate is referred to
+as {\bf q} and the $\phi$ as {\bf p}.
+
+This thorn knows how to handle `physical' and `static conformal'
+metric types.
+
+% Automatically created from the ccl files by using gmake thorndoc
+\include{interface}
+\include{param}
+\include{schedule}
+
+\end{document}
diff --git a/interface.ccl b/interface.ccl
new file mode 100644
index 0000000..5a92716
--- /dev/null
+++ b/interface.ccl
@@ -0,0 +1,28 @@
+# Interface definition for thorn ADMAnalysis
+# $Header$
+
+implements: ADMAnalysis
+
+inherits: ADMBase, StaticConformal, Grid
+
+# For evaltrK
+REAL trace_of_K TYPE = GF
+{
+  trK 
+} "trace of extrinsic curvature"
+
+REAL detofg TYPE = GF
+{
+  detg 
+} "determinant of the conformal metric"
+
+# For carttoshpere (p=phi, q=theta)
+REAL spherical_metric TYPE = GF
+{
+  grr,gqq,gpp,grq,grp,gqp
+} "Metric in spherical coordinates"
+
+REAL spherical_curv TYPE = GF
+{
+  Krr,Kqq,Kpp,Krq,Krp,Kqp
+} "extrinsic curvature in spherical coordinates"
diff --git a/param.ccl b/param.ccl
new file mode 100644
index 0000000..3a9e0bd
--- /dev/null
+++ b/param.ccl
@@ -0,0 +1,12 @@
+# Parameter definitions for thorn ADMAnalysis
+# $Header$
+
+shares: ADMBase
+
+USES KEYWORD metric_type
+
+private:
+
+BOOLEAN rsquared_in_sphm "Use r^2 when calculating spherical components ?"
+{
+} "no"
diff --git a/schedule.ccl b/schedule.ccl
new file mode 100644
index 0000000..0bc95ec
--- /dev/null
+++ b/schedule.ccl
@@ -0,0 +1,29 @@
+# Schedule definitions for thorn ADMAnalysis
+# $Header$
+
+SCHEDULE ADMAnalysis_ParamCheck AT CCTK_PARAMCHECK
+{
+  LANG:C
+} "Check that the metric_type is recognised"
+
+
+SCHEDULE ADMAnalysis_EvaltrK AT CCTK_ANALYSIS
+{
+  STORAGE: trace_of_K,detofg
+  LANG: C
+  TRIGGERS: trace_of_K,detofg
+} "Compute the trace of the extrinsic curvature and the determinant of the metric"
+
+SCHEDULE ADMAnalysis_MetricCartToSphere AT CCTK_ANALYSIS
+{
+  STORAGE: spherical_metric
+  LANG: C
+  TRIGGERS: spherical_metric
+} "Calculate the spherical metric in r,theta(q), phi(p)"
+
+SCHEDULE ADMAnalysis_CurvCartToSphere AT CCTK_ANALYSIS
+{
+  STORAGE: spherical_curv
+  LANG: C
+  TRIGGERS: spherical_curv
+} "Calculate the spherical ex. curvature in r, theta(q), phi(p)"
diff --git a/src/ADMAnalysis.h b/src/ADMAnalysis.h
new file mode 100644
index 0000000..d11aeb7
--- /dev/null
+++ b/src/ADMAnalysis.h
@@ -0,0 +1,58 @@
+ /*@@
+   @header    ADMAnalysis.h
+   @date      Thu Apr 25 16:36:20 2002
+   @author    Tom Goodale
+   @desc 
+   
+   @enddesc
+   @version $Header$
+ @@*/
+
+#ifndef _ADMAnalysis_H_
+#define _ADMAnalysis_H_ 1
+
+#ifdef __cplusplus
+extern "C" 
+{
+#endif
+
+void ADMAnalysis_CartToSphere(const CCTK_INT *lsh,
+                              int r2norm,
+                              const CCTK_REAL *x,
+                              const CCTK_REAL *y,
+                              const CCTK_REAL *z,
+                              const CCTK_REAL *r,
+                              const CCTK_REAL *cart_xx,
+                              const CCTK_REAL *cart_xy,
+                              const CCTK_REAL *cart_xz,
+                              const CCTK_REAL *cart_yy,
+                              const CCTK_REAL *cart_yz,
+                              const CCTK_REAL *cart_zz,
+                              CCTK_REAL *sphere_rr,
+                              CCTK_REAL *sphere_rq,
+                              CCTK_REAL *sphere_rp,
+                              CCTK_REAL *sphere_qq,
+                              CCTK_REAL *sphere_qp,
+                              CCTK_REAL *sphere_pp);
+
+void ADMAnalysis_Trace(const CCTK_INT *lsh,
+                       const CCTK_REAL *g11,
+                       const CCTK_REAL *g12,
+                       const CCTK_REAL *g13,
+                       const CCTK_REAL *g22,
+                       const CCTK_REAL *g23,
+                       const CCTK_REAL *g33,
+                       CCTK_REAL *tensor11,
+                       CCTK_REAL *tensor12,
+                       CCTK_REAL *tensor13,
+                       CCTK_REAL *tensor22,
+                       CCTK_REAL *tensor23,
+                       CCTK_REAL *tensor33,
+                       CCTK_REAL *trace,
+                       CCTK_REAL *detg);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* _ADMAnalysis_H_ */
diff --git a/src/Analysis.c b/src/Analysis.c
new file mode 100644
index 0000000..f154e76
--- /dev/null
+++ b/src/Analysis.c
@@ -0,0 +1,183 @@
+ /*@@
+   @file      Analysis.c
+   @date      Thu Apr 25 18:46:27 2002
+   @author    Tom Goodale
+   @desc 
+   Routines to do the various analysis tasks.
+   @enddesc
+   @version $Header$
+ @@*/
+
+#include "cctk.h"
+
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+
+#include "ADMAnalysis.h"
+
+static const char *rcsid = "$Header$";
+
+CCTK_FILEVERSION(CactusEinstein_ADMAnalysis_Analysis_c)
+
+/********************************************************************
+ *********************     Local Data Types   ***********************
+ ********************************************************************/
+
+#define SQR(a) ((a)*(a))
+
+/********************************************************************
+ ********************* Local Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ ***************** Scheduled Routine Prototypes *********************
+ ********************************************************************/
+
+void ADMAnalysis_EvaltrK(CCTK_ARGUMENTS);
+
+void ADMAnalysis_MetricCartToSphere(CCTK_ARGUMENTS);
+
+void ADMAnalysis_CurvCartToSphere(CCTK_ARGUMENTS);
+
+/********************************************************************
+ ********************* Other Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     Local Data   *****************************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     External Routines   **********************
+ ********************************************************************/
+
+ /*@@
+   @routine    ADMAnalysis_EvaltrK
+   @date       Thu Apr 25 19:11:32 2002
+   @author     Tom Goodale
+   @desc 
+   Scheduled routine to evaluate the trace of the extrinsic curvature
+   and the determinant of the metric.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+ 
+@@*/
+void ADMAnalysis_EvaltrK(CCTK_ARGUMENTS)
+{
+  int i;
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+
+  ADMAnalysis_Trace(cctk_lsh,
+                    gxx,
+                    gxy,
+                    gxz,
+                    gyy,
+                    gyz,
+                    gzz,
+                    kxx,
+                    kxy,
+                    kxz,
+                    kyy,
+                    kyz,
+                    kzz,
+                    trK,
+                    detg);
+
+  /* If we have a conformal metric, must convert trK to non-conformal form */
+  if(CCTK_EQUALS(metric_type, "static conformal") && *conformal_state > 0)
+  {
+    for(i = 0; i< cctk_lsh[0]*cctk_lsh[1]*cctk_lsh[2]; i++)  
+    {
+      trK[i] = trK[i] /(SQR(psi[i])*SQR(psi[i]));
+    }
+  }
+}
+
+/*@@
+  @routine    ADMAnalysis_MetricCartToSphere
+  @date       Thu Apr 25 19:11:32 2002
+  @author     Tom Goodale
+  @desc 
+  Scheduled routine to evaluate the components of the metric
+  in spherical coordinates
+  @enddesc 
+  @calls     
+  @calledby   
+  @history 
+  
+  @endhistory 
+  
+@@*/
+void ADMAnalysis_MetricCartToSphere(CCTK_ARGUMENTS)
+{
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+
+  ADMAnalysis_CartToSphere(cctk_lsh,
+                           rsquared_in_sphm,
+                           x,
+                           y,
+                           z,
+                           r,
+                           gxx,
+                           gxy,
+                           gxz,
+                           gyy,
+                           gyz,
+                           gzz,
+                           grr,
+                           grq,
+                           grp,
+                           gqq,
+                           gqp,
+                           gpp);
+}
+
+/*@@
+  @routine    ADMAnalysis_MetricCartToSphere
+  @date       Thu Apr 25 19:11:32 2002
+  @author     Tom Goodale
+  @desc 
+  Scheduled routine to evaluate the components of the extrinsic curvature
+  in spherical coordinates
+  @enddesc 
+  @calls     
+  @calledby   
+  @history 
+  
+  @endhistory 
+  
+@@*/
+void ADMAnalysis_CurvCartToSphere(CCTK_ARGUMENTS)
+{
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+
+  ADMAnalysis_CartToSphere(cctk_lsh,
+                           rsquared_in_sphm,
+                           x,
+                           y,
+                           z,
+                           r,
+                           kxx,
+                           kxy,
+                           kxz,
+                           kyy,
+                           kyz,
+                           kzz,
+                           Krr,
+                           Krq,
+                           Krp,
+                           Kqq,
+                           Kqp,
+                           Kpp);
+}
+/********************************************************************
+ *********************     Local Routines   *************************
+ ********************************************************************/
+
diff --git a/src/CartToSphere.c b/src/CartToSphere.c
new file mode 100644
index 0000000..00d074a
--- /dev/null
+++ b/src/CartToSphere.c
@@ -0,0 +1,273 @@
+ /*@@
+   @file      CartToSphere.c
+   @date      Thu Apr 25 16:04:53 2002
+   @author    Tom Goodale
+   @desc 
+   
+   @enddesc
+   @version $Header
+ @@*/
+
+#include "cctk.h"
+
+#include <math.h>
+
+#include "ADMAnalysis.h"
+
+static const char *rcsid = "$Header$";
+
+CCTK_FILEVERSION(CactusEinstein_ADMAnalysis_CartToSphere_c)
+
+/********************************************************************
+ *********************     Local Data Types   ***********************
+ ********************************************************************/
+
+#define SQR(a) ((a)*(a))
+
+/********************************************************************
+ ********************* Local Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ ***************** Scheduled Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ ********************* Other Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     Local Data   *****************************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     External Routines   **********************
+ ********************************************************************/
+
+ /*@@
+   @routine    ADMAnalysis_CartToSphere
+   @date       Thu Apr 25 16:17:58 2002
+   @author     Tom Goodale
+   @desc 
+   Calculates the spherical components (r,q,p) corresponding to 
+   the cartesian components (xyz) of a symmetric spatial tensor.   
+   q corresponds to theta, and p to phi.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+   @hdate Thu Apr 25 16:20:07 2002 @hauthor Tom Goodale
+   @hdesc Took old metric_ and curv_cartosphere routines and
+          made this generic one.
+   @endhistory 
+   @var     lsh
+   @vdesc   the local shape of the grid
+   @vtype   const int *
+   @vio     in
+   @var     r2norm
+   @vdesc   a flag which is set if the radius variable shouldbe squared
+   @vtype   int
+   @vio     in
+   @var     x
+   @vdesc   the x coordinate
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     y
+   @vdesc   the y coordinate
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     z
+   @vdesc   the z coordinate
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     r
+   @vdesc   the r coordinate
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     cart_xx
+   @vdesc   the xx component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     cart_xy
+   @vdesc   the xy component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     cart_xz
+   @vdesc   the xz component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     cart_yy
+   @vdesc   the yy component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     cart_yz
+   @vdesc   the yz component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     cart_zz
+   @vdesc   the zz component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     sphere_rr
+   @vdesc   the rr component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     sphere_rq
+   @vdesc   the rq component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     sphere_rp
+   @vdesc   the rp component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     sphere_qq
+   @vdesc   the qq component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     sphere_qp
+   @vdesc   the qp component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     sphere_pp
+   @vdesc   the pp component of the tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+
+ @@*/
+void ADMAnalysis_CartToSphere(const CCTK_INT *lsh,
+                              int r2norm,
+                              const CCTK_REAL *x,
+                              const CCTK_REAL *y,
+                              const CCTK_REAL *z,
+                              const CCTK_REAL *r,
+                              const CCTK_REAL *cart_xx,
+                              const CCTK_REAL *cart_xy,
+                              const CCTK_REAL *cart_xz,
+                              const CCTK_REAL *cart_yy,
+                              const CCTK_REAL *cart_yz,
+                              const CCTK_REAL *cart_zz,
+                              CCTK_REAL *sphere_rr,
+                              CCTK_REAL *sphere_rq,
+                              CCTK_REAL *sphere_rp,
+                              CCTK_REAL *sphere_qq,
+                              CCTK_REAL *sphere_qp,
+                              CCTK_REAL *sphere_pp)
+{
+  int i;
+  CCTK_REAL cost;
+  CCTK_REAL sint;
+  CCTK_REAL sinp;
+  CCTK_REAL cosp;
+  CCTK_REAL r_squared;
+  CCTK_REAL sxy;
+
+  CCTK_REAL    txx,txy,txz,tyy,tyz,tzz;
+    
+  /* loop over all the grid */
+  for(i = 0; i < lsh[0]*lsh[1]*lsh[2]; i++)  
+  {
+    
+     txx = cart_xx[i];
+     txy = cart_xy[i];
+     txz = cart_xz[i];
+     tyy = cart_yy[i];
+     tyz = cart_yz[i];
+     tzz = cart_zz[i];
+     r_squared  = r[i];
+     sxy = sqrt( SQR(x[i]) + SQR(y[i]));
+       
+  /* be careful with r=0 and xy plane */
+     if (r_squared==0.0) 
+     {
+       cost = 1.0;
+       sint = 0.0; 
+       sinp = 0.0;
+       cosp = 1.0;
+     } 
+     else if (sxy==0) 
+     {
+       cost = 1.0; 
+       sint = 0.0;
+       sinp = 0.0;
+       cosp = 1.0;
+     } 
+     else 
+     {
+       cost = z[i]/r_squared;
+       sint = sxy/r_squared;
+       sinp = y[i]/sxy;
+       cosp = x[i]/sxy;
+     }
+     
+     sphere_rr[i]= 
+         tyy*SQR(sinp)*SQR(sint)+
+         2*cosp*txy*sinp*SQR(sint)+
+         SQR(cosp)*txx*SQR(sint)+
+         2*cost*tyz*sinp*sint+
+         2*cosp*cost*txz*sint+
+         SQR(cost)*tzz;
+               
+     sphere_qq[i] = 
+         (tzz*SQR(sint)+
+          (-2*cost*tyz*sinp-
+           2*cosp*cost*txz)*sint+
+          SQR(cost)*tyy*SQR(sinp)+
+          2*cosp*SQR(cost)*txy*sinp
+          +SQR(cosp)*SQR(cost)*txx);
+     
+     if (r2norm) 
+     {
+       sphere_qq[i] *= SQR(r[i]);
+     }
+     
+     sphere_pp[i] = 
+         (txx*SQR(sinp)-
+          2*cosp*txy*sinp+
+          SQR(cosp)*tyy);
+     
+     if (r2norm)  
+     {
+       sphere_pp[i] *= SQR(r[i]) * SQR(sint);
+     }
+               
+     sphere_rq[i] = 
+         (cost*tyy*SQR(sinp)*sint+
+          2*cosp*cost*txy*sinp*sint-
+          cost*tzz*sint+
+          SQR(cosp)*cost*txx*sint+
+          2*SQR(cost)*tyz*sinp-
+          tyz*sinp+
+          2*cosp*SQR(cost)*txz-
+          cosp*txz)*r[i];
+     
+     if (r2norm) 
+     {
+       sphere_rq[i] *= r[i];
+     }
+
+     sphere_rp[i] = 
+         ((-txy*SQR(sinp)+
+           (cosp*tyy-cosp*txx)*sinp+
+           SQR(cosp)*txy)*sint-
+          cost*txz*sinp+cosp*cost*tyz);
+     
+     if (r2norm)  
+     {
+       sphere_rp[i] *= r[i] * sint;
+     }
+     
+     sphere_qp[i] =
+         ((txz*sinp-cosp*tyz)*sint+
+          cost*(-txy*SQR(sinp)+
+                cosp*(tyy-txx)*sinp+SQR(cosp)*txy));
+
+     if (r2norm) 
+     {
+       sphere_qp[i] *= SQR(r[i]) * sint;
+     }
+  }
+}
+/********************************************************************
+ *********************     Local Routines   *************************
+ ********************************************************************/
+
diff --git a/src/ParamCheck.c b/src/ParamCheck.c
new file mode 100644
index 0000000..918bcee
--- /dev/null
+++ b/src/ParamCheck.c
@@ -0,0 +1,73 @@
+ /*@@
+   @file      ParamCheck.c
+   @date      Thu Apr 25 19:02:51 2002
+   @author    Tom Goodale
+   @desc 
+   Parameter checking stuff for ADMAnalysis
+   @enddesc
+   @version $Header$
+ @@*/
+
+#include "cctk.h"
+
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+
+static const char *rcsid = "$Header$";
+
+CCTK_FILEVERSION(CactusEinstein_ADMAnalysis_ParamCheck_c)
+
+/********************************************************************
+ *********************     Local Data Types   ***********************
+ ********************************************************************/
+
+/********************************************************************
+ ********************* Local Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ ***************** Scheduled Routine Prototypes *********************
+ ********************************************************************/
+
+
+/********************************************************************
+ ********************* Other Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     Local Data   *****************************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     External Routines   **********************
+ ********************************************************************/
+
+ /*@@
+   @routine    ADMAnalysis_ParamCheck
+   @date       Thu Apr 25 19:04:06 2002
+   @author     Tom Goodale
+   @desc 
+   Scheduled routine to detect invalid parameter settings.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+ 
+ @@*/
+void ADMAnalysis_ParamCheck(CCTK_ARGUMENTS)
+{
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+
+  if(! CCTK_EQUALS(metric_type, "physical") &&
+     ! CCTK_EQUALS(metric_type, "static conformal"))
+  {
+    CCTK_PARAMWARN("Unknown ADMBase::metric_type - known types are \"physical\" and \"static conformal\"");
+  }
+}
+
+/********************************************************************
+ *********************     Local Routines   *************************
+ ********************************************************************/
diff --git a/src/Trace.c b/src/Trace.c
new file mode 100644
index 0000000..241900b
--- /dev/null
+++ b/src/Trace.c
@@ -0,0 +1,180 @@
+ /*@@
+   @file      Trace.c
+   @date      Thu Apr 25 16:35:04 2002
+   @author    Tom Goodale
+   @desc 
+   Calculates the trace of things
+   @enddesc 
+   @version $Header$
+ @@*/
+
+#include "cctk.h"
+
+#include "ADMAnalysis.h"
+
+static const char *rcsid = "$Header$";
+
+CCTK_FILEVERSION(CactusEinstein_ADMAnalysis_Trace_c)
+
+/********************************************************************
+ *********************     Local Data Types   ***********************
+ ********************************************************************/
+
+#define SQR(a) ((a)*(a))
+
+/********************************************************************
+ ********************* Local Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ ***************** Scheduled Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ ********************* Other Routine Prototypes *********************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     Local Data   *****************************
+ ********************************************************************/
+
+/********************************************************************
+ *********************     External Routines   **********************
+ ********************************************************************/
+
+ /*@@
+   @routine    ADMAnalysis_Trace
+   @date       Thu Apr 25 16:49:57 2002
+   @author     Tom Goodale
+   @desc 
+   Calculates the trace of a symmetric 2 index 3-tensor.
+   Optionally returns the trace of the metric at each point.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+   @hdate Thu Apr 25 16:50:38 2002 @hauthor Tom Goodale
+   @hdesc Took old evlatrK routine and made general 
+   @endhistory 
+    @var     lsh
+   @vdesc   grid size
+   @vtype   const int *
+   @vio     in
+   @var     g11
+   @vdesc   the 11 component of the metric tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     g12
+   @vdesc   the 12 component of the metric tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     g13
+   @vdesc   the 13 component of the metric tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     g22
+   @vdesc   the 22 component of the metric tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     g23
+   @vdesc   the 23 component of the metric tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     g33
+   @vdesc   the 33 component of the metric tensor
+   @vtype   const CCTK_REAL *
+   @vio     in
+   @var     tensor11
+   @vdesc   the 11 component of the tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @var     tensor12
+   @vdesc   the 12 component of the tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @var     tensor13
+   @vdesc   the 13 component of the tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @var     tensor22
+   @vdesc   the 22 component of the tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @var     tensor23
+   @vdesc   the 23 component of the tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @var     tensor33
+   @vdesc   the 33 component of the tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @var     detg
+   @vdesc   the determinant of the metric tensor
+   @vtype   CCTK_REAL *
+   @vio     out
+   @vcomment 
+   Will be ignored if NULL.
+   @endvar 
+
+ @@*/
+void ADMAnalysis_Trace(const CCTK_INT *lsh,
+                       const CCTK_REAL *g11,
+                       const CCTK_REAL *g12,
+                       const CCTK_REAL *g13,
+                       const CCTK_REAL *g22,
+                       const CCTK_REAL *g23,
+                       const CCTK_REAL *g33,
+                       CCTK_REAL *tensor11,
+                       CCTK_REAL *tensor12,
+                       CCTK_REAL *tensor13,
+                       CCTK_REAL *tensor22,
+                       CCTK_REAL *tensor23,
+                       CCTK_REAL *tensor33,
+                       CCTK_REAL *trace,
+                       CCTK_REAL *detg)
+{
+  
+  int  i;
+  CCTK_REAL det, u11, u12, u22, u13, u23, u33, two;
+  CCTK_REAL g11p, g12p, g13p, g22p, g23p, g33p;
+  
+  two = 2.0;
+
+  /* loop over all the gridpoints */
+  for(i = 0; i< lsh[0]*lsh[1]*lsh[2];i++)  
+  {
+       
+    /* get the metric */
+    g11p = g11[i];
+    g12p = g12[i];
+    g13p = g13[i];
+    g22p = g22[i];
+    g23p = g23[i];
+    g33p = g33[i];
+      
+    /* compute determinant */
+    det=-(g13p*g13p*g22p)+2.*g12p*g13p*g23p-g11p*g23p*
+      g23p-g12p*g12p*g33p+g11p*g22p*g33p;
+    
+    /* invert metric. This is the conformal upper metric */
+    u11=(-SQR(g23p) + g22p*g33p)/det;
+    u12=(g13p*g23p - g12p*g33p)/det;
+    u22=(-SQR(g13p) + g11p*g33p)/det;
+    u13=(-g13p*g22p + g12p*g23p)/det;
+    u23=(g12p*g13p - g11p*g23p)/det;
+    u33=(-SQR(g12p) + g11p*g22p)/det;
+    
+    /* Calculate trK */
+    trace[i] = (u11*tensor11[i] + u22*tensor22[i] +
+                u33*tensor33[i]+ two*u12*tensor12[i] +
+                two*u13*tensor13[i] + two*u23*tensor23[i]);
+    if(detg)
+    {
+      detg[i]= det;
+    }
+  }
+}
+/********************************************************************
+ *********************     Local Routines   *************************
+ ********************************************************************/
+
diff --git a/src/make.code.defn b/src/make.code.defn
new file mode 100644
index 0000000..6b92303
--- /dev/null
+++ b/src/make.code.defn
@@ -0,0 +1,9 @@
+# Main make.code.defn file for thorn ADMAnalysis
+# $Header$
+
+# Source files in this directory
+SRCS = ParamCheck.c Analysis.c CartToSphere.c Trace.c
+
+# Subdirectories containing source files
+SUBDIRS = 
+