Improve xTensor support to the point where it can generate an xTensor version of the EM example. There are still several very hack parts which will need to be improved.

1 files changed, 156 insertions, 0 deletions

diff --git a/Examples/EM-xTensor.m b/Examples/EM-xTensor.m
new file mode 100644
index 0000000..944bed3
--- /dev/null
+++ b/Examples/EM-xTensor.m
@@ -0,0 +1,156 @@
+(* ::Package:: *)
+
+(*  Copyright 2004 Sascha Husa, Ian Hinder, Christiane Lechner, Barry Wardell
+
+    This file is part of Kranc.
+
+    Kranc is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    Kranc is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with Foobar; if not, write to the Free Software
+    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+*)
+
+$KrancTensorPackage="xTensorKranc`";
+Get["KrancTensor`"];
+Get["GeneralRelativity`LoadMetric`"];
+$DefInfoQ=False;
+$CVVerbose=False;
+$PrePrint=ScreenDollarIndices;
+
+(**************************************************************************************)
+(* Derivatives *)
+(**************************************************************************************)
+
+derivatives =
+{
+  PDstandard2nd[i_] -> StandardCenteredDifferenceOperator[1,1,i],
+  PDstandard2nd[i_, i_] -> StandardCenteredDifferenceOperator[2,1,i],
+  PDstandard2nd[i_, j_] -> StandardCenteredDifferenceOperator[1,1,i] *
+    StandardCenteredDifferenceOperator[1,1,j],
+
+  PDstandard4th[i_] -> StandardCenteredDifferenceOperator[1,2,i],
+  PDstandard4th[i_, i_] -> StandardCenteredDifferenceOperator[2,2,i],
+  PDstandard4th[i_, j_] -> StandardCenteredDifferenceOperator[1,2,i] *
+    StandardCenteredDifferenceOperator[1,2,j]
+};
+
+(**************************************************************************************)
+(* Tensors *)
+(**************************************************************************************)
+(* Load Euclidean metric and set components of g, epsilon and gamma *)
+LoadMetric["Euclidean"];
+MetricCompute[g, Euclidean, "Christoffel"[1,-1,-1], CVSimplify -> Simplify];
+AllComponentValues[ToBasis[Euclidean][epsilong[-a,e,f]],
+  Table[Signature[{i,j,k}],{i,3},{j,3},{k,3}]];
+RuleToSet[ChristoffelCDPDEuclidean];
+RuleToSet[g];
+RuleToSet[epsilong];
+
+(* Register the tensor quantities with the TensorTools package *)
+DefTensor[El[-a],EuclideanM];
+DefTensor[B[-a],EuclideanM];
+
+(**************************************************************************************)
+(* Groups *)
+(**************************************************************************************)
+
+(* NB This will give B the symmetries of a VECTOR, which is not correct *)
+evolvedGroups = Map[CreateGroupFromTensor, ToBasis[Euclidean][{El[-a], B[-a]}]];
+evaluatedGroups = 
+ {{"constraints", {CEl, CB}},
+  {"endens",{rho}}};
+
+declaredGroups = Join[evolvedGroups, evaluatedGroups];
+declaredGroupNames = Map[First, declaredGroups];
+
+groups = Join[declaredGroups];
+
+(**************************************************************************************)
+(* Initial data *)
+(**************************************************************************************)
+
+initialCalc =
+{
+  Name -> "EM_initial",
+  Schedule -> {"at CCTK_INITIAL"},
+  Equations -> 
+  {
+    El1 -> sigma*Cos[2 Pi (x + y)] ,
+    El2 -> - (1 - sigma)*Cos[2 Pi x] - sigma*Cos[2 Pi (x + y)],
+    El3 -> 0,
+    B1  -> 0,
+    B2  -> 0,
+    B3  -> (1 - sigma)*Cos[2 Pi x] + sigma*Cos[2 Pi (x + y)]
+  }
+};
+
+(**************************************************************************************)
+(* Evolution equations *)
+(**************************************************************************************)
+
+evolCalc =
+{
+  Name -> "EM_evol",
+  Schedule -> {"in MoL_CalcRHS"},
+  Equations -> 
+  {
+    dot[ToBasis[Euclidean][El[-a]]] -> ToBasis[Euclidean][(epsilong[-a,e,f]  CD[-e][B[-f]])],
+    dot[ToBasis[Euclidean][ B[-a]]] -> ToBasis[Euclidean][-(epsilong[-a,e,f] CD[-e][El[-f]])]
+  }
+};
+
+(**************************************************************************************)
+(* Constraint equations *)
+(**************************************************************************************)
+
+constraintsCalc =
+{
+  Name -> "EM_constraints",
+  Equations -> 
+  {
+    CEl -> ToBasis[Euclidean][g[a,b]CD[-b][El[-a]]],
+    CB -> ToBasis[Euclidean][ g[a,b]CD[-b][B[-a]]] 
+  }
+};
+
+(**************************************************************************************)
+(* Energy equation *)
+(**************************************************************************************)
+
+energyCalc =
+{
+  Name -> "EM_energy",
+  Equations -> 
+  {
+    rho -> ToBasis[Euclidean][El[-a] El[a]/2 + B[-a] B[a]/2]
+  }
+};
+
+realParameters = {sigma};
+
+(**************************************************************************************)
+(* Construct the thorn *)
+(**************************************************************************************)
+
+calculations = 
+{
+  initialCalc,
+  evolCalc, 
+  constraintsCalc, 
+  energyCalc
+};
+
+CreateKrancThornTT[groups, ".", "EM-xTensor", 
+  Calculations -> calculations,
+  DeclaredGroups -> declaredGroupNames,
+  PartialDerivatives -> derivatives,
+  RealParameters -> realParameters];