aboutsummaryrefslogtreecommitdiff
path: root/ML_ADM
diff options
context:
space:
mode:
authorErik Schnetter <eschnett@redshift.lan>2007-11-04 16:02:46 -0600
committerErik Schnetter <eschnett@redshift.lan>2007-11-04 16:02:46 -0600
commit257a60a12cd1dfc1c7d46e24b06ba2f320bc426e (patch)
treee2ccac629c4e77df3a0a816fedb862f647d94284 /ML_ADM
parenta48e299f483c7edc4de972a7ecdce4bac9a501f2 (diff)
Add ML_ prefix.
Add README and COPYRIGHT statements.
Diffstat (limited to 'ML_ADM')
-rw-r--r--ML_ADM/configuration.ccl6
-rw-r--r--ML_ADM/interface.ccl117
-rw-r--r--ML_ADM/param.ccl710
-rw-r--r--ML_ADM/schedule.ccl119
-rw-r--r--ML_ADM/src/ADM_Minkowski.c180
-rw-r--r--ML_ADM/src/ADM_RHS.c536
-rw-r--r--ML_ADM/src/ADM_constraints.c445
-rw-r--r--ML_ADM/src/ADM_convertFromADMBase.c212
-rw-r--r--ML_ADM/src/ADM_convertToADMBase.c228
-rw-r--r--ML_ADM/src/Boundaries.c1044
-rw-r--r--ML_ADM/src/Differencing.h12
-rw-r--r--ML_ADM/src/ML_ADM_Minkowski.c180
-rw-r--r--ML_ADM/src/ML_ADM_RHS.c536
-rw-r--r--ML_ADM/src/ML_ADM_constraints.c445
-rw-r--r--ML_ADM/src/ML_ADM_convertFromADMBase.c212
-rw-r--r--ML_ADM/src/ML_ADM_convertToADMBase.c228
-rw-r--r--ML_ADM/src/RegisterMoL.c35
-rw-r--r--ML_ADM/src/RegisterSymmetries.c122
-rw-r--r--ML_ADM/src/Startup.c13
-rw-r--r--ML_ADM/src/make.code.defn6
20 files changed, 5386 insertions, 0 deletions
diff --git a/ML_ADM/configuration.ccl b/ML_ADM/configuration.ccl
new file mode 100644
index 0000000..0374af4
--- /dev/null
+++ b/ML_ADM/configuration.ccl
@@ -0,0 +1,6 @@
+# File produced by user eschnett
+# Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007)
+
+# Mathematica script written by Ian Hinder and Sascha Husa
+
+REQUIRES LoopControl
diff --git a/ML_ADM/interface.ccl b/ML_ADM/interface.ccl
new file mode 100644
index 0000000..c3fae86
--- /dev/null
+++ b/ML_ADM/interface.ccl
@@ -0,0 +1,117 @@
+# File produced by user eschnett
+# Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007)
+
+# Mathematica script written by Ian Hinder and Sascha Husa
+
+implements: ML_ADM
+
+inherits: ADMBase Grid GenericFD Boundary
+
+
+
+USES INCLUDE: GenericFD.h
+USES INCLUDE: Symmetry.h
+USES INCLUDE: sbp_calc_coeffs.h
+USES INCLUDE: Boundary.h
+USES INCLUDE: loopcontrol.h
+
+CCTK_INT FUNCTION MoLRegisterEvolved(CCTK_INT IN EvolvedIndex, CCTK_INT IN RHSIndex)
+USES FUNCTION MoLRegisterEvolved
+
+CCTK_INT FUNCTION MoLRegisterConstrained(CCTK_INT IN ConstrainedIndex)
+USES FUNCTION MoLRegisterConstrained
+
+SUBROUTINE Diff_coeff(CCTK_POINTER_TO_CONST IN cctkGH, CCTK_INT IN dir, CCTK_INT IN nsize, CCTK_INT OUT ARRAY imin, CCTK_INT OUT ARRAY imax, CCTK_REAL OUT ARRAY q, CCTK_INT IN table_handle)
+USES FUNCTION Diff_coeff
+
+CCTK_INT FUNCTION Boundary_SelectGroupForBC(CCTK_POINTER_TO_CONST IN GH, CCTK_INT IN faces, CCTK_INT IN boundary_width, CCTK_INT IN table_handle, CCTK_STRING IN group_name, CCTK_STRING IN bc_name)
+USES FUNCTION Boundary_SelectGroupForBC
+
+CCTK_INT FUNCTION Boundary_SelectVarForBC(CCTK_POINTER_TO_CONST IN GH, CCTK_INT IN faces, CCTK_INT IN boundary_width, CCTK_INT IN table_handle, CCTK_STRING IN var_name, CCTK_STRING IN bc_name)
+USES FUNCTION Boundary_SelectVarForBC
+
+public:
+CCTK_REAL Ham type=GF timelevels=1 tags='tensortypealias="Scalar" tensorweight=1.0000000000000000000'
+{
+ H
+} "Ham"
+
+public:
+CCTK_REAL mom type=GF timelevels=1 tags='tensortypealias="D" tensorweight=1.0000000000000000000'
+{
+ M1,
+ M2,
+ M3
+} "mom"
+
+public:
+CCTK_REAL curvrhs type=GF timelevels=1 tags='tensortypealias="DD_sym" tensorweight=1.0000000000000000000'
+{
+ K11rhs,
+ K21rhs,
+ K31rhs,
+ K22rhs,
+ K32rhs,
+ K33rhs
+} "curvrhs"
+
+public:
+CCTK_REAL lapserhs type=GF timelevels=1 tags='tensortypealias="Scalar" tensorweight=1.0000000000000000000'
+{
+ alpharhs
+} "lapserhs"
+
+public:
+CCTK_REAL metricrhs type=GF timelevels=1 tags='tensortypealias="DD_sym" tensorweight=1.0000000000000000000'
+{
+ g11rhs,
+ g21rhs,
+ g31rhs,
+ g22rhs,
+ g32rhs,
+ g33rhs
+} "metricrhs"
+
+public:
+CCTK_REAL shiftrhs type=GF timelevels=1 tags='tensortypealias="U" tensorweight=1.0000000000000000000'
+{
+ beta1rhs,
+ beta2rhs,
+ beta3rhs
+} "shiftrhs"
+
+public:
+CCTK_REAL curv type=GF timelevels=3 tags='tensortypealias="DD_sym" tensorweight=1.0000000000000000000'
+{
+ K11,
+ K21,
+ K31,
+ K22,
+ K32,
+ K33
+} "curv"
+
+public:
+CCTK_REAL lapse type=GF timelevels=3 tags='tensortypealias="Scalar" tensorweight=1.0000000000000000000'
+{
+ alpha
+} "lapse"
+
+public:
+CCTK_REAL metric type=GF timelevels=3 tags='tensortypealias="DD_sym" tensorweight=1.0000000000000000000'
+{
+ g11,
+ g21,
+ g31,
+ g22,
+ g32,
+ g33
+} "metric"
+
+public:
+CCTK_REAL shift type=GF timelevels=3 tags='tensortypealias="U" tensorweight=1.0000000000000000000'
+{
+ beta1,
+ beta2,
+ beta3
+} "shift"
diff --git a/ML_ADM/param.ccl b/ML_ADM/param.ccl
new file mode 100644
index 0000000..71c2edb
--- /dev/null
+++ b/ML_ADM/param.ccl
@@ -0,0 +1,710 @@
+# File produced by user eschnett
+# Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007)
+
+# Mathematica script written by Ian Hinder and Sascha Husa
+
+
+shares: ADMBase
+
+USES CCTK_KEYWORD initial_data
+
+
+shares: GenericFD
+
+USES CCTK_INT stencil_width
+USES CCTK_INT stencil_width_x
+USES CCTK_INT stencil_width_y
+USES CCTK_INT stencil_width_z
+USES CCTK_INT boundary_width
+
+
+shares: MethodOfLines
+
+USES CCTK_INT MoL_Num_Evolved_Vars
+USES CCTK_INT MoL_Num_Constrained_Vars
+
+restricted:
+CCTK_REAL verbose "verbose"
+{
+ "*:*" :: "no restrictions"
+} 0
+
+private:
+KEYWORD my_initial_data "my_initial_data"
+{
+ "ADMBase" :: "ADMBase"
+ "Minkowski" :: "Minkowski"
+} "ADMBase"
+
+restricted:
+CCTK_INT ML_ADM_MaxNumEvolvedVars "Number of evolved variables used by this thorn" ACCUMULATOR-BASE=MethodofLines::MoL_Num_Evolved_Vars
+{
+ 16:16 :: "Number of evolved variables used by this thorn"
+} 16
+
+restricted:
+CCTK_INT ML_ADM_MaxNumConstrainedVars "Number of constrained variables used by this thorn" ACCUMULATOR-BASE=MethodofLines::MoL_Num_Constrained_Vars
+{
+ 28:28 :: "Number of constrained variables used by this thorn"
+} 28
+
+restricted:
+CCTK_INT ML_ADM_Minkowski_calc_every "ML_ADM_Minkowski_calc_every"
+{
+ *:* :: "no restrictions"
+} 1
+
+restricted:
+CCTK_INT ML_ADM_convertFromADMBase_calc_every "ML_ADM_convertFromADMBase_calc_every"
+{
+ *:* :: "no restrictions"
+} 1
+
+restricted:
+CCTK_INT ML_ADM_RHS_calc_every "ML_ADM_RHS_calc_every"
+{
+ *:* :: "no restrictions"
+} 1
+
+restricted:
+CCTK_INT ML_ADM_convertToADMBase_calc_every "ML_ADM_convertToADMBase_calc_every"
+{
+ *:* :: "no restrictions"
+} 1
+
+restricted:
+CCTK_INT ML_ADM_constraints_calc_every "ML_ADM_constraints_calc_every"
+{
+ *:* :: "no restrictions"
+} 1
+
+restricted:
+CCTK_INT ML_ADM_Minkowski_calc_offset "ML_ADM_Minkowski_calc_offset"
+{
+ *:* :: "no restrictions"
+} 0
+
+restricted:
+CCTK_INT ML_ADM_convertFromADMBase_calc_offset "ML_ADM_convertFromADMBase_calc_offset"
+{
+ *:* :: "no restrictions"
+} 0
+
+restricted:
+CCTK_INT ML_ADM_RHS_calc_offset "ML_ADM_RHS_calc_offset"
+{
+ *:* :: "no restrictions"
+} 0
+
+restricted:
+CCTK_INT ML_ADM_convertToADMBase_calc_offset "ML_ADM_convertToADMBase_calc_offset"
+{
+ *:* :: "no restrictions"
+} 0
+
+restricted:
+CCTK_INT ML_ADM_constraints_calc_offset "ML_ADM_constraints_calc_offset"
+{
+ *:* :: "no restrictions"
+} 0
+
+private:
+KEYWORD K11_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD K21_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD K31_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD K22_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD K32_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD K33_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD alpha_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD g11_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD g21_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD g31_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD g22_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD g32_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD g33_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD beta1_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD beta2_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD beta3_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD curv_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD lapse_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD metric_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+KEYWORD shift_bound "Boundary condition to implement"
+{
+ "flat" :: "Flat boundary condition"
+ "none" :: "No boundary condition"
+ "static" :: "Boundaries held fixed"
+ "radiative" :: "Radiation boundary condition"
+ "scalar" :: "Dirichlet boundary condition"
+ "newrad" :: "Improved radiative boundary condition"
+ "skip" :: "skip boundary condition code"
+} "skip"
+
+private:
+CCTK_REAL K11_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL K21_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL K31_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL K22_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL K32_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL K33_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL alpha_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL g11_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL g21_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL g31_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL g22_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL g32_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL g33_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL beta1_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL beta2_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL beta3_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL curv_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL lapse_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL metric_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL shift_bound_speed "characteristic speed at boundary"
+{
+ "0:*" :: "outgoing characteristic speed > 0"
+} 1.
+
+private:
+CCTK_REAL K11_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL K21_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL K31_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL K22_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL K32_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL K33_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL alpha_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL g11_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL g21_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL g31_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL g22_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL g32_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL g33_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL beta1_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL beta2_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL beta3_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL curv_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL lapse_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL metric_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL shift_bound_limit "limit value for r -> infinity"
+{
+ "*:*" :: "value of limit value is unrestricted"
+} 0.
+
+private:
+CCTK_REAL K11_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL K21_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL K31_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL K22_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL K32_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL K33_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL alpha_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL g11_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL g21_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL g31_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL g22_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL g32_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL g33_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL beta1_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL beta2_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL beta3_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL curv_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL lapse_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL metric_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
+private:
+CCTK_REAL shift_bound_scalar "Dirichlet boundary value"
+{
+ "*:*" :: "unrestricted"
+} 0.
+
diff --git a/ML_ADM/schedule.ccl b/ML_ADM/schedule.ccl
new file mode 100644
index 0000000..418de03
--- /dev/null
+++ b/ML_ADM/schedule.ccl
@@ -0,0 +1,119 @@
+# File produced by user eschnett
+# Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007)
+
+# Mathematica script written by Ian Hinder and Sascha Husa
+
+
+STORAGE: Ham[1]
+
+STORAGE: mom[1]
+
+STORAGE: curvrhs[1]
+
+STORAGE: lapserhs[1]
+
+STORAGE: metricrhs[1]
+
+STORAGE: shiftrhs[1]
+
+STORAGE: curv[3]
+
+STORAGE: lapse[3]
+
+STORAGE: metric[3]
+
+STORAGE: shift[3]
+
+schedule ML_ADM_Startup at STARTUP
+{
+ LANG: C
+ OPTIONS: meta
+} "create banner"
+
+schedule ML_ADM_RegisterVars in MoL_Register
+{
+ LANG: C
+ OPTIONS: meta
+} "Register Variables for MoL"
+
+schedule ML_ADM_RegisterSymmetries at BASEGRID
+{
+ LANG: C
+ OPTIONS: meta
+} "register symmetries"
+
+
+if (CCTK_EQUALS(my_initial_data, "Minkowski"))
+{
+ schedule ML_ADM_Minkowski IN ADMBase_InitialData
+ {
+ LANG: C
+
+ } "ML_ADM_Minkowski"
+}
+
+
+if (CCTK_EQUALS(my_initial_data, "ADMBase"))
+{
+ schedule ML_ADM_convertFromADMBase AT initial
+ {
+ LANG: C
+
+ } "ML_ADM_convertFromADMBase"
+}
+
+
+if (CCTK_EQUALS(my_initial_data, "ADMBase"))
+{
+ schedule ML_ADM_convertFromADMBase AFTER ADMBase_PostInitial
+ {
+ LANG: C
+
+ } "ML_ADM_convertFromADMBase"
+}
+
+schedule ML_ADM_RHS IN MoL_CalcRHS
+{
+ LANG: C
+
+} "ML_ADM_RHS"
+
+schedule ML_ADM_convertToADMBase IN MoL_PostStep
+{
+ LANG: C
+
+} "ML_ADM_convertToADMBase"
+
+schedule ML_ADM_convertToADMBase AFTER ADM_ApplyBoundConds
+{
+ LANG: C
+
+} "ML_ADM_convertToADMBase"
+
+schedule ML_ADM_constraints AT analysis
+{
+ LANG: C
+
+} "ML_ADM_constraints"
+
+schedule ML_ADM_ApplyBoundConds in MoL_PostStep
+{
+ LANG: C
+ OPTIONS: level
+ SYNC: curv
+ SYNC: lapse
+ SYNC: metric
+ SYNC: shift
+} "apply boundary conditions"
+
+schedule ML_ADM_CheckBoundaries at BASEGRID
+{
+ LANG: C
+ OPTIONS: meta
+} "check boundaries treatment"
+
+schedule group ApplyBCs as ML_ADM_ApplyBCs in MoL_PostStep after ML_ADM_ApplyBoundConds
+{
+ # no language specified
+
+} "Apply boundary conditions controlled by thorn Boundary"
diff --git a/ML_ADM/src/ADM_Minkowski.c b/ML_ADM/src/ADM_Minkowski.c
new file mode 100644
index 0000000..c21c98e
--- /dev/null
+++ b/ML_ADM/src/ADM_Minkowski.c
@@ -0,0 +1,180 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ADM_Minkowski_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ADM_Minkowski_Body");
+ }
+
+ if (cctk_iteration % ADM_Minkowski_calc_every != ADM_Minkowski_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alphaL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+
+ /* Assign local copies of grid functions */
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ g11L = 1;
+
+ g21L = 0;
+
+ g31L = 0;
+
+ g22L = 1;
+
+ g32L = 0;
+
+ g33L = 1;
+
+ K11L = 0;
+
+ K21L = 0;
+
+ K31L = 0;
+
+ K22L = 0;
+
+ K32L = 0;
+
+ K33L = 0;
+
+ alphaL = 1;
+
+ beta1L = 0;
+
+ beta2L = 0;
+
+ beta3L = 0;
+
+
+ /* Copy local copies back to grid functions */
+ alpha[index] = alphaL;
+ beta1[index] = beta1L;
+ beta2[index] = beta2L;
+ beta3[index] = beta3L;
+ g11[index] = g11L;
+ g21[index] = g21L;
+ g22[index] = g22L;
+ g31[index] = g31L;
+ g32[index] = g32L;
+ g33[index] = g33L;
+ K11[index] = K11L;
+ K21[index] = K21L;
+ K22[index] = K22L;
+ K31[index] = K31L;
+ K32[index] = K32L;
+ K33[index] = K33L;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ADM_Minkowski(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ADM_Minkowski_Body);
+}
diff --git a/ML_ADM/src/ADM_RHS.c b/ML_ADM/src/ADM_RHS.c
new file mode 100644
index 0000000..1a1a7ed
--- /dev/null
+++ b/ML_ADM/src/ADM_RHS.c
@@ -0,0 +1,536 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ADM_RHS_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ADM_RHS_Body");
+ }
+
+ if (cctk_iteration % ADM_RHS_calc_every != ADM_RHS_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+ CCTK_REAL detg = INITVALUE;
+ CCTK_REAL G111 = INITVALUE, G121 = INITVALUE, G122 = INITVALUE, G131 = INITVALUE, G132 = INITVALUE, G133 = INITVALUE;
+ CCTK_REAL G211 = INITVALUE, G221 = INITVALUE, G222 = INITVALUE, G231 = INITVALUE, G232 = INITVALUE, G233 = INITVALUE;
+ CCTK_REAL G311 = INITVALUE, G321 = INITVALUE, G322 = INITVALUE, G331 = INITVALUE, G332 = INITVALUE, G333 = INITVALUE;
+ CCTK_REAL gu11 = INITVALUE, gu12 = INITVALUE, gu13 = INITVALUE, gu22 = INITVALUE, gu23 = INITVALUE, gu33 = INITVALUE;
+ CCTK_REAL Km11 = INITVALUE, Km12 = INITVALUE, Km13 = INITVALUE, Km21 = INITVALUE, Km22 = INITVALUE, Km23 = INITVALUE;
+ CCTK_REAL Km31 = INITVALUE, Km32 = INITVALUE, Km33 = INITVALUE;
+ CCTK_REAL R11 = INITVALUE, R21 = INITVALUE, R22 = INITVALUE, R31 = INITVALUE, R32 = INITVALUE, R33 = INITVALUE;
+ CCTK_REAL trK = INITVALUE;
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alphaL = INITVALUE, alpharhsL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta1rhsL = INITVALUE, beta2L = INITVALUE, beta2rhsL = INITVALUE, beta3L = INITVALUE, beta3rhsL = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g11rhsL = INITVALUE, g21L = INITVALUE, g21rhsL = INITVALUE, g22L = INITVALUE, g22rhsL = INITVALUE;
+ CCTK_REAL g31L = INITVALUE, g31rhsL = INITVALUE, g32L = INITVALUE, g32rhsL = INITVALUE, g33L = INITVALUE, g33rhsL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K11rhsL = INITVALUE, K21L = INITVALUE, K21rhsL = INITVALUE, K22L = INITVALUE, K22rhsL = INITVALUE;
+ CCTK_REAL K31L = INITVALUE, K31rhsL = INITVALUE, K32L = INITVALUE, K32rhsL = INITVALUE, K33L = INITVALUE, K33rhsL = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+ CCTK_REAL PDstandard4th1alpha = INITVALUE;
+ CCTK_REAL PDstandard4th2alpha = INITVALUE;
+ CCTK_REAL PDstandard4th3alpha = INITVALUE;
+ CCTK_REAL PDstandard4th11alpha = INITVALUE;
+ CCTK_REAL PDstandard4th22alpha = INITVALUE;
+ CCTK_REAL PDstandard4th33alpha = INITVALUE;
+ CCTK_REAL PDstandard4th12alpha = INITVALUE;
+ CCTK_REAL PDstandard4th13alpha = INITVALUE;
+ CCTK_REAL PDstandard4th23alpha = INITVALUE;
+ CCTK_REAL PDstandard4th1beta1 = INITVALUE;
+ CCTK_REAL PDstandard4th2beta1 = INITVALUE;
+ CCTK_REAL PDstandard4th3beta1 = INITVALUE;
+ CCTK_REAL PDstandard4th1beta2 = INITVALUE;
+ CCTK_REAL PDstandard4th2beta2 = INITVALUE;
+ CCTK_REAL PDstandard4th3beta2 = INITVALUE;
+ CCTK_REAL PDstandard4th1beta3 = INITVALUE;
+ CCTK_REAL PDstandard4th2beta3 = INITVALUE;
+ CCTK_REAL PDstandard4th3beta3 = INITVALUE;
+ CCTK_REAL PDstandard4th1g11 = INITVALUE;
+ CCTK_REAL PDstandard4th2g11 = INITVALUE;
+ CCTK_REAL PDstandard4th3g11 = INITVALUE;
+ CCTK_REAL PDstandard4th22g11 = INITVALUE;
+ CCTK_REAL PDstandard4th33g11 = INITVALUE;
+ CCTK_REAL PDstandard4th12g11 = INITVALUE;
+ CCTK_REAL PDstandard4th13g11 = INITVALUE;
+ CCTK_REAL PDstandard4th21g11 = INITVALUE;
+ CCTK_REAL PDstandard4th23g11 = INITVALUE;
+ CCTK_REAL PDstandard4th31g11 = INITVALUE;
+ CCTK_REAL PDstandard4th32g11 = INITVALUE;
+ CCTK_REAL PDstandard4th1g21 = INITVALUE;
+ CCTK_REAL PDstandard4th2g21 = INITVALUE;
+ CCTK_REAL PDstandard4th3g21 = INITVALUE;
+ CCTK_REAL PDstandard4th33g21 = INITVALUE;
+ CCTK_REAL PDstandard4th12g21 = INITVALUE;
+ CCTK_REAL PDstandard4th13g21 = INITVALUE;
+ CCTK_REAL PDstandard4th21g21 = INITVALUE;
+ CCTK_REAL PDstandard4th23g21 = INITVALUE;
+ CCTK_REAL PDstandard4th31g21 = INITVALUE;
+ CCTK_REAL PDstandard4th32g21 = INITVALUE;
+ CCTK_REAL PDstandard4th1g22 = INITVALUE;
+ CCTK_REAL PDstandard4th2g22 = INITVALUE;
+ CCTK_REAL PDstandard4th3g22 = INITVALUE;
+ CCTK_REAL PDstandard4th11g22 = INITVALUE;
+ CCTK_REAL PDstandard4th33g22 = INITVALUE;
+ CCTK_REAL PDstandard4th12g22 = INITVALUE;
+ CCTK_REAL PDstandard4th13g22 = INITVALUE;
+ CCTK_REAL PDstandard4th21g22 = INITVALUE;
+ CCTK_REAL PDstandard4th23g22 = INITVALUE;
+ CCTK_REAL PDstandard4th31g22 = INITVALUE;
+ CCTK_REAL PDstandard4th32g22 = INITVALUE;
+ CCTK_REAL PDstandard4th1g31 = INITVALUE;
+ CCTK_REAL PDstandard4th2g31 = INITVALUE;
+ CCTK_REAL PDstandard4th3g31 = INITVALUE;
+ CCTK_REAL PDstandard4th22g31 = INITVALUE;
+ CCTK_REAL PDstandard4th12g31 = INITVALUE;
+ CCTK_REAL PDstandard4th13g31 = INITVALUE;
+ CCTK_REAL PDstandard4th21g31 = INITVALUE;
+ CCTK_REAL PDstandard4th23g31 = INITVALUE;
+ CCTK_REAL PDstandard4th31g31 = INITVALUE;
+ CCTK_REAL PDstandard4th32g31 = INITVALUE;
+ CCTK_REAL PDstandard4th1g32 = INITVALUE;
+ CCTK_REAL PDstandard4th2g32 = INITVALUE;
+ CCTK_REAL PDstandard4th3g32 = INITVALUE;
+ CCTK_REAL PDstandard4th11g32 = INITVALUE;
+ CCTK_REAL PDstandard4th12g32 = INITVALUE;
+ CCTK_REAL PDstandard4th13g32 = INITVALUE;
+ CCTK_REAL PDstandard4th21g32 = INITVALUE;
+ CCTK_REAL PDstandard4th23g32 = INITVALUE;
+ CCTK_REAL PDstandard4th31g32 = INITVALUE;
+ CCTK_REAL PDstandard4th32g32 = INITVALUE;
+ CCTK_REAL PDstandard4th1g33 = INITVALUE;
+ CCTK_REAL PDstandard4th2g33 = INITVALUE;
+ CCTK_REAL PDstandard4th3g33 = INITVALUE;
+ CCTK_REAL PDstandard4th11g33 = INITVALUE;
+ CCTK_REAL PDstandard4th22g33 = INITVALUE;
+ CCTK_REAL PDstandard4th12g33 = INITVALUE;
+ CCTK_REAL PDstandard4th13g33 = INITVALUE;
+ CCTK_REAL PDstandard4th21g33 = INITVALUE;
+ CCTK_REAL PDstandard4th23g33 = INITVALUE;
+ CCTK_REAL PDstandard4th31g33 = INITVALUE;
+ CCTK_REAL PDstandard4th32g33 = INITVALUE;
+ CCTK_REAL PDstandard4th1K11 = INITVALUE;
+ CCTK_REAL PDstandard4th2K11 = INITVALUE;
+ CCTK_REAL PDstandard4th3K11 = INITVALUE;
+ CCTK_REAL PDstandard4th1K21 = INITVALUE;
+ CCTK_REAL PDstandard4th2K21 = INITVALUE;
+ CCTK_REAL PDstandard4th3K21 = INITVALUE;
+ CCTK_REAL PDstandard4th1K22 = INITVALUE;
+ CCTK_REAL PDstandard4th2K22 = INITVALUE;
+ CCTK_REAL PDstandard4th3K22 = INITVALUE;
+ CCTK_REAL PDstandard4th1K31 = INITVALUE;
+ CCTK_REAL PDstandard4th2K31 = INITVALUE;
+ CCTK_REAL PDstandard4th3K31 = INITVALUE;
+ CCTK_REAL PDstandard4th1K32 = INITVALUE;
+ CCTK_REAL PDstandard4th2K32 = INITVALUE;
+ CCTK_REAL PDstandard4th3K32 = INITVALUE;
+ CCTK_REAL PDstandard4th1K33 = INITVALUE;
+ CCTK_REAL PDstandard4th2K33 = INITVALUE;
+ CCTK_REAL PDstandard4th3K33 = INITVALUE;
+
+ /* Assign local copies of grid functions */
+ alphaL = alpha[index];
+ beta1L = beta1[index];
+ beta2L = beta2[index];
+ beta3L = beta3[index];
+ g11L = g11[index];
+ g21L = g21[index];
+ g22L = g22[index];
+ g31L = g31[index];
+ g32L = g32[index];
+ g33L = g33[index];
+ K11L = K11[index];
+ K21L = K21[index];
+ K22L = K22[index];
+ K31L = K31[index];
+ K32L = K32[index];
+ K33L = K33[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+ PDstandard4th1alpha = PDstandard4th1(alpha, i, j, k);
+ PDstandard4th2alpha = PDstandard4th2(alpha, i, j, k);
+ PDstandard4th3alpha = PDstandard4th3(alpha, i, j, k);
+ PDstandard4th11alpha = PDstandard4th11(alpha, i, j, k);
+ PDstandard4th22alpha = PDstandard4th22(alpha, i, j, k);
+ PDstandard4th33alpha = PDstandard4th33(alpha, i, j, k);
+ PDstandard4th12alpha = PDstandard4th12(alpha, i, j, k);
+ PDstandard4th13alpha = PDstandard4th13(alpha, i, j, k);
+ PDstandard4th23alpha = PDstandard4th23(alpha, i, j, k);
+ PDstandard4th1beta1 = PDstandard4th1(beta1, i, j, k);
+ PDstandard4th2beta1 = PDstandard4th2(beta1, i, j, k);
+ PDstandard4th3beta1 = PDstandard4th3(beta1, i, j, k);
+ PDstandard4th1beta2 = PDstandard4th1(beta2, i, j, k);
+ PDstandard4th2beta2 = PDstandard4th2(beta2, i, j, k);
+ PDstandard4th3beta2 = PDstandard4th3(beta2, i, j, k);
+ PDstandard4th1beta3 = PDstandard4th1(beta3, i, j, k);
+ PDstandard4th2beta3 = PDstandard4th2(beta3, i, j, k);
+ PDstandard4th3beta3 = PDstandard4th3(beta3, i, j, k);
+ PDstandard4th1g11 = PDstandard4th1(g11, i, j, k);
+ PDstandard4th2g11 = PDstandard4th2(g11, i, j, k);
+ PDstandard4th3g11 = PDstandard4th3(g11, i, j, k);
+ PDstandard4th22g11 = PDstandard4th22(g11, i, j, k);
+ PDstandard4th33g11 = PDstandard4th33(g11, i, j, k);
+ PDstandard4th23g11 = PDstandard4th23(g11, i, j, k);
+ PDstandard4th1g21 = PDstandard4th1(g21, i, j, k);
+ PDstandard4th2g21 = PDstandard4th2(g21, i, j, k);
+ PDstandard4th3g21 = PDstandard4th3(g21, i, j, k);
+ PDstandard4th33g21 = PDstandard4th33(g21, i, j, k);
+ PDstandard4th12g21 = PDstandard4th12(g21, i, j, k);
+ PDstandard4th13g21 = PDstandard4th13(g21, i, j, k);
+ PDstandard4th23g21 = PDstandard4th23(g21, i, j, k);
+ PDstandard4th1g22 = PDstandard4th1(g22, i, j, k);
+ PDstandard4th2g22 = PDstandard4th2(g22, i, j, k);
+ PDstandard4th3g22 = PDstandard4th3(g22, i, j, k);
+ PDstandard4th11g22 = PDstandard4th11(g22, i, j, k);
+ PDstandard4th33g22 = PDstandard4th33(g22, i, j, k);
+ PDstandard4th13g22 = PDstandard4th13(g22, i, j, k);
+ PDstandard4th1g31 = PDstandard4th1(g31, i, j, k);
+ PDstandard4th2g31 = PDstandard4th2(g31, i, j, k);
+ PDstandard4th3g31 = PDstandard4th3(g31, i, j, k);
+ PDstandard4th22g31 = PDstandard4th22(g31, i, j, k);
+ PDstandard4th12g31 = PDstandard4th12(g31, i, j, k);
+ PDstandard4th13g31 = PDstandard4th13(g31, i, j, k);
+ PDstandard4th23g31 = PDstandard4th23(g31, i, j, k);
+ PDstandard4th1g32 = PDstandard4th1(g32, i, j, k);
+ PDstandard4th2g32 = PDstandard4th2(g32, i, j, k);
+ PDstandard4th3g32 = PDstandard4th3(g32, i, j, k);
+ PDstandard4th11g32 = PDstandard4th11(g32, i, j, k);
+ PDstandard4th12g32 = PDstandard4th12(g32, i, j, k);
+ PDstandard4th13g32 = PDstandard4th13(g32, i, j, k);
+ PDstandard4th23g32 = PDstandard4th23(g32, i, j, k);
+ PDstandard4th1g33 = PDstandard4th1(g33, i, j, k);
+ PDstandard4th2g33 = PDstandard4th2(g33, i, j, k);
+ PDstandard4th3g33 = PDstandard4th3(g33, i, j, k);
+ PDstandard4th11g33 = PDstandard4th11(g33, i, j, k);
+ PDstandard4th22g33 = PDstandard4th22(g33, i, j, k);
+ PDstandard4th12g33 = PDstandard4th12(g33, i, j, k);
+ PDstandard4th1K11 = PDstandard4th1(K11, i, j, k);
+ PDstandard4th2K11 = PDstandard4th2(K11, i, j, k);
+ PDstandard4th3K11 = PDstandard4th3(K11, i, j, k);
+ PDstandard4th1K21 = PDstandard4th1(K21, i, j, k);
+ PDstandard4th2K21 = PDstandard4th2(K21, i, j, k);
+ PDstandard4th3K21 = PDstandard4th3(K21, i, j, k);
+ PDstandard4th1K22 = PDstandard4th1(K22, i, j, k);
+ PDstandard4th2K22 = PDstandard4th2(K22, i, j, k);
+ PDstandard4th3K22 = PDstandard4th3(K22, i, j, k);
+ PDstandard4th1K31 = PDstandard4th1(K31, i, j, k);
+ PDstandard4th2K31 = PDstandard4th2(K31, i, j, k);
+ PDstandard4th3K31 = PDstandard4th3(K31, i, j, k);
+ PDstandard4th1K32 = PDstandard4th1(K32, i, j, k);
+ PDstandard4th2K32 = PDstandard4th2(K32, i, j, k);
+ PDstandard4th3K32 = PDstandard4th3(K32, i, j, k);
+ PDstandard4th1K33 = PDstandard4th1(K33, i, j, k);
+ PDstandard4th2K33 = PDstandard4th2(K33, i, j, k);
+ PDstandard4th3K33 = PDstandard4th3(K33, i, j, k);
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ detg = 2*g21L*g31L*g32L + g33L*(g11L*g22L - SQR(g21L)) - g22L*SQR(g31L) - g11L*SQR(g32L);
+
+ gu11 = INV(detg)*(g22L*g33L - SQR(g32L));
+
+ gu12 = (g31L*g32L - g21L*g33L)*INV(detg);
+
+ gu13 = (-(g22L*g31L) + g21L*g32L)*INV(detg);
+
+ gu22 = INV(detg)*(g11L*g33L - SQR(g31L));
+
+ gu23 = (g21L*g31L - g11L*g32L)*INV(detg);
+
+ gu33 = INV(detg)*(g11L*g22L - SQR(g21L));
+
+ G111 = khalf*(gu11*PDstandard4th1g11 + 2*(gu12*PDstandard4th1g21 + gu13*PDstandard4th1g31) - gu12*PDstandard4th2g11 -
+ gu13*PDstandard4th3g11);
+
+ G211 = khalf*(gu12*PDstandard4th1g11 + 2*(gu22*PDstandard4th1g21 + gu23*PDstandard4th1g31) - gu22*PDstandard4th2g11 -
+ gu23*PDstandard4th3g11);
+
+ G311 = khalf*(gu13*PDstandard4th1g11 + 2*(gu23*PDstandard4th1g21 + gu33*PDstandard4th1g31) - gu23*PDstandard4th2g11 -
+ gu33*PDstandard4th3g11);
+
+ G121 = khalf*(gu12*PDstandard4th1g22 + gu11*PDstandard4th2g11 +
+ gu13*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G221 = khalf*(gu22*PDstandard4th1g22 + gu12*PDstandard4th2g11 +
+ gu23*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G321 = khalf*(gu23*PDstandard4th1g22 + gu13*PDstandard4th2g11 +
+ gu33*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G131 = khalf*(gu13*PDstandard4th1g33 + gu11*PDstandard4th3g11 +
+ gu12*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G231 = khalf*(gu23*PDstandard4th1g33 + gu12*PDstandard4th3g11 +
+ gu22*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G331 = khalf*(gu33*PDstandard4th1g33 + gu13*PDstandard4th3g11 +
+ gu23*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G122 = khalf*(gu11*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu12*PDstandard4th2g22 +
+ gu13*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G222 = khalf*(gu12*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu22*PDstandard4th2g22 +
+ gu23*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G322 = khalf*(gu13*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu23*PDstandard4th2g22 +
+ gu33*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G132 = khalf*(gu13*PDstandard4th2g33 + gu11*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu12*PDstandard4th3g22);
+
+ G232 = khalf*(gu23*PDstandard4th2g33 + gu12*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu22*PDstandard4th3g22);
+
+ G332 = khalf*(gu33*PDstandard4th2g33 + gu13*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu23*PDstandard4th3g22);
+
+ G133 = khalf*(-(gu11*PDstandard4th1g33) - gu12*PDstandard4th2g33 + 2*gu11*PDstandard4th3g31 +
+ 2*gu12*PDstandard4th3g32 + gu13*PDstandard4th3g33);
+
+ G233 = khalf*(-(gu12*PDstandard4th1g33) - gu22*PDstandard4th2g33 + 2*gu12*PDstandard4th3g31 +
+ 2*gu22*PDstandard4th3g32 + gu23*PDstandard4th3g33);
+
+ G333 = khalf*(-(gu13*PDstandard4th1g33) - gu23*PDstandard4th2g33 + 2*gu13*PDstandard4th3g31 +
+ 2*gu23*PDstandard4th3g32 + gu33*PDstandard4th3g33);
+
+ R11 = khalf*(4*G231*G321 + G211*(2*G121 - 2*G222 - 2*G332) + G311*(2*G131 - 2*G333) - gu22*PDstandard4th11g22 -
+ 2*(G111*G221 + G232*G311 + G111*G331 + gu23*PDstandard4th11g32) - gu33*PDstandard4th11g33 +
+ 2*gu22*PDstandard4th12g21 + 2*gu23*PDstandard4th12g31 + 2*gu23*PDstandard4th13g21 + 2*gu33*PDstandard4th13g31 -
+ gu22*PDstandard4th22g11 - 2*gu23*PDstandard4th23g11 - gu33*PDstandard4th33g11 + 2*SQR(G221) + 2*SQR(G331));
+
+ R21 = khalf*(2*(G122*G211 + G132*G311 + G231*G322 + G331*G332) -
+ 2*(G121*G221 + G121*G331 + G221*G332 + G321*G333 + gu12*PDstandard4th12g21) - gu23*PDstandard4th12g32 -
+ gu33*PDstandard4th12g33 + gu13*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21) +
+ gu23*PDstandard4th13g22 + gu33*PDstandard4th13g32 + gu12*(PDstandard4th11g22 + PDstandard4th22g11) +
+ gu23*PDstandard4th22g31 + gu13*PDstandard4th23g11 - gu23*PDstandard4th23g21 + gu33*PDstandard4th23g31 -
+ gu33*PDstandard4th33g21);
+
+ R31 = khalf*(2*(G132*G211 + G221*G232 + G133*G311 + G233*G321) -
+ 2*(G222*G231 + G232*G331 + G131*(G221 + G331) + gu13*PDstandard4th13g31) +
+ gu12*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21 + PDstandard4th23g11) +
+ gu22*(PDstandard4th12g32 - PDstandard4th13g22 - PDstandard4th22g31 + PDstandard4th23g21) +
+ gu13*(PDstandard4th11g33 + PDstandard4th33g11) +
+ gu23*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21));
+
+ R22 = khalf*(4*G132*G321 + G122*(2*G221 - 2*(G111 + G331)) + G322*(2*G232 - 2*G333) -
+ 2*(G131*G322 + G222*(G121 + G332) + gu13*PDstandard4th13g22) +
+ gu11*(-PDstandard4th11g22 + 2*PDstandard4th12g21 - PDstandard4th22g11) +
+ gu13*(-2*PDstandard4th22g31 + 2*(PDstandard4th12g32 + PDstandard4th23g21)) +
+ gu33*(-PDstandard4th22g33 + 2*PDstandard4th23g32 - PDstandard4th33g22) + 2*(SQR(G121) + SQR(G332)));
+
+ R32 = khalf*(2*(G121*G131 + G122*G231 + G133*G321 + G233*G322) +
+ gu11*(-PDstandard4th11g32 + PDstandard4th12g31 + PDstandard4th13g21 - PDstandard4th23g11) +
+ gu12*(-PDstandard4th12g32 + PDstandard4th13g22 + PDstandard4th22g31 - PDstandard4th23g21) -
+ 2*(G111*G132 + G131*G332 + G232*(G121 + G332) + gu23*PDstandard4th23g32) +
+ gu13*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21) +
+ gu23*(PDstandard4th22g33 + PDstandard4th33g22));
+
+ R33 = khalf*(4*G132*G231 - gu11*PDstandard4th11g33 -
+ 2*(G111*G133 + G133*G221 + G121*G233 + G222*G233 + G131*G333 + G232*G333 + gu12*PDstandard4th12g33) +
+ 2*(G133*G331 + G233*G332 + gu11*PDstandard4th13g31) + 2*gu12*PDstandard4th13g32 - gu22*PDstandard4th22g33 +
+ 2*gu12*PDstandard4th23g31 + 2*gu22*PDstandard4th23g32 - gu11*PDstandard4th33g11 - 2*gu12*PDstandard4th33g21 -
+ gu22*PDstandard4th33g22 + 2*SQR(G131) + 2*SQR(G232));
+
+ Km11 = gu11*K11L + gu12*K21L + gu13*K31L;
+
+ Km21 = gu12*K11L + gu22*K21L + gu23*K31L;
+
+ Km31 = gu13*K11L + gu23*K21L + gu33*K31L;
+
+ Km12 = gu11*K21L + gu12*K22L + gu13*K32L;
+
+ Km22 = gu12*K21L + gu22*K22L + gu23*K32L;
+
+ Km32 = gu13*K21L + gu23*K22L + gu33*K32L;
+
+ Km13 = gu11*K31L + gu12*K32L + gu13*K33L;
+
+ Km23 = gu12*K31L + gu22*K32L + gu23*K33L;
+
+ Km33 = gu13*K31L + gu23*K32L + gu33*K33L;
+
+ trK = Km11 + Km22 + Km33;
+
+ g11rhsL = -2*alphaL*K11L + 2*(g11L*PDstandard4th1beta1 + g21L*PDstandard4th1beta2 + g31L*PDstandard4th1beta3) +
+ beta1L*PDstandard4th1g11 + beta2L*PDstandard4th2g11 + beta3L*PDstandard4th3g11;
+
+ g21rhsL = -2*alphaL*K21L + g22L*PDstandard4th1beta2 + g32L*PDstandard4th1beta3 + beta1L*PDstandard4th1g21 +
+ g11L*PDstandard4th2beta1 + g21L*(PDstandard4th1beta1 + PDstandard4th2beta2) + g31L*PDstandard4th2beta3 +
+ beta2L*PDstandard4th2g21 + beta3L*PDstandard4th3g21;
+
+ g31rhsL = -2*alphaL*K31L + g32L*PDstandard4th1beta2 + g33L*PDstandard4th1beta3 + beta1L*PDstandard4th1g31 +
+ beta2L*PDstandard4th2g31 + g11L*PDstandard4th3beta1 + g21L*PDstandard4th3beta2 +
+ g31L*(PDstandard4th1beta1 + PDstandard4th3beta3) + beta3L*PDstandard4th3g31;
+
+ g22rhsL = -2*alphaL*K22L + beta1L*PDstandard4th1g22 +
+ 2*(g21L*PDstandard4th2beta1 + g22L*PDstandard4th2beta2 + g32L*PDstandard4th2beta3) + beta2L*PDstandard4th2g22 +
+ beta3L*PDstandard4th3g22;
+
+ g32rhsL = -2*alphaL*K32L + beta1L*PDstandard4th1g32 + g31L*PDstandard4th2beta1 + g33L*PDstandard4th2beta3 +
+ beta2L*PDstandard4th2g32 + g21L*PDstandard4th3beta1 + g22L*PDstandard4th3beta2 +
+ g32L*(PDstandard4th2beta2 + PDstandard4th3beta3) + beta3L*PDstandard4th3g32;
+
+ g33rhsL = -2*alphaL*K33L + beta1L*PDstandard4th1g33 + beta2L*PDstandard4th2g33 +
+ 2*(g31L*PDstandard4th3beta1 + g32L*PDstandard4th3beta2 + g33L*PDstandard4th3beta3) + beta3L*PDstandard4th3g33;
+
+ K11rhsL = -PDstandard4th11alpha + G111*PDstandard4th1alpha +
+ 2*(K11L*PDstandard4th1beta1 + K21L*PDstandard4th1beta2 + K31L*PDstandard4th1beta3) + beta1L*PDstandard4th1K11 +
+ G211*PDstandard4th2alpha + beta2L*PDstandard4th2K11 + G311*PDstandard4th3alpha + beta3L*PDstandard4th3K11 +
+ alphaL*(-2*(K11L*Km11 + K21L*Km21 + K31L*Km31) + R11 + K11L*trK);
+
+ K21rhsL = -PDstandard4th12alpha + G121*PDstandard4th1alpha + K22L*PDstandard4th1beta2 + K32L*PDstandard4th1beta3 +
+ beta1L*PDstandard4th1K21 + G221*PDstandard4th2alpha + K11L*PDstandard4th2beta1 +
+ K21L*(PDstandard4th1beta1 + PDstandard4th2beta2) + K31L*PDstandard4th2beta3 + beta2L*PDstandard4th2K21 +
+ G321*PDstandard4th3alpha + beta3L*PDstandard4th3K21 +
+ alphaL*(-2*(K11L*Km12 + K21L*Km22 + K31L*Km32) + R21 + K21L*trK);
+
+ K31rhsL = -PDstandard4th13alpha + G131*PDstandard4th1alpha + K32L*PDstandard4th1beta2 + K33L*PDstandard4th1beta3 +
+ beta1L*PDstandard4th1K31 + G231*PDstandard4th2alpha + beta2L*PDstandard4th2K31 + G331*PDstandard4th3alpha +
+ K11L*PDstandard4th3beta1 + K21L*PDstandard4th3beta2 + K31L*(PDstandard4th1beta1 + PDstandard4th3beta3) +
+ beta3L*PDstandard4th3K31 + alphaL*(-2*(K11L*Km13 + K21L*Km23 + K31L*Km33) + R31 + K31L*trK);
+
+ K22rhsL = G122*PDstandard4th1alpha + beta1L*PDstandard4th1K22 - PDstandard4th22alpha + G222*PDstandard4th2alpha +
+ 2*(K21L*PDstandard4th2beta1 + K22L*PDstandard4th2beta2 + K32L*PDstandard4th2beta3) + beta2L*PDstandard4th2K22 +
+ G322*PDstandard4th3alpha + beta3L*PDstandard4th3K22 +
+ alphaL*(-2*(K21L*Km12 + K22L*Km22 + K32L*Km32) + R22 + K22L*trK);
+
+ K32rhsL = G132*PDstandard4th1alpha + beta1L*PDstandard4th1K32 - PDstandard4th23alpha + G232*PDstandard4th2alpha +
+ K31L*PDstandard4th2beta1 + K33L*PDstandard4th2beta3 + beta2L*PDstandard4th2K32 + G332*PDstandard4th3alpha +
+ K21L*PDstandard4th3beta1 + K22L*PDstandard4th3beta2 + K32L*(PDstandard4th2beta2 + PDstandard4th3beta3) +
+ beta3L*PDstandard4th3K32 + alphaL*(-2*(K21L*Km13 + K22L*Km23 + K32L*Km33) + R32 + K32L*trK);
+
+ K33rhsL = G133*PDstandard4th1alpha + beta1L*PDstandard4th1K33 + G233*PDstandard4th2alpha + beta2L*PDstandard4th2K33 -
+ PDstandard4th33alpha + G333*PDstandard4th3alpha +
+ 2*(K31L*PDstandard4th3beta1 + K32L*PDstandard4th3beta2 + K33L*PDstandard4th3beta3) + beta3L*PDstandard4th3K33 +
+ alphaL*(-2*(K31L*Km13 + K32L*Km23 + K33L*Km33) + R33 + K33L*trK);
+
+ alpharhsL = 0;
+
+ beta1rhsL = 0;
+
+ beta2rhsL = 0;
+
+ beta3rhsL = 0;
+
+
+ /* Copy local copies back to grid functions */
+ alpharhs[index] = alpharhsL;
+ beta1rhs[index] = beta1rhsL;
+ beta2rhs[index] = beta2rhsL;
+ beta3rhs[index] = beta3rhsL;
+ g11rhs[index] = g11rhsL;
+ g21rhs[index] = g21rhsL;
+ g22rhs[index] = g22rhsL;
+ g31rhs[index] = g31rhsL;
+ g32rhs[index] = g32rhsL;
+ g33rhs[index] = g33rhsL;
+ K11rhs[index] = K11rhsL;
+ K21rhs[index] = K21rhsL;
+ K22rhs[index] = K22rhsL;
+ K31rhs[index] = K31rhsL;
+ K32rhs[index] = K32rhsL;
+ K33rhs[index] = K33rhsL;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ADM_RHS(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ADM_RHS_Body);
+}
diff --git a/ML_ADM/src/ADM_constraints.c b/ML_ADM/src/ADM_constraints.c
new file mode 100644
index 0000000..aadfe22
--- /dev/null
+++ b/ML_ADM/src/ADM_constraints.c
@@ -0,0 +1,445 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ADM_constraints_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ADM_constraints_Body");
+ }
+
+ if (cctk_iteration % ADM_constraints_calc_every != ADM_constraints_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+ CCTK_REAL detg = INITVALUE;
+ CCTK_REAL G111 = INITVALUE, G121 = INITVALUE, G122 = INITVALUE, G131 = INITVALUE, G132 = INITVALUE, G133 = INITVALUE;
+ CCTK_REAL G211 = INITVALUE, G221 = INITVALUE, G222 = INITVALUE, G231 = INITVALUE, G232 = INITVALUE, G233 = INITVALUE;
+ CCTK_REAL G311 = INITVALUE, G321 = INITVALUE, G322 = INITVALUE, G331 = INITVALUE, G332 = INITVALUE, G333 = INITVALUE;
+ CCTK_REAL gu11 = INITVALUE, gu12 = INITVALUE, gu13 = INITVALUE, gu22 = INITVALUE, gu23 = INITVALUE, gu33 = INITVALUE;
+ CCTK_REAL Km11 = INITVALUE, Km12 = INITVALUE, Km13 = INITVALUE, Km21 = INITVALUE, Km22 = INITVALUE, Km23 = INITVALUE;
+ CCTK_REAL Km31 = INITVALUE, Km32 = INITVALUE, Km33 = INITVALUE;
+ CCTK_REAL R11 = INITVALUE, R21 = INITVALUE, R22 = INITVALUE, R31 = INITVALUE, R32 = INITVALUE, R33 = INITVALUE;
+ CCTK_REAL trK = INITVALUE;
+ CCTK_REAL trR = INITVALUE;
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL HL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ CCTK_REAL M1L = INITVALUE, M2L = INITVALUE, M3L = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+ CCTK_REAL PDstandard4th1g11 = INITVALUE;
+ CCTK_REAL PDstandard4th2g11 = INITVALUE;
+ CCTK_REAL PDstandard4th3g11 = INITVALUE;
+ CCTK_REAL PDstandard4th22g11 = INITVALUE;
+ CCTK_REAL PDstandard4th33g11 = INITVALUE;
+ CCTK_REAL PDstandard4th12g11 = INITVALUE;
+ CCTK_REAL PDstandard4th13g11 = INITVALUE;
+ CCTK_REAL PDstandard4th21g11 = INITVALUE;
+ CCTK_REAL PDstandard4th23g11 = INITVALUE;
+ CCTK_REAL PDstandard4th31g11 = INITVALUE;
+ CCTK_REAL PDstandard4th32g11 = INITVALUE;
+ CCTK_REAL PDstandard4th1g21 = INITVALUE;
+ CCTK_REAL PDstandard4th2g21 = INITVALUE;
+ CCTK_REAL PDstandard4th3g21 = INITVALUE;
+ CCTK_REAL PDstandard4th33g21 = INITVALUE;
+ CCTK_REAL PDstandard4th12g21 = INITVALUE;
+ CCTK_REAL PDstandard4th13g21 = INITVALUE;
+ CCTK_REAL PDstandard4th21g21 = INITVALUE;
+ CCTK_REAL PDstandard4th23g21 = INITVALUE;
+ CCTK_REAL PDstandard4th31g21 = INITVALUE;
+ CCTK_REAL PDstandard4th32g21 = INITVALUE;
+ CCTK_REAL PDstandard4th1g22 = INITVALUE;
+ CCTK_REAL PDstandard4th2g22 = INITVALUE;
+ CCTK_REAL PDstandard4th3g22 = INITVALUE;
+ CCTK_REAL PDstandard4th11g22 = INITVALUE;
+ CCTK_REAL PDstandard4th33g22 = INITVALUE;
+ CCTK_REAL PDstandard4th12g22 = INITVALUE;
+ CCTK_REAL PDstandard4th13g22 = INITVALUE;
+ CCTK_REAL PDstandard4th21g22 = INITVALUE;
+ CCTK_REAL PDstandard4th23g22 = INITVALUE;
+ CCTK_REAL PDstandard4th31g22 = INITVALUE;
+ CCTK_REAL PDstandard4th32g22 = INITVALUE;
+ CCTK_REAL PDstandard4th1g31 = INITVALUE;
+ CCTK_REAL PDstandard4th2g31 = INITVALUE;
+ CCTK_REAL PDstandard4th3g31 = INITVALUE;
+ CCTK_REAL PDstandard4th22g31 = INITVALUE;
+ CCTK_REAL PDstandard4th12g31 = INITVALUE;
+ CCTK_REAL PDstandard4th13g31 = INITVALUE;
+ CCTK_REAL PDstandard4th21g31 = INITVALUE;
+ CCTK_REAL PDstandard4th23g31 = INITVALUE;
+ CCTK_REAL PDstandard4th31g31 = INITVALUE;
+ CCTK_REAL PDstandard4th32g31 = INITVALUE;
+ CCTK_REAL PDstandard4th1g32 = INITVALUE;
+ CCTK_REAL PDstandard4th2g32 = INITVALUE;
+ CCTK_REAL PDstandard4th3g32 = INITVALUE;
+ CCTK_REAL PDstandard4th11g32 = INITVALUE;
+ CCTK_REAL PDstandard4th12g32 = INITVALUE;
+ CCTK_REAL PDstandard4th13g32 = INITVALUE;
+ CCTK_REAL PDstandard4th21g32 = INITVALUE;
+ CCTK_REAL PDstandard4th23g32 = INITVALUE;
+ CCTK_REAL PDstandard4th31g32 = INITVALUE;
+ CCTK_REAL PDstandard4th32g32 = INITVALUE;
+ CCTK_REAL PDstandard4th1g33 = INITVALUE;
+ CCTK_REAL PDstandard4th2g33 = INITVALUE;
+ CCTK_REAL PDstandard4th3g33 = INITVALUE;
+ CCTK_REAL PDstandard4th11g33 = INITVALUE;
+ CCTK_REAL PDstandard4th22g33 = INITVALUE;
+ CCTK_REAL PDstandard4th12g33 = INITVALUE;
+ CCTK_REAL PDstandard4th13g33 = INITVALUE;
+ CCTK_REAL PDstandard4th21g33 = INITVALUE;
+ CCTK_REAL PDstandard4th23g33 = INITVALUE;
+ CCTK_REAL PDstandard4th31g33 = INITVALUE;
+ CCTK_REAL PDstandard4th32g33 = INITVALUE;
+ CCTK_REAL PDstandard4th2K11 = INITVALUE;
+ CCTK_REAL PDstandard4th3K11 = INITVALUE;
+ CCTK_REAL PDstandard4th1K21 = INITVALUE;
+ CCTK_REAL PDstandard4th2K21 = INITVALUE;
+ CCTK_REAL PDstandard4th3K21 = INITVALUE;
+ CCTK_REAL PDstandard4th1K22 = INITVALUE;
+ CCTK_REAL PDstandard4th3K22 = INITVALUE;
+ CCTK_REAL PDstandard4th1K31 = INITVALUE;
+ CCTK_REAL PDstandard4th2K31 = INITVALUE;
+ CCTK_REAL PDstandard4th3K31 = INITVALUE;
+ CCTK_REAL PDstandard4th1K32 = INITVALUE;
+ CCTK_REAL PDstandard4th2K32 = INITVALUE;
+ CCTK_REAL PDstandard4th3K32 = INITVALUE;
+ CCTK_REAL PDstandard4th1K33 = INITVALUE;
+ CCTK_REAL PDstandard4th2K33 = INITVALUE;
+
+ /* Assign local copies of grid functions */
+ g11L = g11[index];
+ g21L = g21[index];
+ g22L = g22[index];
+ g31L = g31[index];
+ g32L = g32[index];
+ g33L = g33[index];
+ K11L = K11[index];
+ K21L = K21[index];
+ K22L = K22[index];
+ K31L = K31[index];
+ K32L = K32[index];
+ K33L = K33[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+ PDstandard4th1g11 = PDstandard4th1(g11, i, j, k);
+ PDstandard4th2g11 = PDstandard4th2(g11, i, j, k);
+ PDstandard4th3g11 = PDstandard4th3(g11, i, j, k);
+ PDstandard4th22g11 = PDstandard4th22(g11, i, j, k);
+ PDstandard4th33g11 = PDstandard4th33(g11, i, j, k);
+ PDstandard4th23g11 = PDstandard4th23(g11, i, j, k);
+ PDstandard4th1g21 = PDstandard4th1(g21, i, j, k);
+ PDstandard4th2g21 = PDstandard4th2(g21, i, j, k);
+ PDstandard4th3g21 = PDstandard4th3(g21, i, j, k);
+ PDstandard4th33g21 = PDstandard4th33(g21, i, j, k);
+ PDstandard4th12g21 = PDstandard4th12(g21, i, j, k);
+ PDstandard4th13g21 = PDstandard4th13(g21, i, j, k);
+ PDstandard4th23g21 = PDstandard4th23(g21, i, j, k);
+ PDstandard4th1g22 = PDstandard4th1(g22, i, j, k);
+ PDstandard4th2g22 = PDstandard4th2(g22, i, j, k);
+ PDstandard4th3g22 = PDstandard4th3(g22, i, j, k);
+ PDstandard4th11g22 = PDstandard4th11(g22, i, j, k);
+ PDstandard4th33g22 = PDstandard4th33(g22, i, j, k);
+ PDstandard4th13g22 = PDstandard4th13(g22, i, j, k);
+ PDstandard4th1g31 = PDstandard4th1(g31, i, j, k);
+ PDstandard4th2g31 = PDstandard4th2(g31, i, j, k);
+ PDstandard4th3g31 = PDstandard4th3(g31, i, j, k);
+ PDstandard4th22g31 = PDstandard4th22(g31, i, j, k);
+ PDstandard4th12g31 = PDstandard4th12(g31, i, j, k);
+ PDstandard4th13g31 = PDstandard4th13(g31, i, j, k);
+ PDstandard4th23g31 = PDstandard4th23(g31, i, j, k);
+ PDstandard4th1g32 = PDstandard4th1(g32, i, j, k);
+ PDstandard4th2g32 = PDstandard4th2(g32, i, j, k);
+ PDstandard4th3g32 = PDstandard4th3(g32, i, j, k);
+ PDstandard4th11g32 = PDstandard4th11(g32, i, j, k);
+ PDstandard4th12g32 = PDstandard4th12(g32, i, j, k);
+ PDstandard4th13g32 = PDstandard4th13(g32, i, j, k);
+ PDstandard4th23g32 = PDstandard4th23(g32, i, j, k);
+ PDstandard4th1g33 = PDstandard4th1(g33, i, j, k);
+ PDstandard4th2g33 = PDstandard4th2(g33, i, j, k);
+ PDstandard4th3g33 = PDstandard4th3(g33, i, j, k);
+ PDstandard4th11g33 = PDstandard4th11(g33, i, j, k);
+ PDstandard4th22g33 = PDstandard4th22(g33, i, j, k);
+ PDstandard4th12g33 = PDstandard4th12(g33, i, j, k);
+ PDstandard4th2K11 = PDstandard4th2(K11, i, j, k);
+ PDstandard4th3K11 = PDstandard4th3(K11, i, j, k);
+ PDstandard4th1K21 = PDstandard4th1(K21, i, j, k);
+ PDstandard4th2K21 = PDstandard4th2(K21, i, j, k);
+ PDstandard4th3K21 = PDstandard4th3(K21, i, j, k);
+ PDstandard4th1K22 = PDstandard4th1(K22, i, j, k);
+ PDstandard4th3K22 = PDstandard4th3(K22, i, j, k);
+ PDstandard4th1K31 = PDstandard4th1(K31, i, j, k);
+ PDstandard4th2K31 = PDstandard4th2(K31, i, j, k);
+ PDstandard4th3K31 = PDstandard4th3(K31, i, j, k);
+ PDstandard4th1K32 = PDstandard4th1(K32, i, j, k);
+ PDstandard4th2K32 = PDstandard4th2(K32, i, j, k);
+ PDstandard4th3K32 = PDstandard4th3(K32, i, j, k);
+ PDstandard4th1K33 = PDstandard4th1(K33, i, j, k);
+ PDstandard4th2K33 = PDstandard4th2(K33, i, j, k);
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ detg = 2*g21L*g31L*g32L + g33L*(g11L*g22L - SQR(g21L)) - g22L*SQR(g31L) - g11L*SQR(g32L);
+
+ gu11 = INV(detg)*(g22L*g33L - SQR(g32L));
+
+ gu12 = (g31L*g32L - g21L*g33L)*INV(detg);
+
+ gu13 = (-(g22L*g31L) + g21L*g32L)*INV(detg);
+
+ gu22 = INV(detg)*(g11L*g33L - SQR(g31L));
+
+ gu23 = (g21L*g31L - g11L*g32L)*INV(detg);
+
+ gu33 = INV(detg)*(g11L*g22L - SQR(g21L));
+
+ G111 = khalf*(gu11*PDstandard4th1g11 + 2*(gu12*PDstandard4th1g21 + gu13*PDstandard4th1g31) - gu12*PDstandard4th2g11 -
+ gu13*PDstandard4th3g11);
+
+ G211 = khalf*(gu12*PDstandard4th1g11 + 2*(gu22*PDstandard4th1g21 + gu23*PDstandard4th1g31) - gu22*PDstandard4th2g11 -
+ gu23*PDstandard4th3g11);
+
+ G311 = khalf*(gu13*PDstandard4th1g11 + 2*(gu23*PDstandard4th1g21 + gu33*PDstandard4th1g31) - gu23*PDstandard4th2g11 -
+ gu33*PDstandard4th3g11);
+
+ G121 = khalf*(gu12*PDstandard4th1g22 + gu11*PDstandard4th2g11 +
+ gu13*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G221 = khalf*(gu22*PDstandard4th1g22 + gu12*PDstandard4th2g11 +
+ gu23*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G321 = khalf*(gu23*PDstandard4th1g22 + gu13*PDstandard4th2g11 +
+ gu33*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G131 = khalf*(gu13*PDstandard4th1g33 + gu11*PDstandard4th3g11 +
+ gu12*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G231 = khalf*(gu23*PDstandard4th1g33 + gu12*PDstandard4th3g11 +
+ gu22*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G331 = khalf*(gu33*PDstandard4th1g33 + gu13*PDstandard4th3g11 +
+ gu23*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G122 = khalf*(gu11*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu12*PDstandard4th2g22 +
+ gu13*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G222 = khalf*(gu12*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu22*PDstandard4th2g22 +
+ gu23*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G322 = khalf*(gu13*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu23*PDstandard4th2g22 +
+ gu33*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G132 = khalf*(gu13*PDstandard4th2g33 + gu11*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu12*PDstandard4th3g22);
+
+ G232 = khalf*(gu23*PDstandard4th2g33 + gu12*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu22*PDstandard4th3g22);
+
+ G332 = khalf*(gu33*PDstandard4th2g33 + gu13*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu23*PDstandard4th3g22);
+
+ G133 = khalf*(-(gu11*PDstandard4th1g33) - gu12*PDstandard4th2g33 + 2*gu11*PDstandard4th3g31 +
+ 2*gu12*PDstandard4th3g32 + gu13*PDstandard4th3g33);
+
+ G233 = khalf*(-(gu12*PDstandard4th1g33) - gu22*PDstandard4th2g33 + 2*gu12*PDstandard4th3g31 +
+ 2*gu22*PDstandard4th3g32 + gu23*PDstandard4th3g33);
+
+ G333 = khalf*(-(gu13*PDstandard4th1g33) - gu23*PDstandard4th2g33 + 2*gu13*PDstandard4th3g31 +
+ 2*gu23*PDstandard4th3g32 + gu33*PDstandard4th3g33);
+
+ R11 = khalf*(-(gu22*PDstandard4th11g22) - 2*(G111*G122 + G111*G133 + G211*G222 + G211*G233 + G311*G322 + G311*G333 +
+ gu23*PDstandard4th11g32) - gu33*PDstandard4th11g33 + 2*gu22*PDstandard4th12g21 + 2*gu23*PDstandard4th12g31 +
+ 2*gu23*PDstandard4th13g21 + 2*gu33*PDstandard4th13g31 - gu22*PDstandard4th22g11 - 2*gu23*PDstandard4th23g11 -
+ gu33*PDstandard4th33g11 + 2*SQR(G121) + 2*SQR(G131) + 2*SQR(G221) + 2*SQR(G231) + 2*SQR(G321) + 2*SQR(G331));
+
+ R21 = khalf*(2*(G131*G132 + G231*G232 + G331*G332) - 2*(G121*G133 + G221*G233 + G321*G333 + gu12*PDstandard4th12g21) -
+ gu23*PDstandard4th12g32 - gu33*PDstandard4th12g33 +
+ gu13*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21) + gu23*PDstandard4th13g22 +
+ gu33*PDstandard4th13g32 + gu12*(PDstandard4th11g22 + PDstandard4th22g11) + gu23*PDstandard4th22g31 +
+ gu13*PDstandard4th23g11 - gu23*PDstandard4th23g21 + gu33*PDstandard4th23g31 - gu33*PDstandard4th33g21);
+
+ R31 = khalf*(2*(G121*G132 + G221*G232 + G321*G332) - 2*(G122*G131 + G222*G231 + G322*G331 + gu13*PDstandard4th13g31) +
+ gu12*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21 + PDstandard4th23g11) +
+ gu22*(PDstandard4th12g32 - PDstandard4th13g22 - PDstandard4th22g31 + PDstandard4th23g21) +
+ gu13*(PDstandard4th11g33 + PDstandard4th33g11) +
+ gu23*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21));
+
+ R22 = khalf*(-2*(G122*(G111 + G133) + G222*(G211 + G233) + G322*(G311 + G333) + gu13*PDstandard4th13g22) +
+ gu11*(-PDstandard4th11g22 + 2*PDstandard4th12g21 - PDstandard4th22g11) +
+ gu13*(-2*PDstandard4th22g31 + 2*(PDstandard4th12g32 + PDstandard4th23g21)) +
+ gu33*(-PDstandard4th22g33 + 2*PDstandard4th23g32 - PDstandard4th33g22) +
+ 2*(SQR(G121) + SQR(G132) + SQR(G221) + SQR(G232) + SQR(G321) + SQR(G332)));
+
+ R32 = khalf*(2*(G121*G131 + G221*G231 + G321*G331) +
+ gu11*(-PDstandard4th11g32 + PDstandard4th12g31 + PDstandard4th13g21 - PDstandard4th23g11) +
+ gu12*(-PDstandard4th12g32 + PDstandard4th13g22 + PDstandard4th22g31 - PDstandard4th23g21) -
+ 2*(G111*G132 + G211*G232 + G311*G332 + gu23*PDstandard4th23g32) +
+ gu13*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21) +
+ gu23*(PDstandard4th22g33 + PDstandard4th33g22));
+
+ R33 = khalf*(gu11*(-PDstandard4th11g33 + 2*PDstandard4th13g31 - PDstandard4th33g11) -
+ 2*((G111 + G122)*G133 + (G211 + G222)*G233 + (G311 + G322)*G333 +
+ gu12*(PDstandard4th12g33 + PDstandard4th33g21)) +
+ gu22*(-PDstandard4th22g33 + 2*PDstandard4th23g32 - PDstandard4th33g22) +
+ 2*(gu12*(PDstandard4th13g32 + PDstandard4th23g31) + SQR(G131) + SQR(G132) + SQR(G231) + SQR(G232) + SQR(G331) +
+ SQR(G332)));
+
+ trR = gu11*R11 + gu22*R22 + 2*(gu12*R21 + gu13*R31 + gu23*R32) + gu33*R33;
+
+ Km11 = gu11*K11L + gu12*K21L + gu13*K31L;
+
+ Km21 = gu12*K11L + gu22*K21L + gu23*K31L;
+
+ Km31 = gu13*K11L + gu23*K21L + gu33*K31L;
+
+ Km12 = gu11*K21L + gu12*K22L + gu13*K32L;
+
+ Km22 = gu12*K21L + gu22*K22L + gu23*K32L;
+
+ Km32 = gu13*K21L + gu23*K22L + gu33*K32L;
+
+ Km13 = gu11*K31L + gu12*K32L + gu13*K33L;
+
+ Km23 = gu12*K31L + gu22*K32L + gu23*K33L;
+
+ Km33 = gu13*K31L + gu23*K32L + gu33*K33L;
+
+ trK = Km11 + Km22 + Km33;
+
+ HL = -2*(Km12*Km21 + Km13*Km31 + Km23*Km32) + trR - SQR(Km11) - SQR(Km22) - SQR(Km33) + SQR(trK);
+
+ M1L = gu12*(-(G121*K11L) + G111*K21L - G221*K21L + G211*K22L - G321*K31L + G311*K32L - PDstandard4th1K21 +
+ PDstandard4th2K11) + gu22*(-(G122*K11L) + G121*K21L - G222*K21L + G221*K22L - G322*K31L + G321*K32L -
+ PDstandard4th1K22 + PDstandard4th2K21) + gu13*
+ (-(G131*K11L) - G231*K21L + G111*K31L - G331*K31L + G211*K32L + G311*K33L - PDstandard4th1K31 + PDstandard4th3K11)\
+ + gu23*(G131*K21L + G231*K22L + G121*K31L + (G221 + G331)*K32L + G321*K33L -
+ 2*(G132*K11L + G232*K21L + G332*K31L + PDstandard4th1K32) + PDstandard4th2K31 + PDstandard4th3K21) +
+ gu33*(-(G133*K11L) - G233*K21L + G131*K31L - G333*K31L + G231*K32L + G331*K33L - PDstandard4th1K33 +
+ PDstandard4th3K31);
+
+ M2L = gu11*(G121*K11L + (-G111 + G221)*K21L - G211*K22L + G321*K31L - G311*K32L + PDstandard4th1K21 -
+ PDstandard4th2K11) + gu12*(G122*K11L + (-G121 + G222)*K21L - G221*K22L + G322*K31L - G321*K32L +
+ PDstandard4th1K22 - PDstandard4th2K21) + gu13*
+ (G132*K11L + (-2*G131 + G232)*K21L + (G121 + G332)*K31L + G221*K32L + G321*K33L + PDstandard4th1K32 -
+ 2*(G231*K22L + G331*K32L + PDstandard4th2K31) + PDstandard4th3K21) +
+ gu23*(-(G132*K21L) - G232*K22L + G122*K31L + G222*K32L - G332*K32L + G322*K33L - PDstandard4th2K32 +
+ PDstandard4th3K22) + gu33*(-(G133*K21L) - G233*K22L + G132*K31L + G232*K32L - G333*K32L + G332*K33L -
+ PDstandard4th2K33 + PDstandard4th3K32);
+
+ M3L = gu11*(G131*K11L + G231*K21L + (-G111 + G331)*K31L - G211*K32L - G311*K33L + PDstandard4th1K31 -
+ PDstandard4th3K11) + gu12*(G132*K11L + (G131 + G232)*K21L + G231*K22L + (-2*G121 + G332)*K31L +
+ (-2*G221 + G331)*K32L + PDstandard4th1K32 + PDstandard4th2K31 - 2*(G321*K33L + PDstandard4th3K21)) +
+ gu22*(G132*K21L + G232*K22L - G122*K31L - G222*K32L + G332*K32L - G322*K33L + PDstandard4th2K32 -
+ PDstandard4th3K22) + gu13*(G133*K11L + G233*K21L + (-G131 + G333)*K31L - G231*K32L - G331*K33L +
+ PDstandard4th1K33 - PDstandard4th3K31) + gu23*
+ (G133*K21L + G233*K22L - G132*K31L - G232*K32L + G333*K32L - G332*K33L + PDstandard4th2K33 - PDstandard4th3K32);
+
+
+ /* Copy local copies back to grid functions */
+ H[index] = HL;
+ M1[index] = M1L;
+ M2[index] = M2L;
+ M3[index] = M3L;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ADM_constraints(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ADM_constraints_Body);
+}
diff --git a/ML_ADM/src/ADM_convertFromADMBase.c b/ML_ADM/src/ADM_convertFromADMBase.c
new file mode 100644
index 0000000..c3275a4
--- /dev/null
+++ b/ML_ADM/src/ADM_convertFromADMBase.c
@@ -0,0 +1,212 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ADM_convertFromADMBase_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ADM_convertFromADMBase_Body");
+ }
+
+ if (cctk_iteration % ADM_convertFromADMBase_calc_every != ADM_convertFromADMBase_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alpL = INITVALUE;
+ CCTK_REAL alphaL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE;
+ CCTK_REAL betaxL = INITVALUE;
+ CCTK_REAL betayL = INITVALUE;
+ CCTK_REAL betazL = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL gxxL = INITVALUE;
+ CCTK_REAL gxyL = INITVALUE;
+ CCTK_REAL gxzL = INITVALUE;
+ CCTK_REAL gyyL = INITVALUE;
+ CCTK_REAL gyzL = INITVALUE;
+ CCTK_REAL gzzL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ CCTK_REAL kxxL = INITVALUE;
+ CCTK_REAL kxyL = INITVALUE;
+ CCTK_REAL kxzL = INITVALUE;
+ CCTK_REAL kyyL = INITVALUE;
+ CCTK_REAL kyzL = INITVALUE;
+ CCTK_REAL kzzL = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+
+ /* Assign local copies of grid functions */
+ alpL = alp[index];
+ betaxL = betax[index];
+ betayL = betay[index];
+ betazL = betaz[index];
+ gxxL = gxx[index];
+ gxyL = gxy[index];
+ gxzL = gxz[index];
+ gyyL = gyy[index];
+ gyzL = gyz[index];
+ gzzL = gzz[index];
+ kxxL = kxx[index];
+ kxyL = kxy[index];
+ kxzL = kxz[index];
+ kyyL = kyy[index];
+ kyzL = kyz[index];
+ kzzL = kzz[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ g11L = gxxL;
+
+ g21L = gxyL;
+
+ g31L = gxzL;
+
+ g22L = gyyL;
+
+ g32L = gyzL;
+
+ g33L = gzzL;
+
+ K11L = kxxL;
+
+ K21L = kxyL;
+
+ K31L = kxzL;
+
+ K22L = kyyL;
+
+ K32L = kyzL;
+
+ K33L = kzzL;
+
+ alphaL = alpL;
+
+ beta1L = betaxL;
+
+ beta2L = betayL;
+
+ beta3L = betazL;
+
+
+ /* Copy local copies back to grid functions */
+ alpha[index] = alphaL;
+ beta1[index] = beta1L;
+ beta2[index] = beta2L;
+ beta3[index] = beta3L;
+ g11[index] = g11L;
+ g21[index] = g21L;
+ g22[index] = g22L;
+ g31[index] = g31L;
+ g32[index] = g32L;
+ g33[index] = g33L;
+ K11[index] = K11L;
+ K21[index] = K21L;
+ K22[index] = K22L;
+ K31[index] = K31L;
+ K32[index] = K32L;
+ K33[index] = K33L;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ADM_convertFromADMBase(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ADM_convertFromADMBase_Body);
+}
diff --git a/ML_ADM/src/ADM_convertToADMBase.c b/ML_ADM/src/ADM_convertToADMBase.c
new file mode 100644
index 0000000..77d9852
--- /dev/null
+++ b/ML_ADM/src/ADM_convertToADMBase.c
@@ -0,0 +1,228 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ADM_convertToADMBase_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ADM_convertToADMBase_Body");
+ }
+
+ if (cctk_iteration % ADM_convertToADMBase_calc_every != ADM_convertToADMBase_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alpL = INITVALUE;
+ CCTK_REAL alphaL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE;
+ CCTK_REAL betaxL = INITVALUE;
+ CCTK_REAL betayL = INITVALUE;
+ CCTK_REAL betazL = INITVALUE;
+ CCTK_REAL dtalpL = INITVALUE;
+ CCTK_REAL dtbetaxL = INITVALUE;
+ CCTK_REAL dtbetayL = INITVALUE;
+ CCTK_REAL dtbetazL = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL gxxL = INITVALUE;
+ CCTK_REAL gxyL = INITVALUE;
+ CCTK_REAL gxzL = INITVALUE;
+ CCTK_REAL gyyL = INITVALUE;
+ CCTK_REAL gyzL = INITVALUE;
+ CCTK_REAL gzzL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ CCTK_REAL kxxL = INITVALUE;
+ CCTK_REAL kxyL = INITVALUE;
+ CCTK_REAL kxzL = INITVALUE;
+ CCTK_REAL kyyL = INITVALUE;
+ CCTK_REAL kyzL = INITVALUE;
+ CCTK_REAL kzzL = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+
+ /* Assign local copies of grid functions */
+ alphaL = alpha[index];
+ beta1L = beta1[index];
+ beta2L = beta2[index];
+ beta3L = beta3[index];
+ g11L = g11[index];
+ g21L = g21[index];
+ g22L = g22[index];
+ g31L = g31[index];
+ g32L = g32[index];
+ g33L = g33[index];
+ K11L = K11[index];
+ K21L = K21[index];
+ K22L = K22[index];
+ K31L = K31[index];
+ K32L = K32[index];
+ K33L = K33[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ gxxL = g11L;
+
+ gxyL = g21L;
+
+ gxzL = g31L;
+
+ gyyL = g22L;
+
+ gyzL = g32L;
+
+ gzzL = g33L;
+
+ kxxL = K11L;
+
+ kxyL = K21L;
+
+ kxzL = K31L;
+
+ kyyL = K22L;
+
+ kyzL = K32L;
+
+ kzzL = K33L;
+
+ alpL = alphaL;
+
+ dtalpL = 0;
+
+ betaxL = beta1L;
+
+ betayL = beta2L;
+
+ betazL = beta3L;
+
+ dtbetaxL = 0;
+
+ dtbetayL = 0;
+
+ dtbetazL = 0;
+
+
+ /* Copy local copies back to grid functions */
+ alp[index] = alpL;
+ betax[index] = betaxL;
+ betay[index] = betayL;
+ betaz[index] = betazL;
+ dtalp[index] = dtalpL;
+ dtbetax[index] = dtbetaxL;
+ dtbetay[index] = dtbetayL;
+ dtbetaz[index] = dtbetazL;
+ gxx[index] = gxxL;
+ gxy[index] = gxyL;
+ gxz[index] = gxzL;
+ gyy[index] = gyyL;
+ gyz[index] = gyzL;
+ gzz[index] = gzzL;
+ kxx[index] = kxxL;
+ kxy[index] = kxyL;
+ kxz[index] = kxzL;
+ kyy[index] = kyyL;
+ kyz[index] = kyzL;
+ kzz[index] = kzzL;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ADM_convertToADMBase(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ADM_convertToADMBase_Body);
+}
diff --git a/ML_ADM/src/Boundaries.c b/ML_ADM/src/Boundaries.c
new file mode 100644
index 0000000..a79a204
--- /dev/null
+++ b/ML_ADM/src/Boundaries.c
@@ -0,0 +1,1044 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "cctk_Faces.h"
+#include "util_Table.h"
+#include "Symmetry.h"
+
+
+/* the boundary treatment is split into 3 steps: */
+/* 1. excision */
+/* 2. symmetries */
+/* 3. "other" boundary conditions, e.g. radiative */
+
+/* to simplify scheduling and testing, the 3 steps */
+/* are currently applied in separate functions */
+
+
+void ML_ADM_CheckBoundaries(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ return;
+}
+
+void ML_ADM_ApplyBoundConds(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ CCTK_INT ierr = 0;
+
+ if (CCTK_EQUALS(curv_bound, "none" ) ||
+ CCTK_EQUALS(curv_bound, "static") ||
+ CCTK_EQUALS(curv_bound, "flat" ) ||
+ CCTK_EQUALS(curv_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::curv", curv_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register curv_bound BC for ML_ADM::curv!");
+ }
+
+ if (CCTK_EQUALS(lapse_bound, "none" ) ||
+ CCTK_EQUALS(lapse_bound, "static") ||
+ CCTK_EQUALS(lapse_bound, "flat" ) ||
+ CCTK_EQUALS(lapse_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::lapse", lapse_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register lapse_bound BC for ML_ADM::lapse!");
+ }
+
+ if (CCTK_EQUALS(metric_bound, "none" ) ||
+ CCTK_EQUALS(metric_bound, "static") ||
+ CCTK_EQUALS(metric_bound, "flat" ) ||
+ CCTK_EQUALS(metric_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::metric", metric_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register metric_bound BC for ML_ADM::metric!");
+ }
+
+ if (CCTK_EQUALS(shift_bound, "none" ) ||
+ CCTK_EQUALS(shift_bound, "static") ||
+ CCTK_EQUALS(shift_bound, "flat" ) ||
+ CCTK_EQUALS(shift_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::shift", shift_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register shift_bound BC for ML_ADM::shift!");
+ }
+
+ if (CCTK_EQUALS(K11_bound, "none" ) ||
+ CCTK_EQUALS(K11_bound, "static") ||
+ CCTK_EQUALS(K11_bound, "flat" ) ||
+ CCTK_EQUALS(K11_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::K11", K11_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register K11_bound BC for ML_ADM::K11!");
+ }
+
+ if (CCTK_EQUALS(K21_bound, "none" ) ||
+ CCTK_EQUALS(K21_bound, "static") ||
+ CCTK_EQUALS(K21_bound, "flat" ) ||
+ CCTK_EQUALS(K21_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::K21", K21_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register K21_bound BC for ML_ADM::K21!");
+ }
+
+ if (CCTK_EQUALS(K31_bound, "none" ) ||
+ CCTK_EQUALS(K31_bound, "static") ||
+ CCTK_EQUALS(K31_bound, "flat" ) ||
+ CCTK_EQUALS(K31_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::K31", K31_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register K31_bound BC for ML_ADM::K31!");
+ }
+
+ if (CCTK_EQUALS(K22_bound, "none" ) ||
+ CCTK_EQUALS(K22_bound, "static") ||
+ CCTK_EQUALS(K22_bound, "flat" ) ||
+ CCTK_EQUALS(K22_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::K22", K22_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register K22_bound BC for ML_ADM::K22!");
+ }
+
+ if (CCTK_EQUALS(K32_bound, "none" ) ||
+ CCTK_EQUALS(K32_bound, "static") ||
+ CCTK_EQUALS(K32_bound, "flat" ) ||
+ CCTK_EQUALS(K32_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::K32", K32_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register K32_bound BC for ML_ADM::K32!");
+ }
+
+ if (CCTK_EQUALS(K33_bound, "none" ) ||
+ CCTK_EQUALS(K33_bound, "static") ||
+ CCTK_EQUALS(K33_bound, "flat" ) ||
+ CCTK_EQUALS(K33_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::K33", K33_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register K33_bound BC for ML_ADM::K33!");
+ }
+
+ if (CCTK_EQUALS(alpha_bound, "none" ) ||
+ CCTK_EQUALS(alpha_bound, "static") ||
+ CCTK_EQUALS(alpha_bound, "flat" ) ||
+ CCTK_EQUALS(alpha_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::alpha", alpha_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register alpha_bound BC for ML_ADM::alpha!");
+ }
+
+ if (CCTK_EQUALS(g11_bound, "none" ) ||
+ CCTK_EQUALS(g11_bound, "static") ||
+ CCTK_EQUALS(g11_bound, "flat" ) ||
+ CCTK_EQUALS(g11_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::g11", g11_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register g11_bound BC for ML_ADM::g11!");
+ }
+
+ if (CCTK_EQUALS(g21_bound, "none" ) ||
+ CCTK_EQUALS(g21_bound, "static") ||
+ CCTK_EQUALS(g21_bound, "flat" ) ||
+ CCTK_EQUALS(g21_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::g21", g21_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register g21_bound BC for ML_ADM::g21!");
+ }
+
+ if (CCTK_EQUALS(g31_bound, "none" ) ||
+ CCTK_EQUALS(g31_bound, "static") ||
+ CCTK_EQUALS(g31_bound, "flat" ) ||
+ CCTK_EQUALS(g31_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::g31", g31_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register g31_bound BC for ML_ADM::g31!");
+ }
+
+ if (CCTK_EQUALS(g22_bound, "none" ) ||
+ CCTK_EQUALS(g22_bound, "static") ||
+ CCTK_EQUALS(g22_bound, "flat" ) ||
+ CCTK_EQUALS(g22_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::g22", g22_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register g22_bound BC for ML_ADM::g22!");
+ }
+
+ if (CCTK_EQUALS(g32_bound, "none" ) ||
+ CCTK_EQUALS(g32_bound, "static") ||
+ CCTK_EQUALS(g32_bound, "flat" ) ||
+ CCTK_EQUALS(g32_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::g32", g32_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register g32_bound BC for ML_ADM::g32!");
+ }
+
+ if (CCTK_EQUALS(g33_bound, "none" ) ||
+ CCTK_EQUALS(g33_bound, "static") ||
+ CCTK_EQUALS(g33_bound, "flat" ) ||
+ CCTK_EQUALS(g33_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::g33", g33_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register g33_bound BC for ML_ADM::g33!");
+ }
+
+ if (CCTK_EQUALS(beta1_bound, "none" ) ||
+ CCTK_EQUALS(beta1_bound, "static") ||
+ CCTK_EQUALS(beta1_bound, "flat" ) ||
+ CCTK_EQUALS(beta1_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::beta1", beta1_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register beta1_bound BC for ML_ADM::beta1!");
+ }
+
+ if (CCTK_EQUALS(beta2_bound, "none" ) ||
+ CCTK_EQUALS(beta2_bound, "static") ||
+ CCTK_EQUALS(beta2_bound, "flat" ) ||
+ CCTK_EQUALS(beta2_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::beta2", beta2_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register beta2_bound BC for ML_ADM::beta2!");
+ }
+
+ if (CCTK_EQUALS(beta3_bound, "none" ) ||
+ CCTK_EQUALS(beta3_bound, "static") ||
+ CCTK_EQUALS(beta3_bound, "flat" ) ||
+ CCTK_EQUALS(beta3_bound, "zero" ) )
+ {
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, -1,
+ "ML_ADM::beta3", beta3_bound);
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register beta3_bound BC for ML_ADM::beta3!");
+ }
+
+ if (CCTK_EQUALS(curv_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_curv_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_curv_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_curv_bound , curv_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_curv_bound ,curv_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_curv_bound,
+ "ML_ADM::curv", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::curv!");
+
+ }
+
+ if (CCTK_EQUALS(lapse_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_lapse_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_lapse_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_lapse_bound , lapse_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_lapse_bound ,lapse_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_lapse_bound,
+ "ML_ADM::lapse", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::lapse!");
+
+ }
+
+ if (CCTK_EQUALS(metric_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_metric_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_metric_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_metric_bound , metric_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_metric_bound ,metric_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_metric_bound,
+ "ML_ADM::metric", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::metric!");
+
+ }
+
+ if (CCTK_EQUALS(shift_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_shift_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_shift_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_shift_bound , shift_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_shift_bound ,shift_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_shift_bound,
+ "ML_ADM::shift", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::shift!");
+
+ }
+
+ if (CCTK_EQUALS(K11_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_K11_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K11_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K11_bound , K11_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_K11_bound ,K11_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K11_bound,
+ "ML_ADM::K11", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::K11!");
+
+ }
+
+ if (CCTK_EQUALS(K21_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_K21_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K21_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K21_bound , K21_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_K21_bound ,K21_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K21_bound,
+ "ML_ADM::K21", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::K21!");
+
+ }
+
+ if (CCTK_EQUALS(K31_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_K31_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K31_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K31_bound , K31_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_K31_bound ,K31_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K31_bound,
+ "ML_ADM::K31", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::K31!");
+
+ }
+
+ if (CCTK_EQUALS(K22_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_K22_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K22_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K22_bound , K22_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_K22_bound ,K22_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K22_bound,
+ "ML_ADM::K22", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::K22!");
+
+ }
+
+ if (CCTK_EQUALS(K32_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_K32_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K32_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K32_bound , K32_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_K32_bound ,K32_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K32_bound,
+ "ML_ADM::K32", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::K32!");
+
+ }
+
+ if (CCTK_EQUALS(K33_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_K33_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K33_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K33_bound , K33_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_K33_bound ,K33_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K33_bound,
+ "ML_ADM::K33", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::K33!");
+
+ }
+
+ if (CCTK_EQUALS(alpha_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_alpha_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_alpha_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_alpha_bound , alpha_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_alpha_bound ,alpha_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_alpha_bound,
+ "ML_ADM::alpha", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::alpha!");
+
+ }
+
+ if (CCTK_EQUALS(g11_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_g11_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g11_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g11_bound , g11_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_g11_bound ,g11_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g11_bound,
+ "ML_ADM::g11", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::g11!");
+
+ }
+
+ if (CCTK_EQUALS(g21_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_g21_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g21_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g21_bound , g21_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_g21_bound ,g21_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g21_bound,
+ "ML_ADM::g21", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::g21!");
+
+ }
+
+ if (CCTK_EQUALS(g31_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_g31_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g31_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g31_bound , g31_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_g31_bound ,g31_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g31_bound,
+ "ML_ADM::g31", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::g31!");
+
+ }
+
+ if (CCTK_EQUALS(g22_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_g22_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g22_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g22_bound , g22_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_g22_bound ,g22_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g22_bound,
+ "ML_ADM::g22", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::g22!");
+
+ }
+
+ if (CCTK_EQUALS(g32_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_g32_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g32_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g32_bound , g32_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_g32_bound ,g32_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g32_bound,
+ "ML_ADM::g32", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::g32!");
+
+ }
+
+ if (CCTK_EQUALS(g33_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_g33_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g33_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g33_bound , g33_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_g33_bound ,g33_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g33_bound,
+ "ML_ADM::g33", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::g33!");
+
+ }
+
+ if (CCTK_EQUALS(beta1_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_beta1_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_beta1_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_beta1_bound , beta1_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_beta1_bound ,beta1_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_beta1_bound,
+ "ML_ADM::beta1", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::beta1!");
+
+ }
+
+ if (CCTK_EQUALS(beta2_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_beta2_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_beta2_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_beta2_bound , beta2_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_beta2_bound ,beta2_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_beta2_bound,
+ "ML_ADM::beta2", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::beta2!");
+
+ }
+
+ if (CCTK_EQUALS(beta3_bound, "radiative"))
+ {
+ /* apply radiation boundary condition */
+ CCTK_INT handle_beta3_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_beta3_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_beta3_bound , beta3_bound_limit, "LIMIT") < 0)
+ CCTK_WARN(-1, "could not set LIMIT value in table!");
+ if (Util_TableSetReal(handle_beta3_bound ,beta3_bound_speed, "SPEED") < 0)
+ CCTK_WARN(-1, "could not set SPEED value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_beta3_bound,
+ "ML_ADM::beta3", "Radiation");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Radiation BC for ML_ADM::beta3!");
+
+ }
+
+ if (CCTK_EQUALS(curv_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_curv_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_curv_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_curv_bound ,curv_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_curv_bound,
+ "ML_ADM::curv", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Scalar BC for ML_ADM::curv!");
+
+ }
+
+ if (CCTK_EQUALS(lapse_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_lapse_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_lapse_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_lapse_bound ,lapse_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_lapse_bound,
+ "ML_ADM::lapse", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Scalar BC for ML_ADM::lapse!");
+
+ }
+
+ if (CCTK_EQUALS(metric_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_metric_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_metric_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_metric_bound ,metric_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_metric_bound,
+ "ML_ADM::metric", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Scalar BC for ML_ADM::metric!");
+
+ }
+
+ if (CCTK_EQUALS(shift_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_shift_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_shift_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_shift_bound ,shift_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectGroupForBC(cctkGH, CCTK_ALL_FACES, 1, handle_shift_bound,
+ "ML_ADM::shift", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Failed to register Scalar BC for ML_ADM::shift!");
+
+ }
+
+ if (CCTK_EQUALS(K11_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_K11_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K11_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K11_bound ,K11_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K11_bound,
+ "ML_ADM::K11", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::K11!");
+
+ }
+
+ if (CCTK_EQUALS(K21_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_K21_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K21_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K21_bound ,K21_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K21_bound,
+ "ML_ADM::K21", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::K21!");
+
+ }
+
+ if (CCTK_EQUALS(K31_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_K31_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K31_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K31_bound ,K31_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K31_bound,
+ "ML_ADM::K31", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::K31!");
+
+ }
+
+ if (CCTK_EQUALS(K22_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_K22_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K22_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K22_bound ,K22_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K22_bound,
+ "ML_ADM::K22", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::K22!");
+
+ }
+
+ if (CCTK_EQUALS(K32_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_K32_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K32_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K32_bound ,K32_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K32_bound,
+ "ML_ADM::K32", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::K32!");
+
+ }
+
+ if (CCTK_EQUALS(K33_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_K33_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_K33_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_K33_bound ,K33_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_K33_bound,
+ "ML_ADM::K33", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::K33!");
+
+ }
+
+ if (CCTK_EQUALS(alpha_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_alpha_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_alpha_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_alpha_bound ,alpha_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_alpha_bound,
+ "ML_ADM::alpha", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::alpha!");
+
+ }
+
+ if (CCTK_EQUALS(g11_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_g11_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g11_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g11_bound ,g11_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g11_bound,
+ "ML_ADM::g11", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::g11!");
+
+ }
+
+ if (CCTK_EQUALS(g21_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_g21_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g21_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g21_bound ,g21_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g21_bound,
+ "ML_ADM::g21", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::g21!");
+
+ }
+
+ if (CCTK_EQUALS(g31_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_g31_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g31_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g31_bound ,g31_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g31_bound,
+ "ML_ADM::g31", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::g31!");
+
+ }
+
+ if (CCTK_EQUALS(g22_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_g22_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g22_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g22_bound ,g22_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g22_bound,
+ "ML_ADM::g22", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::g22!");
+
+ }
+
+ if (CCTK_EQUALS(g32_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_g32_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g32_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g32_bound ,g32_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g32_bound,
+ "ML_ADM::g32", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::g32!");
+
+ }
+
+ if (CCTK_EQUALS(g33_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_g33_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_g33_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_g33_bound ,g33_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_g33_bound,
+ "ML_ADM::g33", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::g33!");
+
+ }
+
+ if (CCTK_EQUALS(beta1_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_beta1_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_beta1_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_beta1_bound ,beta1_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_beta1_bound,
+ "ML_ADM::beta1", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::beta1!");
+
+ }
+
+ if (CCTK_EQUALS(beta2_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_beta2_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_beta2_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_beta2_bound ,beta2_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_beta2_bound,
+ "ML_ADM::beta2", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::beta2!");
+
+ }
+
+ if (CCTK_EQUALS(beta3_bound, "scalar"))
+ {
+ /* apply scalar boundary condition */
+ CCTK_INT handle_beta3_bound = Util_TableCreate(UTIL_TABLE_FLAGS_CASE_INSENSITIVE);
+ if (handle_beta3_bound < 0) CCTK_WARN(-1, "could not create table!");
+ if (Util_TableSetReal(handle_beta3_bound ,beta3_bound_scalar, "SCALAR") < 0)
+ CCTK_WARN(-1, "could not set SCALAR value in table!");
+
+ ierr = Boundary_SelectVarForBC(cctkGH, CCTK_ALL_FACES, 1, handle_beta3_bound,
+ "ML_ADM::beta3", "scalar");
+
+ if (ierr < 0)
+ CCTK_WARN(-1, "Error in registering Scalar BC for ML_ADM::beta3!");
+
+ }
+ return;
+}
+
+
+
+/* template for entries in parameter file:
+#$bound$#ML_ADM::curv_bound = "skip"
+#$bound$#ML_ADM::curv_bound_speed = 1.0
+#$bound$#ML_ADM::curv_bound_limit = 0.0
+#$bound$#ML_ADM::curv_bound_scalar = 0.0
+
+#$bound$#ML_ADM::lapse_bound = "skip"
+#$bound$#ML_ADM::lapse_bound_speed = 1.0
+#$bound$#ML_ADM::lapse_bound_limit = 0.0
+#$bound$#ML_ADM::lapse_bound_scalar = 0.0
+
+#$bound$#ML_ADM::metric_bound = "skip"
+#$bound$#ML_ADM::metric_bound_speed = 1.0
+#$bound$#ML_ADM::metric_bound_limit = 0.0
+#$bound$#ML_ADM::metric_bound_scalar = 0.0
+
+#$bound$#ML_ADM::shift_bound = "skip"
+#$bound$#ML_ADM::shift_bound_speed = 1.0
+#$bound$#ML_ADM::shift_bound_limit = 0.0
+#$bound$#ML_ADM::shift_bound_scalar = 0.0
+
+#$bound$#ML_ADM::K11_bound = "skip"
+#$bound$#ML_ADM::K11_bound_speed = 1.0
+#$bound$#ML_ADM::K11_bound_limit = 0.0
+#$bound$#ML_ADM::K11_bound_scalar = 0.0
+
+#$bound$#ML_ADM::K21_bound = "skip"
+#$bound$#ML_ADM::K21_bound_speed = 1.0
+#$bound$#ML_ADM::K21_bound_limit = 0.0
+#$bound$#ML_ADM::K21_bound_scalar = 0.0
+
+#$bound$#ML_ADM::K31_bound = "skip"
+#$bound$#ML_ADM::K31_bound_speed = 1.0
+#$bound$#ML_ADM::K31_bound_limit = 0.0
+#$bound$#ML_ADM::K31_bound_scalar = 0.0
+
+#$bound$#ML_ADM::K22_bound = "skip"
+#$bound$#ML_ADM::K22_bound_speed = 1.0
+#$bound$#ML_ADM::K22_bound_limit = 0.0
+#$bound$#ML_ADM::K22_bound_scalar = 0.0
+
+#$bound$#ML_ADM::K32_bound = "skip"
+#$bound$#ML_ADM::K32_bound_speed = 1.0
+#$bound$#ML_ADM::K32_bound_limit = 0.0
+#$bound$#ML_ADM::K32_bound_scalar = 0.0
+
+#$bound$#ML_ADM::K33_bound = "skip"
+#$bound$#ML_ADM::K33_bound_speed = 1.0
+#$bound$#ML_ADM::K33_bound_limit = 0.0
+#$bound$#ML_ADM::K33_bound_scalar = 0.0
+
+#$bound$#ML_ADM::alpha_bound = "skip"
+#$bound$#ML_ADM::alpha_bound_speed = 1.0
+#$bound$#ML_ADM::alpha_bound_limit = 0.0
+#$bound$#ML_ADM::alpha_bound_scalar = 0.0
+
+#$bound$#ML_ADM::g11_bound = "skip"
+#$bound$#ML_ADM::g11_bound_speed = 1.0
+#$bound$#ML_ADM::g11_bound_limit = 0.0
+#$bound$#ML_ADM::g11_bound_scalar = 0.0
+
+#$bound$#ML_ADM::g21_bound = "skip"
+#$bound$#ML_ADM::g21_bound_speed = 1.0
+#$bound$#ML_ADM::g21_bound_limit = 0.0
+#$bound$#ML_ADM::g21_bound_scalar = 0.0
+
+#$bound$#ML_ADM::g31_bound = "skip"
+#$bound$#ML_ADM::g31_bound_speed = 1.0
+#$bound$#ML_ADM::g31_bound_limit = 0.0
+#$bound$#ML_ADM::g31_bound_scalar = 0.0
+
+#$bound$#ML_ADM::g22_bound = "skip"
+#$bound$#ML_ADM::g22_bound_speed = 1.0
+#$bound$#ML_ADM::g22_bound_limit = 0.0
+#$bound$#ML_ADM::g22_bound_scalar = 0.0
+
+#$bound$#ML_ADM::g32_bound = "skip"
+#$bound$#ML_ADM::g32_bound_speed = 1.0
+#$bound$#ML_ADM::g32_bound_limit = 0.0
+#$bound$#ML_ADM::g32_bound_scalar = 0.0
+
+#$bound$#ML_ADM::g33_bound = "skip"
+#$bound$#ML_ADM::g33_bound_speed = 1.0
+#$bound$#ML_ADM::g33_bound_limit = 0.0
+#$bound$#ML_ADM::g33_bound_scalar = 0.0
+
+#$bound$#ML_ADM::beta1_bound = "skip"
+#$bound$#ML_ADM::beta1_bound_speed = 1.0
+#$bound$#ML_ADM::beta1_bound_limit = 0.0
+#$bound$#ML_ADM::beta1_bound_scalar = 0.0
+
+#$bound$#ML_ADM::beta2_bound = "skip"
+#$bound$#ML_ADM::beta2_bound_speed = 1.0
+#$bound$#ML_ADM::beta2_bound_limit = 0.0
+#$bound$#ML_ADM::beta2_bound_scalar = 0.0
+
+#$bound$#ML_ADM::beta3_bound = "skip"
+#$bound$#ML_ADM::beta3_bound_speed = 1.0
+#$bound$#ML_ADM::beta3_bound_limit = 0.0
+#$bound$#ML_ADM::beta3_bound_scalar = 0.0
+
+*/
+
diff --git a/ML_ADM/src/Differencing.h b/ML_ADM/src/Differencing.h
new file mode 100644
index 0000000..7ada6a7
--- /dev/null
+++ b/ML_ADM/src/Differencing.h
@@ -0,0 +1,12 @@
+#define PDstandard4th1(u,i,j,k) (p1o12dx*(-8*u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,k)] + 8*u[CCTK_GFINDEX3D(cctkGH,1 + i,j,k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,k)] - u[CCTK_GFINDEX3D(cctkGH,2 + i,j,k)]))
+#define PDstandard4th2(u,i,j,k) (p1o12dy*(-8*u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,k)] + 8*u[CCTK_GFINDEX3D(cctkGH,i,1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,k)] - u[CCTK_GFINDEX3D(cctkGH,i,2 + j,k)]))
+#define PDstandard4th3(u,i,j,k) (p1o12dz*(-8*u[CCTK_GFINDEX3D(cctkGH,i,j,-1 + k)] + 8*u[CCTK_GFINDEX3D(cctkGH,i,j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,j,-2 + k)] - u[CCTK_GFINDEX3D(cctkGH,i,j,2 + k)]))
+#define PDstandard4th11(u,i,j,k) (pm1o12dx2*(-16*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,k)]) + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,k)] + 30*u[CCTK_GFINDEX3D(cctkGH,i,j,k)]))
+#define PDstandard4th22(u,i,j,k) (pm1o12dy2*(-16*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,k)]) + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,k)] + 30*u[CCTK_GFINDEX3D(cctkGH,i,j,k)]))
+#define PDstandard4th33(u,i,j,k) (pm1o12dz2*(-16*(u[CCTK_GFINDEX3D(cctkGH,i,j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,j,1 + k)]) + u[CCTK_GFINDEX3D(cctkGH,i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,j,2 + k)] + 30*u[CCTK_GFINDEX3D(cctkGH,i,j,k)]))
+#define PDstandard4th12(u,i,j,k) (p1o144dxdy*(-64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,-1 + j,k)]) + 64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,-1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,1 + j,k)]) + 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,-1 + j,k)]) - 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,-1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,1 + j,k)]) + u[CCTK_GFINDEX3D(cctkGH,-2 + i,-2 + j,k)] - u[CCTK_GFINDEX3D(cctkGH,-2 + i,2 + j,k)] - u[CCTK_GFINDEX3D(cctkGH,2 + i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,2 + j,k)]))
+#define PDstandard4th13(u,i,j,k) (p1o144dxdz*(-64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,-1 + k)]) + 64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,1 + k)]) + 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,-1 + k)]) - 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,1 + k)]) + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,-2 + k)] - u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,2 + k)] - u[CCTK_GFINDEX3D(cctkGH,2 + i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,2 + k)]))
+#define PDstandard4th21(u,i,j,k) (p1o144dxdy*(-64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,-1 + j,k)]) + 64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,-1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,1 + j,k)]) + 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,-1 + j,k)]) - 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,-1 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,1 + j,k)]) + u[CCTK_GFINDEX3D(cctkGH,-2 + i,-2 + j,k)] - u[CCTK_GFINDEX3D(cctkGH,-2 + i,2 + j,k)] - u[CCTK_GFINDEX3D(cctkGH,2 + i,-2 + j,k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,2 + j,k)]))
+#define PDstandard4th23(u,i,j,k) (p1o144dydz*(-64*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,-1 + k)]) + 64*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,1 + k)]) + 8*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,-1 + k)]) - 8*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,1 + k)]) + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,-2 + k)] - u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,2 + k)] - u[CCTK_GFINDEX3D(cctkGH,i,2 + j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,2 + k)]))
+#define PDstandard4th31(u,i,j,k) (p1o144dxdz*(-64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,-1 + k)]) + 64*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,1 + k)]) + 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,-1 + k)]) - 8*(u[CCTK_GFINDEX3D(cctkGH,-1 + i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,1 + i,j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,1 + k)]) + u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,-2 + k)] - u[CCTK_GFINDEX3D(cctkGH,-2 + i,j,2 + k)] - u[CCTK_GFINDEX3D(cctkGH,2 + i,j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,2 + i,j,2 + k)]))
+#define PDstandard4th32(u,i,j,k) (p1o144dydz*(-64*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,-1 + k)]) + 64*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,1 + k)]) + 8*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,-1 + k)]) - 8*(u[CCTK_GFINDEX3D(cctkGH,i,-1 + j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,1 + j,2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,-1 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,1 + k)]) + u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,-2 + k)] - u[CCTK_GFINDEX3D(cctkGH,i,-2 + j,2 + k)] - u[CCTK_GFINDEX3D(cctkGH,i,2 + j,-2 + k)] + u[CCTK_GFINDEX3D(cctkGH,i,2 + j,2 + k)]))
diff --git a/ML_ADM/src/ML_ADM_Minkowski.c b/ML_ADM/src/ML_ADM_Minkowski.c
new file mode 100644
index 0000000..25d37cf
--- /dev/null
+++ b/ML_ADM/src/ML_ADM_Minkowski.c
@@ -0,0 +1,180 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ML_ADM_Minkowski_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_ADM_Minkowski_Body");
+ }
+
+ if (cctk_iteration % ML_ADM_Minkowski_calc_every != ML_ADM_Minkowski_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alphaL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+
+ /* Assign local copies of grid functions */
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ g11L = 1;
+
+ g21L = 0;
+
+ g31L = 0;
+
+ g22L = 1;
+
+ g32L = 0;
+
+ g33L = 1;
+
+ K11L = 0;
+
+ K21L = 0;
+
+ K31L = 0;
+
+ K22L = 0;
+
+ K32L = 0;
+
+ K33L = 0;
+
+ alphaL = 1;
+
+ beta1L = 0;
+
+ beta2L = 0;
+
+ beta3L = 0;
+
+
+ /* Copy local copies back to grid functions */
+ alpha[index] = alphaL;
+ beta1[index] = beta1L;
+ beta2[index] = beta2L;
+ beta3[index] = beta3L;
+ g11[index] = g11L;
+ g21[index] = g21L;
+ g22[index] = g22L;
+ g31[index] = g31L;
+ g32[index] = g32L;
+ g33[index] = g33L;
+ K11[index] = K11L;
+ K21[index] = K21L;
+ K22[index] = K22L;
+ K31[index] = K31L;
+ K32[index] = K32L;
+ K33[index] = K33L;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ML_ADM_Minkowski(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ML_ADM_Minkowski_Body);
+}
diff --git a/ML_ADM/src/ML_ADM_RHS.c b/ML_ADM/src/ML_ADM_RHS.c
new file mode 100644
index 0000000..46fecca
--- /dev/null
+++ b/ML_ADM/src/ML_ADM_RHS.c
@@ -0,0 +1,536 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ML_ADM_RHS_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_ADM_RHS_Body");
+ }
+
+ if (cctk_iteration % ML_ADM_RHS_calc_every != ML_ADM_RHS_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+ CCTK_REAL detg = INITVALUE;
+ CCTK_REAL G111 = INITVALUE, G121 = INITVALUE, G122 = INITVALUE, G131 = INITVALUE, G132 = INITVALUE, G133 = INITVALUE;
+ CCTK_REAL G211 = INITVALUE, G221 = INITVALUE, G222 = INITVALUE, G231 = INITVALUE, G232 = INITVALUE, G233 = INITVALUE;
+ CCTK_REAL G311 = INITVALUE, G321 = INITVALUE, G322 = INITVALUE, G331 = INITVALUE, G332 = INITVALUE, G333 = INITVALUE;
+ CCTK_REAL gu11 = INITVALUE, gu12 = INITVALUE, gu13 = INITVALUE, gu22 = INITVALUE, gu23 = INITVALUE, gu33 = INITVALUE;
+ CCTK_REAL Km11 = INITVALUE, Km12 = INITVALUE, Km13 = INITVALUE, Km21 = INITVALUE, Km22 = INITVALUE, Km23 = INITVALUE;
+ CCTK_REAL Km31 = INITVALUE, Km32 = INITVALUE, Km33 = INITVALUE;
+ CCTK_REAL R11 = INITVALUE, R21 = INITVALUE, R22 = INITVALUE, R31 = INITVALUE, R32 = INITVALUE, R33 = INITVALUE;
+ CCTK_REAL trK = INITVALUE;
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alphaL = INITVALUE, alpharhsL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta1rhsL = INITVALUE, beta2L = INITVALUE, beta2rhsL = INITVALUE, beta3L = INITVALUE, beta3rhsL = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g11rhsL = INITVALUE, g21L = INITVALUE, g21rhsL = INITVALUE, g22L = INITVALUE, g22rhsL = INITVALUE;
+ CCTK_REAL g31L = INITVALUE, g31rhsL = INITVALUE, g32L = INITVALUE, g32rhsL = INITVALUE, g33L = INITVALUE, g33rhsL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K11rhsL = INITVALUE, K21L = INITVALUE, K21rhsL = INITVALUE, K22L = INITVALUE, K22rhsL = INITVALUE;
+ CCTK_REAL K31L = INITVALUE, K31rhsL = INITVALUE, K32L = INITVALUE, K32rhsL = INITVALUE, K33L = INITVALUE, K33rhsL = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+ CCTK_REAL PDstandard4th1alpha = INITVALUE;
+ CCTK_REAL PDstandard4th2alpha = INITVALUE;
+ CCTK_REAL PDstandard4th3alpha = INITVALUE;
+ CCTK_REAL PDstandard4th11alpha = INITVALUE;
+ CCTK_REAL PDstandard4th22alpha = INITVALUE;
+ CCTK_REAL PDstandard4th33alpha = INITVALUE;
+ CCTK_REAL PDstandard4th12alpha = INITVALUE;
+ CCTK_REAL PDstandard4th13alpha = INITVALUE;
+ CCTK_REAL PDstandard4th23alpha = INITVALUE;
+ CCTK_REAL PDstandard4th1beta1 = INITVALUE;
+ CCTK_REAL PDstandard4th2beta1 = INITVALUE;
+ CCTK_REAL PDstandard4th3beta1 = INITVALUE;
+ CCTK_REAL PDstandard4th1beta2 = INITVALUE;
+ CCTK_REAL PDstandard4th2beta2 = INITVALUE;
+ CCTK_REAL PDstandard4th3beta2 = INITVALUE;
+ CCTK_REAL PDstandard4th1beta3 = INITVALUE;
+ CCTK_REAL PDstandard4th2beta3 = INITVALUE;
+ CCTK_REAL PDstandard4th3beta3 = INITVALUE;
+ CCTK_REAL PDstandard4th1g11 = INITVALUE;
+ CCTK_REAL PDstandard4th2g11 = INITVALUE;
+ CCTK_REAL PDstandard4th3g11 = INITVALUE;
+ CCTK_REAL PDstandard4th22g11 = INITVALUE;
+ CCTK_REAL PDstandard4th33g11 = INITVALUE;
+ CCTK_REAL PDstandard4th12g11 = INITVALUE;
+ CCTK_REAL PDstandard4th13g11 = INITVALUE;
+ CCTK_REAL PDstandard4th21g11 = INITVALUE;
+ CCTK_REAL PDstandard4th23g11 = INITVALUE;
+ CCTK_REAL PDstandard4th31g11 = INITVALUE;
+ CCTK_REAL PDstandard4th32g11 = INITVALUE;
+ CCTK_REAL PDstandard4th1g21 = INITVALUE;
+ CCTK_REAL PDstandard4th2g21 = INITVALUE;
+ CCTK_REAL PDstandard4th3g21 = INITVALUE;
+ CCTK_REAL PDstandard4th33g21 = INITVALUE;
+ CCTK_REAL PDstandard4th12g21 = INITVALUE;
+ CCTK_REAL PDstandard4th13g21 = INITVALUE;
+ CCTK_REAL PDstandard4th21g21 = INITVALUE;
+ CCTK_REAL PDstandard4th23g21 = INITVALUE;
+ CCTK_REAL PDstandard4th31g21 = INITVALUE;
+ CCTK_REAL PDstandard4th32g21 = INITVALUE;
+ CCTK_REAL PDstandard4th1g22 = INITVALUE;
+ CCTK_REAL PDstandard4th2g22 = INITVALUE;
+ CCTK_REAL PDstandard4th3g22 = INITVALUE;
+ CCTK_REAL PDstandard4th11g22 = INITVALUE;
+ CCTK_REAL PDstandard4th33g22 = INITVALUE;
+ CCTK_REAL PDstandard4th12g22 = INITVALUE;
+ CCTK_REAL PDstandard4th13g22 = INITVALUE;
+ CCTK_REAL PDstandard4th21g22 = INITVALUE;
+ CCTK_REAL PDstandard4th23g22 = INITVALUE;
+ CCTK_REAL PDstandard4th31g22 = INITVALUE;
+ CCTK_REAL PDstandard4th32g22 = INITVALUE;
+ CCTK_REAL PDstandard4th1g31 = INITVALUE;
+ CCTK_REAL PDstandard4th2g31 = INITVALUE;
+ CCTK_REAL PDstandard4th3g31 = INITVALUE;
+ CCTK_REAL PDstandard4th22g31 = INITVALUE;
+ CCTK_REAL PDstandard4th12g31 = INITVALUE;
+ CCTK_REAL PDstandard4th13g31 = INITVALUE;
+ CCTK_REAL PDstandard4th21g31 = INITVALUE;
+ CCTK_REAL PDstandard4th23g31 = INITVALUE;
+ CCTK_REAL PDstandard4th31g31 = INITVALUE;
+ CCTK_REAL PDstandard4th32g31 = INITVALUE;
+ CCTK_REAL PDstandard4th1g32 = INITVALUE;
+ CCTK_REAL PDstandard4th2g32 = INITVALUE;
+ CCTK_REAL PDstandard4th3g32 = INITVALUE;
+ CCTK_REAL PDstandard4th11g32 = INITVALUE;
+ CCTK_REAL PDstandard4th12g32 = INITVALUE;
+ CCTK_REAL PDstandard4th13g32 = INITVALUE;
+ CCTK_REAL PDstandard4th21g32 = INITVALUE;
+ CCTK_REAL PDstandard4th23g32 = INITVALUE;
+ CCTK_REAL PDstandard4th31g32 = INITVALUE;
+ CCTK_REAL PDstandard4th32g32 = INITVALUE;
+ CCTK_REAL PDstandard4th1g33 = INITVALUE;
+ CCTK_REAL PDstandard4th2g33 = INITVALUE;
+ CCTK_REAL PDstandard4th3g33 = INITVALUE;
+ CCTK_REAL PDstandard4th11g33 = INITVALUE;
+ CCTK_REAL PDstandard4th22g33 = INITVALUE;
+ CCTK_REAL PDstandard4th12g33 = INITVALUE;
+ CCTK_REAL PDstandard4th13g33 = INITVALUE;
+ CCTK_REAL PDstandard4th21g33 = INITVALUE;
+ CCTK_REAL PDstandard4th23g33 = INITVALUE;
+ CCTK_REAL PDstandard4th31g33 = INITVALUE;
+ CCTK_REAL PDstandard4th32g33 = INITVALUE;
+ CCTK_REAL PDstandard4th1K11 = INITVALUE;
+ CCTK_REAL PDstandard4th2K11 = INITVALUE;
+ CCTK_REAL PDstandard4th3K11 = INITVALUE;
+ CCTK_REAL PDstandard4th1K21 = INITVALUE;
+ CCTK_REAL PDstandard4th2K21 = INITVALUE;
+ CCTK_REAL PDstandard4th3K21 = INITVALUE;
+ CCTK_REAL PDstandard4th1K22 = INITVALUE;
+ CCTK_REAL PDstandard4th2K22 = INITVALUE;
+ CCTK_REAL PDstandard4th3K22 = INITVALUE;
+ CCTK_REAL PDstandard4th1K31 = INITVALUE;
+ CCTK_REAL PDstandard4th2K31 = INITVALUE;
+ CCTK_REAL PDstandard4th3K31 = INITVALUE;
+ CCTK_REAL PDstandard4th1K32 = INITVALUE;
+ CCTK_REAL PDstandard4th2K32 = INITVALUE;
+ CCTK_REAL PDstandard4th3K32 = INITVALUE;
+ CCTK_REAL PDstandard4th1K33 = INITVALUE;
+ CCTK_REAL PDstandard4th2K33 = INITVALUE;
+ CCTK_REAL PDstandard4th3K33 = INITVALUE;
+
+ /* Assign local copies of grid functions */
+ alphaL = alpha[index];
+ beta1L = beta1[index];
+ beta2L = beta2[index];
+ beta3L = beta3[index];
+ g11L = g11[index];
+ g21L = g21[index];
+ g22L = g22[index];
+ g31L = g31[index];
+ g32L = g32[index];
+ g33L = g33[index];
+ K11L = K11[index];
+ K21L = K21[index];
+ K22L = K22[index];
+ K31L = K31[index];
+ K32L = K32[index];
+ K33L = K33[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+ PDstandard4th1alpha = PDstandard4th1(alpha, i, j, k);
+ PDstandard4th2alpha = PDstandard4th2(alpha, i, j, k);
+ PDstandard4th3alpha = PDstandard4th3(alpha, i, j, k);
+ PDstandard4th11alpha = PDstandard4th11(alpha, i, j, k);
+ PDstandard4th22alpha = PDstandard4th22(alpha, i, j, k);
+ PDstandard4th33alpha = PDstandard4th33(alpha, i, j, k);
+ PDstandard4th12alpha = PDstandard4th12(alpha, i, j, k);
+ PDstandard4th13alpha = PDstandard4th13(alpha, i, j, k);
+ PDstandard4th23alpha = PDstandard4th23(alpha, i, j, k);
+ PDstandard4th1beta1 = PDstandard4th1(beta1, i, j, k);
+ PDstandard4th2beta1 = PDstandard4th2(beta1, i, j, k);
+ PDstandard4th3beta1 = PDstandard4th3(beta1, i, j, k);
+ PDstandard4th1beta2 = PDstandard4th1(beta2, i, j, k);
+ PDstandard4th2beta2 = PDstandard4th2(beta2, i, j, k);
+ PDstandard4th3beta2 = PDstandard4th3(beta2, i, j, k);
+ PDstandard4th1beta3 = PDstandard4th1(beta3, i, j, k);
+ PDstandard4th2beta3 = PDstandard4th2(beta3, i, j, k);
+ PDstandard4th3beta3 = PDstandard4th3(beta3, i, j, k);
+ PDstandard4th1g11 = PDstandard4th1(g11, i, j, k);
+ PDstandard4th2g11 = PDstandard4th2(g11, i, j, k);
+ PDstandard4th3g11 = PDstandard4th3(g11, i, j, k);
+ PDstandard4th22g11 = PDstandard4th22(g11, i, j, k);
+ PDstandard4th33g11 = PDstandard4th33(g11, i, j, k);
+ PDstandard4th23g11 = PDstandard4th23(g11, i, j, k);
+ PDstandard4th1g21 = PDstandard4th1(g21, i, j, k);
+ PDstandard4th2g21 = PDstandard4th2(g21, i, j, k);
+ PDstandard4th3g21 = PDstandard4th3(g21, i, j, k);
+ PDstandard4th33g21 = PDstandard4th33(g21, i, j, k);
+ PDstandard4th12g21 = PDstandard4th12(g21, i, j, k);
+ PDstandard4th13g21 = PDstandard4th13(g21, i, j, k);
+ PDstandard4th23g21 = PDstandard4th23(g21, i, j, k);
+ PDstandard4th1g22 = PDstandard4th1(g22, i, j, k);
+ PDstandard4th2g22 = PDstandard4th2(g22, i, j, k);
+ PDstandard4th3g22 = PDstandard4th3(g22, i, j, k);
+ PDstandard4th11g22 = PDstandard4th11(g22, i, j, k);
+ PDstandard4th33g22 = PDstandard4th33(g22, i, j, k);
+ PDstandard4th13g22 = PDstandard4th13(g22, i, j, k);
+ PDstandard4th1g31 = PDstandard4th1(g31, i, j, k);
+ PDstandard4th2g31 = PDstandard4th2(g31, i, j, k);
+ PDstandard4th3g31 = PDstandard4th3(g31, i, j, k);
+ PDstandard4th22g31 = PDstandard4th22(g31, i, j, k);
+ PDstandard4th12g31 = PDstandard4th12(g31, i, j, k);
+ PDstandard4th13g31 = PDstandard4th13(g31, i, j, k);
+ PDstandard4th23g31 = PDstandard4th23(g31, i, j, k);
+ PDstandard4th1g32 = PDstandard4th1(g32, i, j, k);
+ PDstandard4th2g32 = PDstandard4th2(g32, i, j, k);
+ PDstandard4th3g32 = PDstandard4th3(g32, i, j, k);
+ PDstandard4th11g32 = PDstandard4th11(g32, i, j, k);
+ PDstandard4th12g32 = PDstandard4th12(g32, i, j, k);
+ PDstandard4th13g32 = PDstandard4th13(g32, i, j, k);
+ PDstandard4th23g32 = PDstandard4th23(g32, i, j, k);
+ PDstandard4th1g33 = PDstandard4th1(g33, i, j, k);
+ PDstandard4th2g33 = PDstandard4th2(g33, i, j, k);
+ PDstandard4th3g33 = PDstandard4th3(g33, i, j, k);
+ PDstandard4th11g33 = PDstandard4th11(g33, i, j, k);
+ PDstandard4th22g33 = PDstandard4th22(g33, i, j, k);
+ PDstandard4th12g33 = PDstandard4th12(g33, i, j, k);
+ PDstandard4th1K11 = PDstandard4th1(K11, i, j, k);
+ PDstandard4th2K11 = PDstandard4th2(K11, i, j, k);
+ PDstandard4th3K11 = PDstandard4th3(K11, i, j, k);
+ PDstandard4th1K21 = PDstandard4th1(K21, i, j, k);
+ PDstandard4th2K21 = PDstandard4th2(K21, i, j, k);
+ PDstandard4th3K21 = PDstandard4th3(K21, i, j, k);
+ PDstandard4th1K22 = PDstandard4th1(K22, i, j, k);
+ PDstandard4th2K22 = PDstandard4th2(K22, i, j, k);
+ PDstandard4th3K22 = PDstandard4th3(K22, i, j, k);
+ PDstandard4th1K31 = PDstandard4th1(K31, i, j, k);
+ PDstandard4th2K31 = PDstandard4th2(K31, i, j, k);
+ PDstandard4th3K31 = PDstandard4th3(K31, i, j, k);
+ PDstandard4th1K32 = PDstandard4th1(K32, i, j, k);
+ PDstandard4th2K32 = PDstandard4th2(K32, i, j, k);
+ PDstandard4th3K32 = PDstandard4th3(K32, i, j, k);
+ PDstandard4th1K33 = PDstandard4th1(K33, i, j, k);
+ PDstandard4th2K33 = PDstandard4th2(K33, i, j, k);
+ PDstandard4th3K33 = PDstandard4th3(K33, i, j, k);
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ detg = 2*g21L*g31L*g32L + g33L*(g11L*g22L - SQR(g21L)) - g22L*SQR(g31L) - g11L*SQR(g32L);
+
+ gu11 = INV(detg)*(g22L*g33L - SQR(g32L));
+
+ gu12 = (g31L*g32L - g21L*g33L)*INV(detg);
+
+ gu13 = (-(g22L*g31L) + g21L*g32L)*INV(detg);
+
+ gu22 = INV(detg)*(g11L*g33L - SQR(g31L));
+
+ gu23 = (g21L*g31L - g11L*g32L)*INV(detg);
+
+ gu33 = INV(detg)*(g11L*g22L - SQR(g21L));
+
+ G111 = khalf*(gu11*PDstandard4th1g11 + 2*(gu12*PDstandard4th1g21 + gu13*PDstandard4th1g31) - gu12*PDstandard4th2g11 -
+ gu13*PDstandard4th3g11);
+
+ G211 = khalf*(gu12*PDstandard4th1g11 + 2*(gu22*PDstandard4th1g21 + gu23*PDstandard4th1g31) - gu22*PDstandard4th2g11 -
+ gu23*PDstandard4th3g11);
+
+ G311 = khalf*(gu13*PDstandard4th1g11 + 2*(gu23*PDstandard4th1g21 + gu33*PDstandard4th1g31) - gu23*PDstandard4th2g11 -
+ gu33*PDstandard4th3g11);
+
+ G121 = khalf*(gu12*PDstandard4th1g22 + gu11*PDstandard4th2g11 +
+ gu13*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G221 = khalf*(gu22*PDstandard4th1g22 + gu12*PDstandard4th2g11 +
+ gu23*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G321 = khalf*(gu23*PDstandard4th1g22 + gu13*PDstandard4th2g11 +
+ gu33*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G131 = khalf*(gu13*PDstandard4th1g33 + gu11*PDstandard4th3g11 +
+ gu12*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G231 = khalf*(gu23*PDstandard4th1g33 + gu12*PDstandard4th3g11 +
+ gu22*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G331 = khalf*(gu33*PDstandard4th1g33 + gu13*PDstandard4th3g11 +
+ gu23*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G122 = khalf*(gu11*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu12*PDstandard4th2g22 +
+ gu13*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G222 = khalf*(gu12*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu22*PDstandard4th2g22 +
+ gu23*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G322 = khalf*(gu13*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu23*PDstandard4th2g22 +
+ gu33*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G132 = khalf*(gu13*PDstandard4th2g33 + gu11*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu12*PDstandard4th3g22);
+
+ G232 = khalf*(gu23*PDstandard4th2g33 + gu12*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu22*PDstandard4th3g22);
+
+ G332 = khalf*(gu33*PDstandard4th2g33 + gu13*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu23*PDstandard4th3g22);
+
+ G133 = khalf*(-(gu11*PDstandard4th1g33) - gu12*PDstandard4th2g33 + 2*gu11*PDstandard4th3g31 +
+ 2*gu12*PDstandard4th3g32 + gu13*PDstandard4th3g33);
+
+ G233 = khalf*(-(gu12*PDstandard4th1g33) - gu22*PDstandard4th2g33 + 2*gu12*PDstandard4th3g31 +
+ 2*gu22*PDstandard4th3g32 + gu23*PDstandard4th3g33);
+
+ G333 = khalf*(-(gu13*PDstandard4th1g33) - gu23*PDstandard4th2g33 + 2*gu13*PDstandard4th3g31 +
+ 2*gu23*PDstandard4th3g32 + gu33*PDstandard4th3g33);
+
+ R11 = khalf*(4*G231*G321 + G211*(2*G121 - 2*G222 - 2*G332) + G311*(2*G131 - 2*G333) - gu22*PDstandard4th11g22 -
+ 2*(G111*G221 + G232*G311 + G111*G331 + gu23*PDstandard4th11g32) - gu33*PDstandard4th11g33 +
+ 2*gu22*PDstandard4th12g21 + 2*gu23*PDstandard4th12g31 + 2*gu23*PDstandard4th13g21 + 2*gu33*PDstandard4th13g31 -
+ gu22*PDstandard4th22g11 - 2*gu23*PDstandard4th23g11 - gu33*PDstandard4th33g11 + 2*SQR(G221) + 2*SQR(G331));
+
+ R21 = khalf*(2*(G122*G211 + G132*G311 + G231*G322 + G331*G332) -
+ 2*(G121*G221 + G121*G331 + G221*G332 + G321*G333 + gu12*PDstandard4th12g21) - gu23*PDstandard4th12g32 -
+ gu33*PDstandard4th12g33 + gu13*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21) +
+ gu23*PDstandard4th13g22 + gu33*PDstandard4th13g32 + gu12*(PDstandard4th11g22 + PDstandard4th22g11) +
+ gu23*PDstandard4th22g31 + gu13*PDstandard4th23g11 - gu23*PDstandard4th23g21 + gu33*PDstandard4th23g31 -
+ gu33*PDstandard4th33g21);
+
+ R31 = khalf*(2*(G132*G211 + G221*G232 + G133*G311 + G233*G321) -
+ 2*(G222*G231 + G232*G331 + G131*(G221 + G331) + gu13*PDstandard4th13g31) +
+ gu12*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21 + PDstandard4th23g11) +
+ gu22*(PDstandard4th12g32 - PDstandard4th13g22 - PDstandard4th22g31 + PDstandard4th23g21) +
+ gu13*(PDstandard4th11g33 + PDstandard4th33g11) +
+ gu23*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21));
+
+ R22 = khalf*(4*G132*G321 + G122*(2*G221 - 2*(G111 + G331)) + G322*(2*G232 - 2*G333) -
+ 2*(G131*G322 + G222*(G121 + G332) + gu13*PDstandard4th13g22) +
+ gu11*(-PDstandard4th11g22 + 2*PDstandard4th12g21 - PDstandard4th22g11) +
+ gu13*(-2*PDstandard4th22g31 + 2*(PDstandard4th12g32 + PDstandard4th23g21)) +
+ gu33*(-PDstandard4th22g33 + 2*PDstandard4th23g32 - PDstandard4th33g22) + 2*(SQR(G121) + SQR(G332)));
+
+ R32 = khalf*(2*(G121*G131 + G122*G231 + G133*G321 + G233*G322) +
+ gu11*(-PDstandard4th11g32 + PDstandard4th12g31 + PDstandard4th13g21 - PDstandard4th23g11) +
+ gu12*(-PDstandard4th12g32 + PDstandard4th13g22 + PDstandard4th22g31 - PDstandard4th23g21) -
+ 2*(G111*G132 + G131*G332 + G232*(G121 + G332) + gu23*PDstandard4th23g32) +
+ gu13*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21) +
+ gu23*(PDstandard4th22g33 + PDstandard4th33g22));
+
+ R33 = khalf*(4*G132*G231 - gu11*PDstandard4th11g33 -
+ 2*(G111*G133 + G133*G221 + G121*G233 + G222*G233 + G131*G333 + G232*G333 + gu12*PDstandard4th12g33) +
+ 2*(G133*G331 + G233*G332 + gu11*PDstandard4th13g31) + 2*gu12*PDstandard4th13g32 - gu22*PDstandard4th22g33 +
+ 2*gu12*PDstandard4th23g31 + 2*gu22*PDstandard4th23g32 - gu11*PDstandard4th33g11 - 2*gu12*PDstandard4th33g21 -
+ gu22*PDstandard4th33g22 + 2*SQR(G131) + 2*SQR(G232));
+
+ Km11 = gu11*K11L + gu12*K21L + gu13*K31L;
+
+ Km21 = gu12*K11L + gu22*K21L + gu23*K31L;
+
+ Km31 = gu13*K11L + gu23*K21L + gu33*K31L;
+
+ Km12 = gu11*K21L + gu12*K22L + gu13*K32L;
+
+ Km22 = gu12*K21L + gu22*K22L + gu23*K32L;
+
+ Km32 = gu13*K21L + gu23*K22L + gu33*K32L;
+
+ Km13 = gu11*K31L + gu12*K32L + gu13*K33L;
+
+ Km23 = gu12*K31L + gu22*K32L + gu23*K33L;
+
+ Km33 = gu13*K31L + gu23*K32L + gu33*K33L;
+
+ trK = Km11 + Km22 + Km33;
+
+ g11rhsL = -2*alphaL*K11L + 2*(g11L*PDstandard4th1beta1 + g21L*PDstandard4th1beta2 + g31L*PDstandard4th1beta3) +
+ beta1L*PDstandard4th1g11 + beta2L*PDstandard4th2g11 + beta3L*PDstandard4th3g11;
+
+ g21rhsL = -2*alphaL*K21L + g22L*PDstandard4th1beta2 + g32L*PDstandard4th1beta3 + beta1L*PDstandard4th1g21 +
+ g11L*PDstandard4th2beta1 + g21L*(PDstandard4th1beta1 + PDstandard4th2beta2) + g31L*PDstandard4th2beta3 +
+ beta2L*PDstandard4th2g21 + beta3L*PDstandard4th3g21;
+
+ g31rhsL = -2*alphaL*K31L + g32L*PDstandard4th1beta2 + g33L*PDstandard4th1beta3 + beta1L*PDstandard4th1g31 +
+ beta2L*PDstandard4th2g31 + g11L*PDstandard4th3beta1 + g21L*PDstandard4th3beta2 +
+ g31L*(PDstandard4th1beta1 + PDstandard4th3beta3) + beta3L*PDstandard4th3g31;
+
+ g22rhsL = -2*alphaL*K22L + beta1L*PDstandard4th1g22 +
+ 2*(g21L*PDstandard4th2beta1 + g22L*PDstandard4th2beta2 + g32L*PDstandard4th2beta3) + beta2L*PDstandard4th2g22 +
+ beta3L*PDstandard4th3g22;
+
+ g32rhsL = -2*alphaL*K32L + beta1L*PDstandard4th1g32 + g31L*PDstandard4th2beta1 + g33L*PDstandard4th2beta3 +
+ beta2L*PDstandard4th2g32 + g21L*PDstandard4th3beta1 + g22L*PDstandard4th3beta2 +
+ g32L*(PDstandard4th2beta2 + PDstandard4th3beta3) + beta3L*PDstandard4th3g32;
+
+ g33rhsL = -2*alphaL*K33L + beta1L*PDstandard4th1g33 + beta2L*PDstandard4th2g33 +
+ 2*(g31L*PDstandard4th3beta1 + g32L*PDstandard4th3beta2 + g33L*PDstandard4th3beta3) + beta3L*PDstandard4th3g33;
+
+ K11rhsL = -PDstandard4th11alpha + G111*PDstandard4th1alpha +
+ 2*(K11L*PDstandard4th1beta1 + K21L*PDstandard4th1beta2 + K31L*PDstandard4th1beta3) + beta1L*PDstandard4th1K11 +
+ G211*PDstandard4th2alpha + beta2L*PDstandard4th2K11 + G311*PDstandard4th3alpha + beta3L*PDstandard4th3K11 +
+ alphaL*(-2*(K11L*Km11 + K21L*Km21 + K31L*Km31) + R11 + K11L*trK);
+
+ K21rhsL = -PDstandard4th12alpha + G121*PDstandard4th1alpha + K22L*PDstandard4th1beta2 + K32L*PDstandard4th1beta3 +
+ beta1L*PDstandard4th1K21 + G221*PDstandard4th2alpha + K11L*PDstandard4th2beta1 +
+ K21L*(PDstandard4th1beta1 + PDstandard4th2beta2) + K31L*PDstandard4th2beta3 + beta2L*PDstandard4th2K21 +
+ G321*PDstandard4th3alpha + beta3L*PDstandard4th3K21 +
+ alphaL*(-2*(K11L*Km12 + K21L*Km22 + K31L*Km32) + R21 + K21L*trK);
+
+ K31rhsL = -PDstandard4th13alpha + G131*PDstandard4th1alpha + K32L*PDstandard4th1beta2 + K33L*PDstandard4th1beta3 +
+ beta1L*PDstandard4th1K31 + G231*PDstandard4th2alpha + beta2L*PDstandard4th2K31 + G331*PDstandard4th3alpha +
+ K11L*PDstandard4th3beta1 + K21L*PDstandard4th3beta2 + K31L*(PDstandard4th1beta1 + PDstandard4th3beta3) +
+ beta3L*PDstandard4th3K31 + alphaL*(-2*(K11L*Km13 + K21L*Km23 + K31L*Km33) + R31 + K31L*trK);
+
+ K22rhsL = G122*PDstandard4th1alpha + beta1L*PDstandard4th1K22 - PDstandard4th22alpha + G222*PDstandard4th2alpha +
+ 2*(K21L*PDstandard4th2beta1 + K22L*PDstandard4th2beta2 + K32L*PDstandard4th2beta3) + beta2L*PDstandard4th2K22 +
+ G322*PDstandard4th3alpha + beta3L*PDstandard4th3K22 +
+ alphaL*(-2*(K21L*Km12 + K22L*Km22 + K32L*Km32) + R22 + K22L*trK);
+
+ K32rhsL = G132*PDstandard4th1alpha + beta1L*PDstandard4th1K32 - PDstandard4th23alpha + G232*PDstandard4th2alpha +
+ K31L*PDstandard4th2beta1 + K33L*PDstandard4th2beta3 + beta2L*PDstandard4th2K32 + G332*PDstandard4th3alpha +
+ K21L*PDstandard4th3beta1 + K22L*PDstandard4th3beta2 + K32L*(PDstandard4th2beta2 + PDstandard4th3beta3) +
+ beta3L*PDstandard4th3K32 + alphaL*(-2*(K21L*Km13 + K22L*Km23 + K32L*Km33) + R32 + K32L*trK);
+
+ K33rhsL = G133*PDstandard4th1alpha + beta1L*PDstandard4th1K33 + G233*PDstandard4th2alpha + beta2L*PDstandard4th2K33 -
+ PDstandard4th33alpha + G333*PDstandard4th3alpha +
+ 2*(K31L*PDstandard4th3beta1 + K32L*PDstandard4th3beta2 + K33L*PDstandard4th3beta3) + beta3L*PDstandard4th3K33 +
+ alphaL*(-2*(K31L*Km13 + K32L*Km23 + K33L*Km33) + R33 + K33L*trK);
+
+ alpharhsL = 0;
+
+ beta1rhsL = 0;
+
+ beta2rhsL = 0;
+
+ beta3rhsL = 0;
+
+
+ /* Copy local copies back to grid functions */
+ alpharhs[index] = alpharhsL;
+ beta1rhs[index] = beta1rhsL;
+ beta2rhs[index] = beta2rhsL;
+ beta3rhs[index] = beta3rhsL;
+ g11rhs[index] = g11rhsL;
+ g21rhs[index] = g21rhsL;
+ g22rhs[index] = g22rhsL;
+ g31rhs[index] = g31rhsL;
+ g32rhs[index] = g32rhsL;
+ g33rhs[index] = g33rhsL;
+ K11rhs[index] = K11rhsL;
+ K21rhs[index] = K21rhsL;
+ K22rhs[index] = K22rhsL;
+ K31rhs[index] = K31rhsL;
+ K32rhs[index] = K32rhsL;
+ K33rhs[index] = K33rhsL;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ML_ADM_RHS(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ML_ADM_RHS_Body);
+}
diff --git a/ML_ADM/src/ML_ADM_constraints.c b/ML_ADM/src/ML_ADM_constraints.c
new file mode 100644
index 0000000..f3a9296
--- /dev/null
+++ b/ML_ADM/src/ML_ADM_constraints.c
@@ -0,0 +1,445 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ML_ADM_constraints_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_ADM_constraints_Body");
+ }
+
+ if (cctk_iteration % ML_ADM_constraints_calc_every != ML_ADM_constraints_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+ CCTK_REAL detg = INITVALUE;
+ CCTK_REAL G111 = INITVALUE, G121 = INITVALUE, G122 = INITVALUE, G131 = INITVALUE, G132 = INITVALUE, G133 = INITVALUE;
+ CCTK_REAL G211 = INITVALUE, G221 = INITVALUE, G222 = INITVALUE, G231 = INITVALUE, G232 = INITVALUE, G233 = INITVALUE;
+ CCTK_REAL G311 = INITVALUE, G321 = INITVALUE, G322 = INITVALUE, G331 = INITVALUE, G332 = INITVALUE, G333 = INITVALUE;
+ CCTK_REAL gu11 = INITVALUE, gu12 = INITVALUE, gu13 = INITVALUE, gu22 = INITVALUE, gu23 = INITVALUE, gu33 = INITVALUE;
+ CCTK_REAL Km11 = INITVALUE, Km12 = INITVALUE, Km13 = INITVALUE, Km21 = INITVALUE, Km22 = INITVALUE, Km23 = INITVALUE;
+ CCTK_REAL Km31 = INITVALUE, Km32 = INITVALUE, Km33 = INITVALUE;
+ CCTK_REAL R11 = INITVALUE, R21 = INITVALUE, R22 = INITVALUE, R31 = INITVALUE, R32 = INITVALUE, R33 = INITVALUE;
+ CCTK_REAL trK = INITVALUE;
+ CCTK_REAL trR = INITVALUE;
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL HL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ CCTK_REAL M1L = INITVALUE, M2L = INITVALUE, M3L = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+ CCTK_REAL PDstandard4th1g11 = INITVALUE;
+ CCTK_REAL PDstandard4th2g11 = INITVALUE;
+ CCTK_REAL PDstandard4th3g11 = INITVALUE;
+ CCTK_REAL PDstandard4th22g11 = INITVALUE;
+ CCTK_REAL PDstandard4th33g11 = INITVALUE;
+ CCTK_REAL PDstandard4th12g11 = INITVALUE;
+ CCTK_REAL PDstandard4th13g11 = INITVALUE;
+ CCTK_REAL PDstandard4th21g11 = INITVALUE;
+ CCTK_REAL PDstandard4th23g11 = INITVALUE;
+ CCTK_REAL PDstandard4th31g11 = INITVALUE;
+ CCTK_REAL PDstandard4th32g11 = INITVALUE;
+ CCTK_REAL PDstandard4th1g21 = INITVALUE;
+ CCTK_REAL PDstandard4th2g21 = INITVALUE;
+ CCTK_REAL PDstandard4th3g21 = INITVALUE;
+ CCTK_REAL PDstandard4th33g21 = INITVALUE;
+ CCTK_REAL PDstandard4th12g21 = INITVALUE;
+ CCTK_REAL PDstandard4th13g21 = INITVALUE;
+ CCTK_REAL PDstandard4th21g21 = INITVALUE;
+ CCTK_REAL PDstandard4th23g21 = INITVALUE;
+ CCTK_REAL PDstandard4th31g21 = INITVALUE;
+ CCTK_REAL PDstandard4th32g21 = INITVALUE;
+ CCTK_REAL PDstandard4th1g22 = INITVALUE;
+ CCTK_REAL PDstandard4th2g22 = INITVALUE;
+ CCTK_REAL PDstandard4th3g22 = INITVALUE;
+ CCTK_REAL PDstandard4th11g22 = INITVALUE;
+ CCTK_REAL PDstandard4th33g22 = INITVALUE;
+ CCTK_REAL PDstandard4th12g22 = INITVALUE;
+ CCTK_REAL PDstandard4th13g22 = INITVALUE;
+ CCTK_REAL PDstandard4th21g22 = INITVALUE;
+ CCTK_REAL PDstandard4th23g22 = INITVALUE;
+ CCTK_REAL PDstandard4th31g22 = INITVALUE;
+ CCTK_REAL PDstandard4th32g22 = INITVALUE;
+ CCTK_REAL PDstandard4th1g31 = INITVALUE;
+ CCTK_REAL PDstandard4th2g31 = INITVALUE;
+ CCTK_REAL PDstandard4th3g31 = INITVALUE;
+ CCTK_REAL PDstandard4th22g31 = INITVALUE;
+ CCTK_REAL PDstandard4th12g31 = INITVALUE;
+ CCTK_REAL PDstandard4th13g31 = INITVALUE;
+ CCTK_REAL PDstandard4th21g31 = INITVALUE;
+ CCTK_REAL PDstandard4th23g31 = INITVALUE;
+ CCTK_REAL PDstandard4th31g31 = INITVALUE;
+ CCTK_REAL PDstandard4th32g31 = INITVALUE;
+ CCTK_REAL PDstandard4th1g32 = INITVALUE;
+ CCTK_REAL PDstandard4th2g32 = INITVALUE;
+ CCTK_REAL PDstandard4th3g32 = INITVALUE;
+ CCTK_REAL PDstandard4th11g32 = INITVALUE;
+ CCTK_REAL PDstandard4th12g32 = INITVALUE;
+ CCTK_REAL PDstandard4th13g32 = INITVALUE;
+ CCTK_REAL PDstandard4th21g32 = INITVALUE;
+ CCTK_REAL PDstandard4th23g32 = INITVALUE;
+ CCTK_REAL PDstandard4th31g32 = INITVALUE;
+ CCTK_REAL PDstandard4th32g32 = INITVALUE;
+ CCTK_REAL PDstandard4th1g33 = INITVALUE;
+ CCTK_REAL PDstandard4th2g33 = INITVALUE;
+ CCTK_REAL PDstandard4th3g33 = INITVALUE;
+ CCTK_REAL PDstandard4th11g33 = INITVALUE;
+ CCTK_REAL PDstandard4th22g33 = INITVALUE;
+ CCTK_REAL PDstandard4th12g33 = INITVALUE;
+ CCTK_REAL PDstandard4th13g33 = INITVALUE;
+ CCTK_REAL PDstandard4th21g33 = INITVALUE;
+ CCTK_REAL PDstandard4th23g33 = INITVALUE;
+ CCTK_REAL PDstandard4th31g33 = INITVALUE;
+ CCTK_REAL PDstandard4th32g33 = INITVALUE;
+ CCTK_REAL PDstandard4th2K11 = INITVALUE;
+ CCTK_REAL PDstandard4th3K11 = INITVALUE;
+ CCTK_REAL PDstandard4th1K21 = INITVALUE;
+ CCTK_REAL PDstandard4th2K21 = INITVALUE;
+ CCTK_REAL PDstandard4th3K21 = INITVALUE;
+ CCTK_REAL PDstandard4th1K22 = INITVALUE;
+ CCTK_REAL PDstandard4th3K22 = INITVALUE;
+ CCTK_REAL PDstandard4th1K31 = INITVALUE;
+ CCTK_REAL PDstandard4th2K31 = INITVALUE;
+ CCTK_REAL PDstandard4th3K31 = INITVALUE;
+ CCTK_REAL PDstandard4th1K32 = INITVALUE;
+ CCTK_REAL PDstandard4th2K32 = INITVALUE;
+ CCTK_REAL PDstandard4th3K32 = INITVALUE;
+ CCTK_REAL PDstandard4th1K33 = INITVALUE;
+ CCTK_REAL PDstandard4th2K33 = INITVALUE;
+
+ /* Assign local copies of grid functions */
+ g11L = g11[index];
+ g21L = g21[index];
+ g22L = g22[index];
+ g31L = g31[index];
+ g32L = g32[index];
+ g33L = g33[index];
+ K11L = K11[index];
+ K21L = K21[index];
+ K22L = K22[index];
+ K31L = K31[index];
+ K32L = K32[index];
+ K33L = K33[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+ PDstandard4th1g11 = PDstandard4th1(g11, i, j, k);
+ PDstandard4th2g11 = PDstandard4th2(g11, i, j, k);
+ PDstandard4th3g11 = PDstandard4th3(g11, i, j, k);
+ PDstandard4th22g11 = PDstandard4th22(g11, i, j, k);
+ PDstandard4th33g11 = PDstandard4th33(g11, i, j, k);
+ PDstandard4th23g11 = PDstandard4th23(g11, i, j, k);
+ PDstandard4th1g21 = PDstandard4th1(g21, i, j, k);
+ PDstandard4th2g21 = PDstandard4th2(g21, i, j, k);
+ PDstandard4th3g21 = PDstandard4th3(g21, i, j, k);
+ PDstandard4th33g21 = PDstandard4th33(g21, i, j, k);
+ PDstandard4th12g21 = PDstandard4th12(g21, i, j, k);
+ PDstandard4th13g21 = PDstandard4th13(g21, i, j, k);
+ PDstandard4th23g21 = PDstandard4th23(g21, i, j, k);
+ PDstandard4th1g22 = PDstandard4th1(g22, i, j, k);
+ PDstandard4th2g22 = PDstandard4th2(g22, i, j, k);
+ PDstandard4th3g22 = PDstandard4th3(g22, i, j, k);
+ PDstandard4th11g22 = PDstandard4th11(g22, i, j, k);
+ PDstandard4th33g22 = PDstandard4th33(g22, i, j, k);
+ PDstandard4th13g22 = PDstandard4th13(g22, i, j, k);
+ PDstandard4th1g31 = PDstandard4th1(g31, i, j, k);
+ PDstandard4th2g31 = PDstandard4th2(g31, i, j, k);
+ PDstandard4th3g31 = PDstandard4th3(g31, i, j, k);
+ PDstandard4th22g31 = PDstandard4th22(g31, i, j, k);
+ PDstandard4th12g31 = PDstandard4th12(g31, i, j, k);
+ PDstandard4th13g31 = PDstandard4th13(g31, i, j, k);
+ PDstandard4th23g31 = PDstandard4th23(g31, i, j, k);
+ PDstandard4th1g32 = PDstandard4th1(g32, i, j, k);
+ PDstandard4th2g32 = PDstandard4th2(g32, i, j, k);
+ PDstandard4th3g32 = PDstandard4th3(g32, i, j, k);
+ PDstandard4th11g32 = PDstandard4th11(g32, i, j, k);
+ PDstandard4th12g32 = PDstandard4th12(g32, i, j, k);
+ PDstandard4th13g32 = PDstandard4th13(g32, i, j, k);
+ PDstandard4th23g32 = PDstandard4th23(g32, i, j, k);
+ PDstandard4th1g33 = PDstandard4th1(g33, i, j, k);
+ PDstandard4th2g33 = PDstandard4th2(g33, i, j, k);
+ PDstandard4th3g33 = PDstandard4th3(g33, i, j, k);
+ PDstandard4th11g33 = PDstandard4th11(g33, i, j, k);
+ PDstandard4th22g33 = PDstandard4th22(g33, i, j, k);
+ PDstandard4th12g33 = PDstandard4th12(g33, i, j, k);
+ PDstandard4th2K11 = PDstandard4th2(K11, i, j, k);
+ PDstandard4th3K11 = PDstandard4th3(K11, i, j, k);
+ PDstandard4th1K21 = PDstandard4th1(K21, i, j, k);
+ PDstandard4th2K21 = PDstandard4th2(K21, i, j, k);
+ PDstandard4th3K21 = PDstandard4th3(K21, i, j, k);
+ PDstandard4th1K22 = PDstandard4th1(K22, i, j, k);
+ PDstandard4th3K22 = PDstandard4th3(K22, i, j, k);
+ PDstandard4th1K31 = PDstandard4th1(K31, i, j, k);
+ PDstandard4th2K31 = PDstandard4th2(K31, i, j, k);
+ PDstandard4th3K31 = PDstandard4th3(K31, i, j, k);
+ PDstandard4th1K32 = PDstandard4th1(K32, i, j, k);
+ PDstandard4th2K32 = PDstandard4th2(K32, i, j, k);
+ PDstandard4th3K32 = PDstandard4th3(K32, i, j, k);
+ PDstandard4th1K33 = PDstandard4th1(K33, i, j, k);
+ PDstandard4th2K33 = PDstandard4th2(K33, i, j, k);
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ detg = 2*g21L*g31L*g32L + g33L*(g11L*g22L - SQR(g21L)) - g22L*SQR(g31L) - g11L*SQR(g32L);
+
+ gu11 = INV(detg)*(g22L*g33L - SQR(g32L));
+
+ gu12 = (g31L*g32L - g21L*g33L)*INV(detg);
+
+ gu13 = (-(g22L*g31L) + g21L*g32L)*INV(detg);
+
+ gu22 = INV(detg)*(g11L*g33L - SQR(g31L));
+
+ gu23 = (g21L*g31L - g11L*g32L)*INV(detg);
+
+ gu33 = INV(detg)*(g11L*g22L - SQR(g21L));
+
+ G111 = khalf*(gu11*PDstandard4th1g11 + 2*(gu12*PDstandard4th1g21 + gu13*PDstandard4th1g31) - gu12*PDstandard4th2g11 -
+ gu13*PDstandard4th3g11);
+
+ G211 = khalf*(gu12*PDstandard4th1g11 + 2*(gu22*PDstandard4th1g21 + gu23*PDstandard4th1g31) - gu22*PDstandard4th2g11 -
+ gu23*PDstandard4th3g11);
+
+ G311 = khalf*(gu13*PDstandard4th1g11 + 2*(gu23*PDstandard4th1g21 + gu33*PDstandard4th1g31) - gu23*PDstandard4th2g11 -
+ gu33*PDstandard4th3g11);
+
+ G121 = khalf*(gu12*PDstandard4th1g22 + gu11*PDstandard4th2g11 +
+ gu13*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G221 = khalf*(gu22*PDstandard4th1g22 + gu12*PDstandard4th2g11 +
+ gu23*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G321 = khalf*(gu23*PDstandard4th1g22 + gu13*PDstandard4th2g11 +
+ gu33*(PDstandard4th1g32 + PDstandard4th2g31 - PDstandard4th3g21));
+
+ G131 = khalf*(gu13*PDstandard4th1g33 + gu11*PDstandard4th3g11 +
+ gu12*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G231 = khalf*(gu23*PDstandard4th1g33 + gu12*PDstandard4th3g11 +
+ gu22*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G331 = khalf*(gu33*PDstandard4th1g33 + gu13*PDstandard4th3g11 +
+ gu23*(PDstandard4th1g32 - PDstandard4th2g31 + PDstandard4th3g21));
+
+ G122 = khalf*(gu11*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu12*PDstandard4th2g22 +
+ gu13*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G222 = khalf*(gu12*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu22*PDstandard4th2g22 +
+ gu23*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G322 = khalf*(gu13*(-PDstandard4th1g22 + 2*PDstandard4th2g21) + gu23*PDstandard4th2g22 +
+ gu33*(2*PDstandard4th2g32 - PDstandard4th3g22));
+
+ G132 = khalf*(gu13*PDstandard4th2g33 + gu11*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu12*PDstandard4th3g22);
+
+ G232 = khalf*(gu23*PDstandard4th2g33 + gu12*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu22*PDstandard4th3g22);
+
+ G332 = khalf*(gu33*PDstandard4th2g33 + gu13*(-PDstandard4th1g32 + PDstandard4th2g31 + PDstandard4th3g21) +
+ gu23*PDstandard4th3g22);
+
+ G133 = khalf*(-(gu11*PDstandard4th1g33) - gu12*PDstandard4th2g33 + 2*gu11*PDstandard4th3g31 +
+ 2*gu12*PDstandard4th3g32 + gu13*PDstandard4th3g33);
+
+ G233 = khalf*(-(gu12*PDstandard4th1g33) - gu22*PDstandard4th2g33 + 2*gu12*PDstandard4th3g31 +
+ 2*gu22*PDstandard4th3g32 + gu23*PDstandard4th3g33);
+
+ G333 = khalf*(-(gu13*PDstandard4th1g33) - gu23*PDstandard4th2g33 + 2*gu13*PDstandard4th3g31 +
+ 2*gu23*PDstandard4th3g32 + gu33*PDstandard4th3g33);
+
+ R11 = khalf*(-(gu22*PDstandard4th11g22) - 2*(G111*G122 + G111*G133 + G211*G222 + G211*G233 + G311*G322 + G311*G333 +
+ gu23*PDstandard4th11g32) - gu33*PDstandard4th11g33 + 2*gu22*PDstandard4th12g21 + 2*gu23*PDstandard4th12g31 +
+ 2*gu23*PDstandard4th13g21 + 2*gu33*PDstandard4th13g31 - gu22*PDstandard4th22g11 - 2*gu23*PDstandard4th23g11 -
+ gu33*PDstandard4th33g11 + 2*SQR(G121) + 2*SQR(G131) + 2*SQR(G221) + 2*SQR(G231) + 2*SQR(G321) + 2*SQR(G331));
+
+ R21 = khalf*(2*(G131*G132 + G231*G232 + G331*G332) - 2*(G121*G133 + G221*G233 + G321*G333 + gu12*PDstandard4th12g21) -
+ gu23*PDstandard4th12g32 - gu33*PDstandard4th12g33 +
+ gu13*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21) + gu23*PDstandard4th13g22 +
+ gu33*PDstandard4th13g32 + gu12*(PDstandard4th11g22 + PDstandard4th22g11) + gu23*PDstandard4th22g31 +
+ gu13*PDstandard4th23g11 - gu23*PDstandard4th23g21 + gu33*PDstandard4th23g31 - gu33*PDstandard4th33g21);
+
+ R31 = khalf*(2*(G121*G132 + G221*G232 + G321*G332) - 2*(G122*G131 + G222*G231 + G322*G331 + gu13*PDstandard4th13g31) +
+ gu12*(PDstandard4th11g32 - PDstandard4th12g31 - PDstandard4th13g21 + PDstandard4th23g11) +
+ gu22*(PDstandard4th12g32 - PDstandard4th13g22 - PDstandard4th22g31 + PDstandard4th23g21) +
+ gu13*(PDstandard4th11g33 + PDstandard4th33g11) +
+ gu23*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21));
+
+ R22 = khalf*(-2*(G122*(G111 + G133) + G222*(G211 + G233) + G322*(G311 + G333) + gu13*PDstandard4th13g22) +
+ gu11*(-PDstandard4th11g22 + 2*PDstandard4th12g21 - PDstandard4th22g11) +
+ gu13*(-2*PDstandard4th22g31 + 2*(PDstandard4th12g32 + PDstandard4th23g21)) +
+ gu33*(-PDstandard4th22g33 + 2*PDstandard4th23g32 - PDstandard4th33g22) +
+ 2*(SQR(G121) + SQR(G132) + SQR(G221) + SQR(G232) + SQR(G321) + SQR(G332)));
+
+ R32 = khalf*(2*(G121*G131 + G221*G231 + G321*G331) +
+ gu11*(-PDstandard4th11g32 + PDstandard4th12g31 + PDstandard4th13g21 - PDstandard4th23g11) +
+ gu12*(-PDstandard4th12g32 + PDstandard4th13g22 + PDstandard4th22g31 - PDstandard4th23g21) -
+ 2*(G111*G132 + G211*G232 + G311*G332 + gu23*PDstandard4th23g32) +
+ gu13*(PDstandard4th12g33 - PDstandard4th13g32 - PDstandard4th23g31 + PDstandard4th33g21) +
+ gu23*(PDstandard4th22g33 + PDstandard4th33g22));
+
+ R33 = khalf*(gu11*(-PDstandard4th11g33 + 2*PDstandard4th13g31 - PDstandard4th33g11) -
+ 2*((G111 + G122)*G133 + (G211 + G222)*G233 + (G311 + G322)*G333 +
+ gu12*(PDstandard4th12g33 + PDstandard4th33g21)) +
+ gu22*(-PDstandard4th22g33 + 2*PDstandard4th23g32 - PDstandard4th33g22) +
+ 2*(gu12*(PDstandard4th13g32 + PDstandard4th23g31) + SQR(G131) + SQR(G132) + SQR(G231) + SQR(G232) + SQR(G331) +
+ SQR(G332)));
+
+ trR = gu11*R11 + gu22*R22 + 2*(gu12*R21 + gu13*R31 + gu23*R32) + gu33*R33;
+
+ Km11 = gu11*K11L + gu12*K21L + gu13*K31L;
+
+ Km21 = gu12*K11L + gu22*K21L + gu23*K31L;
+
+ Km31 = gu13*K11L + gu23*K21L + gu33*K31L;
+
+ Km12 = gu11*K21L + gu12*K22L + gu13*K32L;
+
+ Km22 = gu12*K21L + gu22*K22L + gu23*K32L;
+
+ Km32 = gu13*K21L + gu23*K22L + gu33*K32L;
+
+ Km13 = gu11*K31L + gu12*K32L + gu13*K33L;
+
+ Km23 = gu12*K31L + gu22*K32L + gu23*K33L;
+
+ Km33 = gu13*K31L + gu23*K32L + gu33*K33L;
+
+ trK = Km11 + Km22 + Km33;
+
+ HL = -2*(Km12*Km21 + Km13*Km31 + Km23*Km32) + trR - SQR(Km11) - SQR(Km22) - SQR(Km33) + SQR(trK);
+
+ M1L = gu12*(-(G121*K11L) + G111*K21L - G221*K21L + G211*K22L - G321*K31L + G311*K32L - PDstandard4th1K21 +
+ PDstandard4th2K11) + gu22*(-(G122*K11L) + G121*K21L - G222*K21L + G221*K22L - G322*K31L + G321*K32L -
+ PDstandard4th1K22 + PDstandard4th2K21) + gu13*
+ (-(G131*K11L) - G231*K21L + G111*K31L - G331*K31L + G211*K32L + G311*K33L - PDstandard4th1K31 + PDstandard4th3K11)\
+ + gu23*(G131*K21L + G231*K22L + G121*K31L + (G221 + G331)*K32L + G321*K33L -
+ 2*(G132*K11L + G232*K21L + G332*K31L + PDstandard4th1K32) + PDstandard4th2K31 + PDstandard4th3K21) +
+ gu33*(-(G133*K11L) - G233*K21L + G131*K31L - G333*K31L + G231*K32L + G331*K33L - PDstandard4th1K33 +
+ PDstandard4th3K31);
+
+ M2L = gu11*(G121*K11L + (-G111 + G221)*K21L - G211*K22L + G321*K31L - G311*K32L + PDstandard4th1K21 -
+ PDstandard4th2K11) + gu12*(G122*K11L + (-G121 + G222)*K21L - G221*K22L + G322*K31L - G321*K32L +
+ PDstandard4th1K22 - PDstandard4th2K21) + gu13*
+ (G132*K11L + (-2*G131 + G232)*K21L + (G121 + G332)*K31L + G221*K32L + G321*K33L + PDstandard4th1K32 -
+ 2*(G231*K22L + G331*K32L + PDstandard4th2K31) + PDstandard4th3K21) +
+ gu23*(-(G132*K21L) - G232*K22L + G122*K31L + G222*K32L - G332*K32L + G322*K33L - PDstandard4th2K32 +
+ PDstandard4th3K22) + gu33*(-(G133*K21L) - G233*K22L + G132*K31L + G232*K32L - G333*K32L + G332*K33L -
+ PDstandard4th2K33 + PDstandard4th3K32);
+
+ M3L = gu11*(G131*K11L + G231*K21L + (-G111 + G331)*K31L - G211*K32L - G311*K33L + PDstandard4th1K31 -
+ PDstandard4th3K11) + gu12*(G132*K11L + (G131 + G232)*K21L + G231*K22L + (-2*G121 + G332)*K31L +
+ (-2*G221 + G331)*K32L + PDstandard4th1K32 + PDstandard4th2K31 - 2*(G321*K33L + PDstandard4th3K21)) +
+ gu22*(G132*K21L + G232*K22L - G122*K31L - G222*K32L + G332*K32L - G322*K33L + PDstandard4th2K32 -
+ PDstandard4th3K22) + gu13*(G133*K11L + G233*K21L + (-G131 + G333)*K31L - G231*K32L - G331*K33L +
+ PDstandard4th1K33 - PDstandard4th3K31) + gu23*
+ (G133*K21L + G233*K22L - G132*K31L - G232*K32L + G333*K32L - G332*K33L + PDstandard4th2K33 - PDstandard4th3K32);
+
+
+ /* Copy local copies back to grid functions */
+ H[index] = HL;
+ M1[index] = M1L;
+ M2[index] = M2L;
+ M3[index] = M3L;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ML_ADM_constraints(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ML_ADM_constraints_Body);
+}
diff --git a/ML_ADM/src/ML_ADM_convertFromADMBase.c b/ML_ADM/src/ML_ADM_convertFromADMBase.c
new file mode 100644
index 0000000..6f47f79
--- /dev/null
+++ b/ML_ADM/src/ML_ADM_convertFromADMBase.c
@@ -0,0 +1,212 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ML_ADM_convertFromADMBase_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_ADM_convertFromADMBase_Body");
+ }
+
+ if (cctk_iteration % ML_ADM_convertFromADMBase_calc_every != ML_ADM_convertFromADMBase_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alpL = INITVALUE;
+ CCTK_REAL alphaL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE;
+ CCTK_REAL betaxL = INITVALUE;
+ CCTK_REAL betayL = INITVALUE;
+ CCTK_REAL betazL = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL gxxL = INITVALUE;
+ CCTK_REAL gxyL = INITVALUE;
+ CCTK_REAL gxzL = INITVALUE;
+ CCTK_REAL gyyL = INITVALUE;
+ CCTK_REAL gyzL = INITVALUE;
+ CCTK_REAL gzzL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ CCTK_REAL kxxL = INITVALUE;
+ CCTK_REAL kxyL = INITVALUE;
+ CCTK_REAL kxzL = INITVALUE;
+ CCTK_REAL kyyL = INITVALUE;
+ CCTK_REAL kyzL = INITVALUE;
+ CCTK_REAL kzzL = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+
+ /* Assign local copies of grid functions */
+ alpL = alp[index];
+ betaxL = betax[index];
+ betayL = betay[index];
+ betazL = betaz[index];
+ gxxL = gxx[index];
+ gxyL = gxy[index];
+ gxzL = gxz[index];
+ gyyL = gyy[index];
+ gyzL = gyz[index];
+ gzzL = gzz[index];
+ kxxL = kxx[index];
+ kxyL = kxy[index];
+ kxzL = kxz[index];
+ kyyL = kyy[index];
+ kyzL = kyz[index];
+ kzzL = kzz[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ g11L = gxxL;
+
+ g21L = gxyL;
+
+ g31L = gxzL;
+
+ g22L = gyyL;
+
+ g32L = gyzL;
+
+ g33L = gzzL;
+
+ K11L = kxxL;
+
+ K21L = kxyL;
+
+ K31L = kxzL;
+
+ K22L = kyyL;
+
+ K32L = kyzL;
+
+ K33L = kzzL;
+
+ alphaL = alpL;
+
+ beta1L = betaxL;
+
+ beta2L = betayL;
+
+ beta3L = betazL;
+
+
+ /* Copy local copies back to grid functions */
+ alpha[index] = alphaL;
+ beta1[index] = beta1L;
+ beta2[index] = beta2L;
+ beta3[index] = beta3L;
+ g11[index] = g11L;
+ g21[index] = g21L;
+ g22[index] = g22L;
+ g31[index] = g31L;
+ g32[index] = g32L;
+ g33[index] = g33L;
+ K11[index] = K11L;
+ K21[index] = K21L;
+ K22[index] = K22L;
+ K31[index] = K31L;
+ K32[index] = K32L;
+ K33[index] = K33L;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ML_ADM_convertFromADMBase(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ML_ADM_convertFromADMBase_Body);
+}
diff --git a/ML_ADM/src/ML_ADM_convertToADMBase.c b/ML_ADM/src/ML_ADM_convertToADMBase.c
new file mode 100644
index 0000000..0ce9dae
--- /dev/null
+++ b/ML_ADM/src/ML_ADM_convertToADMBase.c
@@ -0,0 +1,228 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#define KRANC_C
+
+#include <math.h>
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "GenericFD.h"
+#include "Differencing.h"
+#include "loopcontrol.h"
+
+/* Define macros used in calculations */
+#define INITVALUE (42)
+#define INV(x) ((1.0) / (x))
+#define SQR(x) ((x) * (x))
+#define CUB(x) ((x) * (x) * (x))
+#define QAD(x) ((x) * (x) * (x) * (x))
+
+void ML_ADM_convertToADMBase_Body(cGH *cctkGH, CCTK_INT dir, CCTK_INT face, CCTK_REAL normal[3], CCTK_REAL tangentA[3], CCTK_REAL tangentB[3], CCTK_INT min[3], CCTK_INT max[3], CCTK_INT n_subblock_gfs, CCTK_REAL *subblock_gfs[])
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* Declare finite differencing variables */
+ CCTK_REAL dx = INITVALUE, dy = INITVALUE, dz = INITVALUE;
+ CCTK_REAL dxi = INITVALUE, dyi = INITVALUE, dzi = INITVALUE;
+ CCTK_REAL khalf = INITVALUE, kthird = INITVALUE, ktwothird = INITVALUE, kfourthird = INITVALUE, keightthird = INITVALUE;
+ CCTK_REAL hdxi = INITVALUE, hdyi = INITVALUE, hdzi = INITVALUE;
+
+
+ /* Declare predefined quantities */
+ CCTK_REAL p1o12dx = INITVALUE;
+ CCTK_REAL p1o12dy = INITVALUE;
+ CCTK_REAL p1o12dz = INITVALUE;
+ CCTK_REAL p1o144dxdy = INITVALUE;
+ CCTK_REAL p1o144dxdz = INITVALUE;
+ CCTK_REAL p1o144dydz = INITVALUE;
+ CCTK_REAL pm1o12dx2 = INITVALUE;
+ CCTK_REAL pm1o12dy2 = INITVALUE;
+ CCTK_REAL pm1o12dz2 = INITVALUE;
+
+ if (verbose > 1)
+ {
+ CCTK_VInfo(CCTK_THORNSTRING,"Entering ML_ADM_convertToADMBase_Body");
+ }
+
+ if (cctk_iteration % ML_ADM_convertToADMBase_calc_every != ML_ADM_convertToADMBase_calc_offset)
+ {
+ return;
+ }
+
+ /* Include user-supplied include files */
+
+ /* Initialise finite differencing variables */
+ dx = CCTK_DELTA_SPACE(0);
+ dy = CCTK_DELTA_SPACE(1);
+ dz = CCTK_DELTA_SPACE(2);
+ dxi = 1.0 / dx;
+ dyi = 1.0 / dy;
+ dzi = 1.0 / dz;
+ khalf = 0.5;
+ kthird = 1/3.0;
+ ktwothird = 2.0/3.0;
+ kfourthird = 4.0/3.0;
+ keightthird = 8.0/3.0;
+ hdxi = 0.5 * dxi;
+ hdyi = 0.5 * dyi;
+ hdzi = 0.5 * dzi;
+
+ /* Initialize predefined quantities */
+ p1o12dx = INV(dx)/12.;
+ p1o12dy = INV(dy)/12.;
+ p1o12dz = INV(dz)/12.;
+ p1o144dxdy = (INV(dx)*INV(dy))/144.;
+ p1o144dxdz = (INV(dx)*INV(dz))/144.;
+ p1o144dydz = (INV(dy)*INV(dz))/144.;
+ pm1o12dx2 = -pow(dx,-2)/12.;
+ pm1o12dy2 = -pow(dy,-2)/12.;
+ pm1o12dz2 = -pow(dz,-2)/12.;
+
+ /* Loop over the grid points */
+ LC_LOOP3 (somename,
+ i,j,k, min[0],min[1],min[2], max[0],max[1],max[2],
+ max[0]-min[0],max[1]-min[1],max[2]-min[2])
+ {
+ int index = INITVALUE;
+ int subblock_index = INITVALUE;
+ index = CCTK_GFINDEX3D(cctkGH,i,j,k);
+ subblock_index = i - min[0] + (max[0] - min[0]) * (j - min[1] + (max[1]-min[1]) * (k - min[2]));
+
+ /* Declare shorthands */
+
+ /* Declare local copies of grid functions */
+ CCTK_REAL alpL = INITVALUE;
+ CCTK_REAL alphaL = INITVALUE;
+ CCTK_REAL beta1L = INITVALUE, beta2L = INITVALUE, beta3L = INITVALUE;
+ CCTK_REAL betaxL = INITVALUE;
+ CCTK_REAL betayL = INITVALUE;
+ CCTK_REAL betazL = INITVALUE;
+ CCTK_REAL dtalpL = INITVALUE;
+ CCTK_REAL dtbetaxL = INITVALUE;
+ CCTK_REAL dtbetayL = INITVALUE;
+ CCTK_REAL dtbetazL = INITVALUE;
+ CCTK_REAL g11L = INITVALUE, g21L = INITVALUE, g22L = INITVALUE, g31L = INITVALUE, g32L = INITVALUE, g33L = INITVALUE;
+ CCTK_REAL gxxL = INITVALUE;
+ CCTK_REAL gxyL = INITVALUE;
+ CCTK_REAL gxzL = INITVALUE;
+ CCTK_REAL gyyL = INITVALUE;
+ CCTK_REAL gyzL = INITVALUE;
+ CCTK_REAL gzzL = INITVALUE;
+ CCTK_REAL K11L = INITVALUE, K21L = INITVALUE, K22L = INITVALUE, K31L = INITVALUE, K32L = INITVALUE, K33L = INITVALUE;
+ CCTK_REAL kxxL = INITVALUE;
+ CCTK_REAL kxyL = INITVALUE;
+ CCTK_REAL kxzL = INITVALUE;
+ CCTK_REAL kyyL = INITVALUE;
+ CCTK_REAL kyzL = INITVALUE;
+ CCTK_REAL kzzL = INITVALUE;
+ /* Declare precomputed derivatives*/
+
+ /* Declare derivatives */
+
+ /* Assign local copies of grid functions */
+ alphaL = alpha[index];
+ beta1L = beta1[index];
+ beta2L = beta2[index];
+ beta3L = beta3[index];
+ g11L = g11[index];
+ g21L = g21[index];
+ g22L = g22[index];
+ g31L = g31[index];
+ g32L = g32[index];
+ g33L = g33[index];
+ K11L = K11[index];
+ K21L = K21[index];
+ K22L = K22[index];
+ K31L = K31[index];
+ K32L = K32[index];
+ K33L = K33[index];
+
+ /* Assign local copies of subblock grid functions */
+
+ /* Include user supplied include files */
+
+ /* Precompute derivatives (new style) */
+
+ /* Precompute derivatives (old style) */
+
+ /* Calculate temporaries and grid functions */
+ gxxL = g11L;
+
+ gxyL = g21L;
+
+ gxzL = g31L;
+
+ gyyL = g22L;
+
+ gyzL = g32L;
+
+ gzzL = g33L;
+
+ kxxL = K11L;
+
+ kxyL = K21L;
+
+ kxzL = K31L;
+
+ kyyL = K22L;
+
+ kyzL = K32L;
+
+ kzzL = K33L;
+
+ alpL = alphaL;
+
+ dtalpL = 0;
+
+ betaxL = beta1L;
+
+ betayL = beta2L;
+
+ betazL = beta3L;
+
+ dtbetaxL = 0;
+
+ dtbetayL = 0;
+
+ dtbetazL = 0;
+
+
+ /* Copy local copies back to grid functions */
+ alp[index] = alpL;
+ betax[index] = betaxL;
+ betay[index] = betayL;
+ betaz[index] = betazL;
+ dtalp[index] = dtalpL;
+ dtbetax[index] = dtbetaxL;
+ dtbetay[index] = dtbetayL;
+ dtbetaz[index] = dtbetazL;
+ gxx[index] = gxxL;
+ gxy[index] = gxyL;
+ gxz[index] = gxzL;
+ gyy[index] = gyyL;
+ gyz[index] = gyzL;
+ gzz[index] = gzzL;
+ kxx[index] = kxxL;
+ kxy[index] = kxyL;
+ kxz[index] = kxzL;
+ kyy[index] = kyyL;
+ kyz[index] = kyzL;
+ kzz[index] = kzzL;
+
+ /* Copy local copies back to subblock grid functions */
+ }
+ LC_ENDLOOP3 (somename);
+}
+
+void ML_ADM_convertToADMBase(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ GenericFD_LoopOverEverything(cctkGH, &ML_ADM_convertToADMBase_Body);
+}
diff --git a/ML_ADM/src/RegisterMoL.c b/ML_ADM/src/RegisterMoL.c
new file mode 100644
index 0000000..a87e176
--- /dev/null
+++ b/ML_ADM/src/RegisterMoL.c
@@ -0,0 +1,35 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+
+void ML_ADM_RegisterVars(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+ CCTK_INT ierr = 0;
+
+ /* Register all the evolved grid functions with MoL */
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::K11"), CCTK_VarIndex("ML_ADM::K11rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::K21"), CCTK_VarIndex("ML_ADM::K21rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::K31"), CCTK_VarIndex("ML_ADM::K31rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::K22"), CCTK_VarIndex("ML_ADM::K22rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::K32"), CCTK_VarIndex("ML_ADM::K32rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::K33"), CCTK_VarIndex("ML_ADM::K33rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::alpha"), CCTK_VarIndex("ML_ADM::alpharhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::g11"), CCTK_VarIndex("ML_ADM::g11rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::g21"), CCTK_VarIndex("ML_ADM::g21rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::g31"), CCTK_VarIndex("ML_ADM::g31rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::g22"), CCTK_VarIndex("ML_ADM::g22rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::g32"), CCTK_VarIndex("ML_ADM::g32rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::g33"), CCTK_VarIndex("ML_ADM::g33rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::beta1"), CCTK_VarIndex("ML_ADM::beta1rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::beta2"), CCTK_VarIndex("ML_ADM::beta2rhs"));
+ ierr += MoLRegisterEvolved(CCTK_VarIndex("ML_ADM::beta3"), CCTK_VarIndex("ML_ADM::beta3rhs"));
+ return;
+}
diff --git a/ML_ADM/src/RegisterSymmetries.c b/ML_ADM/src/RegisterSymmetries.c
new file mode 100644
index 0000000..9a2e722
--- /dev/null
+++ b/ML_ADM/src/RegisterSymmetries.c
@@ -0,0 +1,122 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#include "cctk.h"
+#include "cctk_Arguments.h"
+#include "cctk_Parameters.h"
+#include "Symmetry.h"
+
+void ML_ADM_RegisterSymmetries(CCTK_ARGUMENTS)
+{
+ DECLARE_CCTK_ARGUMENTS
+ DECLARE_CCTK_PARAMETERS
+
+
+ /* array holding symmetry definitions */
+ CCTK_INT sym[3];
+
+
+ /* Register symmetries of grid functions */
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::K11");
+
+ sym[0] = -1;
+ sym[1] = -1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::K21");
+
+ sym[0] = -1;
+ sym[1] = 1;
+ sym[2] = -1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::K31");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::K22");
+
+ sym[0] = 1;
+ sym[1] = -1;
+ sym[2] = -1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::K32");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::K33");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::alpha");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::g11");
+
+ sym[0] = -1;
+ sym[1] = -1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::g21");
+
+ sym[0] = -1;
+ sym[1] = 1;
+ sym[2] = -1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::g31");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::g22");
+
+ sym[0] = 1;
+ sym[1] = -1;
+ sym[2] = -1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::g32");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::g33");
+
+ sym[0] = -1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::beta1");
+
+ sym[0] = 1;
+ sym[1] = -1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::beta2");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = -1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::beta3");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::H");
+
+ sym[0] = -1;
+ sym[1] = 1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::M1");
+
+ sym[0] = 1;
+ sym[1] = -1;
+ sym[2] = 1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::M2");
+
+ sym[0] = 1;
+ sym[1] = 1;
+ sym[2] = -1;
+ SetCartSymVN(cctkGH, sym, "ML_ADM::M3");
+
+}
diff --git a/ML_ADM/src/Startup.c b/ML_ADM/src/Startup.c
new file mode 100644
index 0000000..136b76f
--- /dev/null
+++ b/ML_ADM/src/Startup.c
@@ -0,0 +1,13 @@
+/* File produced by user eschnett */
+/* Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007) */
+
+/* Mathematica script written by Ian Hinder and Sascha Husa */
+
+#include "cctk.h"
+
+int ML_ADM_Startup(void)
+{
+ const char * banner = "ML_ADM";
+ CCTK_RegisterBanner(banner);
+ return 0;
+}
diff --git a/ML_ADM/src/make.code.defn b/ML_ADM/src/make.code.defn
new file mode 100644
index 0000000..6e2a3f8
--- /dev/null
+++ b/ML_ADM/src/make.code.defn
@@ -0,0 +1,6 @@
+# File produced by user eschnett
+# Produced with Mathematica Version 6.0 for Mac OS X x86 (32-bit) (April 20, 2007)
+
+# Mathematica script written by Ian Hinder and Sascha Husa
+
+SRCS = Startup.c RegisterMoL.c RegisterSymmetries.c ML_ADM_Minkowski.c ML_ADM_convertFromADMBase.c ML_ADM_RHS.c ML_ADM_convertToADMBase.c ML_ADM_constraints.c Boundaries.c