Merge branch 'hydro'

Conflicts:
	Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/MathematicaCompat.h
	Examples/kranc.th
	Tools/CodeGen/Kranc.m
	Tools/CodeGen/KrancTensor.m
	Tools/CodeGen/Schedule.m
	Tools/CodeGen/Thorn.m

8 files changed, 303 insertions, 19 deletions

diff --git a/Tools/CodeGen/CalculationFunction.m b/Tools/CodeGen/CalculationFunction.m
index 33bf34c..bdbf040 100644
--- a/Tools/CodeGen/CalculationFunction.m
+++ b/Tools/CodeGen/CalculationFunction.m
@@ -25,6 +25,8 @@ BeginPackage["CalculationFunction`", {"CodeGenCactus`", "CodeGenC`", "CodeGen`",
 
 CreateCalculationFunction::usage = "";
 VerifyCalculation::usage = "";
+calculationSymbols::usage = "";
+GridFunctionsInExpression;
 
 Begin["`Private`"];
 
@@ -811,6 +813,11 @@ DefFn[
 
       If[debugInLoop, Map[InfoVariable[GridName[#]] &, gfsInLHS], ""]}, opts]}]];
 
+(* Unsorted *)
+
+GridFunctionsInExpression[x_, groups_] :=
+  Union[Cases[x, _ ? (MemberQ[allGroupVariables[groups],#] &), Infinity]];
+
 End[];
 
 EndPackage[];
diff --git a/Tools/CodeGen/ConservationCalculation.m b/Tools/CodeGen/ConservationCalculation.m
new file mode 100644
index 0000000..66fc7f3
--- /dev/null
+++ b/Tools/CodeGen/ConservationCalculation.m
@@ -0,0 +1,223 @@
+
+(*  Copyright 2010 Ian Hinder
+
+    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 Kranc; if not, write to the Free Software
+    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
+*)
+
+BeginPackage["ConservationCalculation`", {"Errors`", "Helpers`", "Kranc`",
+  "MapLookup`", "KrancGroups`", "CalculationFunction`", "Differencing`"}];
+
+ProcessConservationCalculation;
+ConservationDifferencingOperators;
+DiffPlus;
+DiffMinus;
+ShiftMinus;
+PDplus;
+ConservationCalculationDeclaredGroups;
+ConservationDifferencingRealParameters;
+hlleAlpha;
+PrimitiveEquations;
+ConservedEquations;
+
+Begin["`Private`"];
+
+ConservationDifferencingOperators[] :=
+{
+  DiffPlus[i_] -> DiffPlusOp[i],
+  DiffMinus[i_] -> DiffMinusOp[i],
+  ShiftMinus[i_] -> 1/shift[i],
+  PDplus[i_] -> DPlus[i]
+};
+
+ConservationDifferencingRealParameters[] :=
+{
+  hlleAlpha
+};
+
+zeroRHSCalc[calc_] :=
+{
+  Name -> lookup[calc,Name] <> "_zero_rhs",
+  Schedule -> {"in MoL_CalcRHS"},
+  Equations ->
+    (Map[First, lookup[calc, Equations]] /. {flux[v_, rest___] :> (dot[v] -> 0)})
+};
+
+minmodVar[v_, i_, vLeft_, vRight_] :=
+{
+  slopeL -> DiffMinus[v, i],
+  slopeR -> DiffPlus[v, i],
+  slope -> MinMod[slopeL, slopeR],
+  vLeft -> v - 0.5 slope,
+  vRight -> v + 0.5 slope
+}
+
+vanLeerVar[v_, i_, vLeft_, vRight_] :=
+{
+  slopeL -> DiffMinus[v, i],
+  slopeR -> DiffPlus[v, i],
+  slope -> VanLeer[slopeL, slopeR],
+  vLeft -> v - 0.5 slope,
+  vRight -> v + 0.5 slope
+}
+
+leftSymbol[v_] :=
+  Symbol["Global`" <> ToString[v] <> "Left"];
+
+rightSymbol[v_] :=
+  Symbol["Global`" <> ToString[v] <> "Right"];
+
+fluxSymbol[v_] :=
+  Symbol["Global`" <> ToString[v] <> "Flux"];
+
+(* Return the list of conserved variables in a calculation *)
+consVars[calc_] :=
+  Union[(Map[First, lookup[calc, Equations]] /. {flux[v_, rest___] :> v})]
+
+(* Return the list of variables to reconstruct in a calculation *)
+primitiveVars[calc_] :=
+  lookup[calc, Primitives];
+
+  (* Module[{allGFs, calcSyms, gfsUsed, conserved, primitive}, *)
+  (*   allGFs = allGroupVariables[lookup[calc, Groups]]; *)
+  (*   calcSyms = calculationSymbols[calc]; *)
+  (*   gfsUsed = Intersection[allGFs, calcSyms]; *)
+  (*   conserved = consVars[calc]; *)
+  (*   primitive = Complement[gfsUsed, conserved]; *)
+  (*   primitive]; *)
+
+(* Return the variables for which Left and Right GFs need to be created *)
+lrGFs[calc_] :=
+  Join[primitiveVars[calc], consVars[calc]];
+
+  (* Module[{allGFs, calcSyms, gfsUsed, conserved, primitive}, *)
+  (*   allGFs = allGroupVariables[lookup[calc, Groups]]; *)
+  (*   calcSyms = calculationSymbols[calc]; *)
+  (*   gfsUsed = Intersection[allGFs, calcSyms]; *)
+  (*   conserved = consVars[calc]; *)
+  (*   primitive = Complement[gfsUsed, conserved]; *)
+  (*   Join[primitive, conserved]]; *)
+
+reconstructCalc[calc_, i_] :=
+{
+  Name -> lookup[calc,Name] <> "_reconstruct_" <> ToString[i],
+  Where -> Interior,
+  Schedule -> {"in MoL_CalcRHS after " <>
+    If[i == 1, lookup[calc,Name] <> "_zero_rhs",
+               lookup[calc,Name] <> "_rhs_" <> ToString[i-1]]},
+  Shorthands -> {slopeL, slopeR, slope},
+  ApplyBCs -> True,
+  Equations ->
+    Flatten[Table[vanLeerVar[v,i, leftSymbol[v], rightSymbol[v]],
+                  {v, primitiveVars[calc]}], 1]
+};
+
+replaceVars[x_, vars_, f_] :=
+  Module[{},
+    x /. Map[(# -> f[#])&, vars]];
+
+hlle[flux[q_, j_] -> frhs_, vars_] :=
+  Module[{},
+  {
+    leftSymbol[fluxSymbol[q]] -> replaceVars[frhs, vars, leftSymbol],
+    rightSymbol[fluxSymbol[q]] -> replaceVars[frhs, vars, Function[v,ShiftMinus[rightSymbol[v],j]]],
+    fluxSymbol[q] ->
+      1/2(leftSymbol[fluxSymbol[q]] + rightSymbol[fluxSymbol[q]] +
+          hlleAlpha(ShiftMinus[rightSymbol[q],j] - leftSymbol[q]))
+  }];
+
+fluxCalc[calc_, i_] :=
+  Module[{fluxes = Select[lookup[calc, Equations], MatchQ[#, flux[_,i]->_] &]},
+  {
+    Name -> lookup[calc,Name] <> "_flux_" <> ToString[i],
+    ApplyBCs -> True,
+    Where -> Interior,
+    Schedule -> {"in MoL_CalcRHS after " <> lookup[calc,Name] <>
+      "_intercell_conserved_" <> ToString[i]},
+    Shorthands -> Join[Map[leftSymbol[fluxSymbol[#]]&, consVars[calc]],
+      Map[rightSymbol[fluxSymbol[#]]&, consVars[calc]]],
+    Equations ->
+      Flatten[Map[hlle[#,GridFunctionsInExpression[#[[2]], lookup[calc, Groups]]] &, fluxes],1]
+  }];
+
+rhs[calc_, i_] :=
+{
+  Name -> lookup[calc,Name] <> "_rhs_" <> ToString[i],
+  Schedule -> {"in MoL_CalcRHS after " <> lookup[calc,Name] <> "_flux_" <> ToString[i]},
+  Where -> Interior,
+  Equations ->
+    Table[dot[v] -> dot[v] - PDplus[fluxSymbol[v], i], {v, consVars[calc]}]
+};
+
+primitivesCalc[calc_, thornName_] :=
+{
+  Name -> lookup[calc, Name] <> "_primitives",
+  Schedule -> {"in MoL_PostStep after " <> thornName <>"_ApplyBCs"},
+  Equations -> lookup[calc, PrimitiveEquations],
+  Shorthands -> lookupDefault[calc, Shorthands, {}]
+};
+
+conservedCalc[calc_] :=
+{
+  Name -> lookup[calc, Name] <> "_conserved",
+  Schedule -> {"at POSTINITIAL"},
+  Equations -> lookup[calc, ConservedEquations],
+  Shorthands -> lookupDefault[calc, Shorthands, {}]
+};
+
+conservedIntercellCalc[calc_, i_] :=
+{
+  Name -> lookup[calc, Name] <> "_intercell_conserved_" <> ToString[i],
+  Schedule -> {"in MoL_CalcRHS after " <> lookup[calc, Name] <> "_reconstruct_" <> ToString[i]},
+
+  Shorthands -> lookupDefault[calc, Shorthands, {}],
+  Equations ->
+    Module[{vars = Join[primitiveVars[calc], consVars[calc]]},
+      Join[lookup[calc, ConservedEquations] /. (Map[# -> leftSymbol[#] &, vars]),
+           lookup[calc, ConservedEquations] /. (Map[# -> rightSymbol[#] &, vars])]]
+};
+
+(* Given a ConservationCalculation structure, return a list of
+   calculations which should be scheduled to solve that conservation
+   law. *)
+ProcessConservationCalculation[calc_, thornName_] :=
+  Module[{},
+  {
+    zeroRHSCalc[calc],
+    conservedCalc[calc],
+    primitivesCalc[calc, thornName],
+    Sequence@@Flatten[
+      Table[
+        {reconstructCalc[calc, i],
+         conservedIntercellCalc[calc, i],
+         fluxCalc[calc, i],
+         rhs[calc, i]}, {i, 1, 3}], 1]
+  }];
+
+(* Return all the new groups which need to be created for this
+   conservation calculation *)
+ConservationCalculationDeclaredGroups[calc_] :=
+  Module[{},
+    Map[CreateGroup[
+      ToString[#]<>"_lr_group",
+      {leftSymbol[#], rightSymbol[#]}, {}] &, lrGFs[calc]] ~Join~
+    Map[CreateGroup[
+      ToString[#]<>"_flux_group",
+      {fluxSymbol[#]}, {}] &, consVars[calc]]];
+
+End[];
+
+EndPackage[];
diff --git a/Tools/CodeGen/Differencing.m b/Tools/CodeGen/Differencing.m
index 499fe81..52ed4b9 100644
--- a/Tools/CodeGen/Differencing.m
+++ b/Tools/CodeGen/Differencing.m
@@ -136,6 +136,8 @@ StandardCenteredDifferenceOperator::usage = "";
 GridFunctionDerivativesInExpression::usage = "";
 DPlus::usage = "";
 DMinus::usage = "";
+DiffPlusOp::usage = "";
+DiffMinusOp::usage = "";
 DZero::usage = "";
 shift::usage = "";
 spacing::usage = "";
@@ -149,6 +151,8 @@ Begin["`Private`"];
 
 DPlus[n_] := (shift[n] - 1)/spacing[n];
 DMinus[n_] := (1 - 1/shift[n])/spacing[n];
+DiffPlusOp[n_] := (shift[n] - 1);
+DiffMinusOp[n_] := (1 - 1/shift[n]);
 DZero[n_] := (DPlus[n] + DMinus[n])/2;
 
 (*************************************************************)
@@ -376,8 +380,8 @@ DefFn[
 
 (*    Print["rhs2 == ", FullForm[rhs]];*)
 
-    pat = Times[spInExpr:(Power[spacing[_],_]..), (Rational[x_,y_])..., rest__];
-(*    Print["pat == ", pat//FullForm];*)
+    pat = (Times[spInExpr:(Power[spacing[_],_]...), (Rational[x_,y_])..., rest__]) | (rest__);
+    (* Print["pat == ", pat//FullForm]; *)
 
     If[MatchQ[rhs, pat],
 (*       Print["matches!"];*)
diff --git a/Tools/CodeGen/Kranc.m b/Tools/CodeGen/Kranc.m
index 4abca9d..3a9cf25 100644
--- a/Tools/CodeGen/Kranc.m
+++ b/Tools/CodeGen/Kranc.m
@@ -39,11 +39,11 @@ dummy;
 LoopPreIncludes, GroupImplementations, PartialDerivatives, NoSimplify,
 Boundary, Interior, InteriorNoSync, Where, AddToStencilWidth,
 Everywhere, normal1, normal2, normal3, INV, SQR, CUB, QAD, dot, pow,
-exp, dt, dx, dy, dz, idx, idy, idz, t}
+exp, dt, dx, dy, dz, idx, idy, idz, t, MinMod, VanLeer}
 
 {ConditionalOnKeyword, ConditionalOnKeywords, CollectList, Interior,
 InteriorNoSync, Boundary, BoundaryWithGhosts, Where, PreDefinitions,
-AllowedSymbols, Parameters, ConditionalOnTextuals};
+AllowedSymbols, Parameters, ConditionalOnTextuals, ApplyBCs};
 
 (* Differencing.m *)
 
@@ -53,6 +53,7 @@ AllowedSymbols, Parameters, ConditionalOnTextuals};
 
 ThornOptions =
  {Calculations -> {},
+  ConservationCalculations -> {},
   DeclaredGroups -> {},
   ODEGroups -> {},
   Implementation -> None,
@@ -93,7 +94,7 @@ Contents, ThornName, BaseImplementation, EvolvedGFs, EvolvedArrays, PrimitiveGFs
 Groups, Calculation, GridFunctions, Shorthands, Equations, Parameter,
 Value, UsesFunctions, ArgString, Conditional, Conditionals, D1, D2, D3, D11, D22,
 D33, D21, D31, D32, Textual, TriggerGroups, Include, RHSGroups, Tags, 
-Steerable, Never, Always, Recover};
+Steerable, Never, Always, Recover, Primitives};
 
 {ExcisionGFs};
 
@@ -101,4 +102,7 @@ Steerable, Never, Always, Recover};
 
 {D1, D2, D3, D11, D22, D33, D21, D31, D32, D12, D13, D23, dot, Eps, Zero3}
 
+(* ConservationCalculation.m *)
+{flux, slopeL, slopeR, slope};
+
 EndPackage[];
diff --git a/Tools/CodeGen/KrancTensor.m b/Tools/CodeGen/KrancTensor.m
index 4480924..ea56d38 100644
--- a/Tools/CodeGen/KrancTensor.m
+++ b/Tools/CodeGen/KrancTensor.m
@@ -23,9 +23,10 @@
 (****************************************************************************)
 (* Wrapper providing tensor support to Kranc (from TensorTools or xTensor)  *)
 (****************************************************************************)
-If[!ValueQ[$KrancTensorPackage], $KrancTensorPackage = "TensorToolsKranc`"];
 
-BeginPackage["KrancTensor`", {"Errors`", "KrancThorn`", "MapLookup`", "KrancGroups`", "Kranc`", $KrancTensorPackage}];
+$KrancTensorPackage = "TensorToolsKranc`";
+
+BeginPackage["KrancTensor`", {"Errors`", "KrancThorn`", "MapLookup`", "KrancGroups`", "Kranc`", $KrancTensorPackage, "ConservationCalculation`", "TensorTools`"}];
 
 CreateKrancThornTT::usage = "Construct a Kranc thorn using tensor expressions.";
 
@@ -38,13 +39,15 @@ Begin["`Private`"];
    -------------------------------------------------------------------------- *)
 
 CreateKrancThornTT[groups_, parentDirectory_, thornName_, opts___] :=
-  Module[{calcs, expCalcs, expGroups, options, derivs, expDerivs, reflectionSymmetries, declaredGroups},
+  Module[{calcs, expCalcs, expGroups, options, derivs, expDerivs, reflectionSymmetries, declaredGroups, consCalcs, expConsCalcs},
     InfoMessage[Terse, "Creating thorn "<>thornName];
     InfoMessage[Terse, "Processing tensorial arguments"];
     calcs = lookup[{opts}, Calculations];
+    consCalcs = lookupDefault[{opts}, ConservationCalculations, {}];
     derivs = lookupDefault[{opts}, PartialDerivatives, {}];
     Map[CheckCalculationTensors, calcs];
     expCalcs = Map[makeCalculationExplicit, calcs];
+    expConsCalcs = Map[makeCalculationExplicit, consCalcs];
 
     InfoMessage[Info, "Group definitions:", groups];
     VerifyGroups[groups];
@@ -52,6 +55,8 @@ CreateKrancThornTT[groups_, parentDirectory_, thornName_, opts___] :=
     expDerivs = Flatten[Map[ExpandComponents,derivs],1];
     expGroups = Map[makeGroupExplicit, groups];
     options = Join[DeleteCases[{opts}, Calculations -> _], {Calculations -> expCalcs}];
+    options = Join[DeleteCases[options, ConservationCalculations -> _],
+      {ConservationCalculations -> expConsCalcs}];
     options = Join[DeleteCases[options, PartialDerivatives -> _], {PartialDerivatives -> expDerivs}];
 
     declaredGroups = lookupDefault[{opts}, DeclaredGroups, {}];
@@ -80,7 +85,10 @@ makeCalculationExplicit[calc_] :=
   mapValueMapMultiple[calc, 
     {Shorthands -> ExpandComponents,
      CollectList -> ExpandComponents,
-     Equations -> ExpandComponents}];
+     Equations -> ExpandComponents,
+     PrimitiveEquations -> MakeExplicit,
+     ConservedEquations -> MakeExplicit,
+     Primitives -> MakeExplicit}];
 
 (* DeleteDuplicates is not available in Mathematica before version 7 *)
 deleteDuplicates[l_] :=
diff --git a/Tools/CodeGen/KrancThorn.m b/Tools/CodeGen/KrancThorn.m
index 16287d7..93760cc 100644
--- a/Tools/CodeGen/KrancThorn.m
+++ b/Tools/CodeGen/KrancThorn.m
@@ -28,7 +28,8 @@
 BeginPackage["KrancThorn`", {"CodeGen`", "Thorn`",
  "MapLookup`", "KrancGroups`", "Differencing`",
  "CalculationFunction`", "Errors`", "Helpers`", "CactusBoundary`",
- "KrancTensor`", "Param`", "Schedule`", "Interface`", "Kranc`", "Jacobian`"}];
+ "KrancTensor`", "Param`", "Schedule`", "Interface`", "Kranc`", "Jacobian`",
+ "ConservationCalculation`"}];
 
 CreateKrancThorn::usage = "Construct a Kranc thorn";
 
@@ -94,7 +95,7 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[
     evolvedODEGroupDefinitions, rhsODEGroupDefinitions, rhsODEGroups,
     allParams, boundarySources, reflectionSymmetries,
     realParamDefs, intParamDefs,
-    pDefs},
+    pDefs, consCalcs, consCalcsIn, consGroups},
 
     (* Parse named arguments *)
 
@@ -112,7 +113,7 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[
     includeFiles = OptionValue[IncludeFiles];
     evolutionTimelevels = OptionValue[EvolutionTimelevels]; (* Redundant *)
     defaultEvolutionTimelevels = lookupDefault[{opts}, DefaultEvolutionTimelevels, evolutionTimelevels];
-    realParams = OptionValue[RealParameters];
+    realParams = OptionValue[RealParameters] ~Join~ ConservationDifferencingRealParameters[];
     intParams = OptionValue[IntParameters];
     realParamDefs = MakeFullParamDefs[realParams];
     intParamDefs = MakeFullParamDefs[intParams];
@@ -123,7 +124,8 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[
     extendedRealParams = OptionValue[ExtendedRealParameters];
     extendedIntParams = OptionValue[ExtendedIntParameters];
     extendedKeywordParams = OptionValue[ExtendedKeywordParameters];
-    partialDerivs = OptionValue[PartialDerivatives];
+    partialDerivs = OptionValue[PartialDerivatives] ~Join~
+      ConservationDifferencingOperators[];
     reflectionSymmetries = OptionValue[ReflectionSymmetries];
 
     coordGroup = {"grid::coordinates", {Kranc`x,Kranc`y,Kranc`z,Kranc`r}};
@@ -151,6 +153,22 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[
     InfoMessage[Terse, "Creating startup file"];
     startup = CreateStartupFile[thornName, thornName];
 
+    consCalcsIn = Append[#,Groups -> groups]& /@
+                    OptionValue[ConservationCalculations];
+
+    (* Add in calculations to solve any conservation laws that have
+       been provided *)
+    calcs = Join[calcs,
+      consCalcs = Flatten[Map[ProcessConservationCalculation[#,thornName] &,
+        consCalcsIn],1]];
+    (* Print["consCalcs = ", consCalcs]; *)
+
+    consGroups = Union@Flatten[
+      Map[ConservationCalculationDeclaredGroups, consCalcsIn],1];
+
+    groups = Join[groups, consGroups];
+    declaredGroups = Join[declaredGroups, Map[groupName, consGroups]];
+
     (* Get the different types of group *)
     evolvedGroups = extractEvolvedGroups[declaredGroups, calcs, groups];
     nonevolvedGroups = extractNonevolvedGroups[declaredGroups, calcs, groups];
diff --git a/Tools/CodeGen/Schedule.m b/Tools/CodeGen/Schedule.m
index c375fd6..87604b1 100644
--- a/Tools/CodeGen/Schedule.m
+++ b/Tools/CodeGen/Schedule.m
@@ -100,8 +100,10 @@ scheduleCalc[calc_, groups_] :=
     relStr = If[before =!= None, " before " <> before, ""]
              <> If[after =!= None, " after " <> after, ""];
 
+    applyBCs = lookupDefault[calc, ApplyBCs, False];
     userSchedule = lookupDefault[calc, Schedule, Automatic];
-    If[userSchedule =!= Automatic,
+
+    If[userSchedule =!= Automatic && !applyBCs,
     Return[Map[
       Join[
       {
@@ -138,7 +140,7 @@ scheduleCalc[calc_, groups_] :=
 
       groupSched = {
         Name               -> "group " <> groupName,
-        SchedulePoint      -> "in MoL_PseudoEvolution" <> relStr,
+        SchedulePoint      -> If[applyBCs, First[userSchedule] <> relStr, "in MoL_PseudoEvolution" <> relStr],
         SynchronizedGroups -> {},
         Language           -> "None",
         Comment            -> lookup[calc, Name]
diff --git a/Tools/CodeGen/Thorn.m b/Tools/CodeGen/Thorn.m
index 5d88067..428d764 100644
--- a/Tools/CodeGen/Thorn.m
+++ b/Tools/CodeGen/Thorn.m
@@ -540,13 +540,16 @@ CreateSetterSource[calcs_, debug_, include_, imp_,
    Symmetries Registration
    ------------------------------------------------------------------------ *)
 
-(* Symmetries registration spec = {{FullName -> "impl::GFname", 
-                                    Sym      -> {symX, symY, symZ}}, ...} *)
+(* Symmetries registration spec = {FullName -> "impl::GFname", 
+                                    Sym      -> {symX, symY, symZ}} *)
 
 SymmetriesBlock[spec_] :=
 
   Module[{i, KrancDim},
 
+  If[!MatchQ[spec, {FullName -> _String, Sym -> {_,_,_}}],
+    ThrowError["SymmetriesBlock: Expecting a symmetry registration spec but got ", spec]];
+
   KrancDim = 3;
 
   sym = lookup[spec, Sym];
@@ -559,8 +562,16 @@ SymmetriesBlock[spec_] :=
 ];
 
 (* syms is a list of rules mapping gridfunctions to their symmetry structures *)
-calcSymmetry[gf_, syms_] := 
-  gf /. syms;
+calcSymmetry[gf_, syms_] :=
+  Module[{},
+    If[mapContains[syms, gf],
+      Return[lookup[syms,gf]],
+      (* FIXME: We are defaulting to scalar symmetries if no information is
+         available.  This shouldn't happen, but I am bypassing this check
+         temporarily. *)
+      Print["WARNING: defaulting to symmetries of a scalar for "<>ToString[gf]];
+      Return[{1,1,1}]]];
+
 
 (* This function guesses the symmetries based on component names as we
    have not been given them *)
@@ -1393,7 +1404,14 @@ CreateThorn[thorn_] :=
                                                lookup[thorn, Sources]];
 
     GenerateFile[sourceDirectory <> "/make.code.defn", lookup[thorn, Makefile]];
+
+    (* Update thorn directory timestamp so that it can be used in makefiles *)
+    GenerateFile[thornDirectory <> "/temp", {}];
+    DeleteFile[thornDirectory <> "/temp"];
+
     Print["Thorn ", thornDirectory, " created successfully"]];
+];
+>>>>>>> hydro
 
 End[];