Added preliminary support for boundary conditions

2 files changed, 59 insertions, 9 deletions

diff --git a/Tools/CodeGen/CalculationFunction.m b/Tools/CodeGen/CalculationFunction.m
index 87b12b8..e324c86 100644
--- a/Tools/CodeGen/CalculationFunction.m
+++ b/Tools/CodeGen/CalculationFunction.m
@@ -23,7 +23,8 @@
 BeginPackage["sym`"];
 
 {GridFunctions, Shorthands, Equations, t, DeclarationIncludes,
-LoopPreIncludes, GroupImplementations, PartialDerivatives, Dplus1, Dplus2, Dplus3}
+LoopPreIncludes, GroupImplementations, PartialDerivatives, Dplus1,
+Dplus2, Dplus3, Boundary, Interior, Where}
 
 {INV, SQR, CUB, QAD, dot, pow, exp} 
 
@@ -385,7 +386,7 @@ GrepSyncGroups[x_, func_] := Module[{pick},
 CreateCalculationFunction[calc_, debug_] :=
   Module[{gfs, allSymbols, knownSymbols,
           shorts, eqs, syncGroups, parameters,
-          functionName, dsUsed, groups, pddefs, cleancalc, numeq, eqLoop, GrepSYNC},
+          functionName, dsUsed, groups, pddefs, cleancalc, numeq, eqLoop, GrepSYNC, where},
 
    cleancalc = cleanCalculation[calc];
    cleancalc = cleanCalculation[cleancalc];
@@ -400,7 +401,7 @@ CreateCalculationFunction[calc_, debug_] :=
           parameters = lookupDefault[cleancalc, Parameters, {}];
           groups = lookup[cleancalc, Groups];
           pddefs = lookupDefault[cleancalc, PartialDerivatives, {}];
-
+          where = lookupDefault[cleancalc, Where, Everywhere];
   Print["number of equations in calculation: ", numeq = Length@eqs];
 
   VerifyCalculation[cleancalc];
@@ -494,7 +495,7 @@ CreateCalculationFunction[calc_, debug_] :=
        (* Have removed ability to include external header files here.
           Can be put back when we need it. *)
 
-	eqLoop = Map[equationLoop[#, gfs, shorts, {}, groups, syncGroups, pddefs] &, eqs]};
+	eqLoop = Map[equationLoop[#, gfs, shorts, {}, groups, syncGroups, pddefs, where] &, eqs]};
 
        (* search for SYNCs *)
        If[numeq <= 1,
@@ -562,9 +563,9 @@ checkShorthandAssignmentOrder[eqs_, shorthand_] :=
 
 
 
-equationLoop[eqs_, gfs_, shorts_, incs_, groups_, syncGroups_, pddefs_] :=
+equationLoop[eqs_, gfs_, shorts_, incs_, groups_, syncGroups_, pddefs_, where_] :=
   Module[{rhss, lhss, gfsInRHS, gfsInLHS, localGFs, localMap, eqs2,
-          derivSwitch, actualSyncGroups, code, syncCode},
+          derivSwitch, actualSyncGroups, code, syncCode, loopFunction},
 
     rhss = Map[#[[2]] &, eqs];
     lhss = Map[#[[1]] &, eqs];
@@ -585,7 +586,11 @@ equationLoop[eqs_, gfs_, shorts_, incs_, groups_, syncGroups_, pddefs_] :=
     checkEquationAssignmentOrder[eqs2, shorts];
    code = {InitialiseGridLoopVariables[derivSwitch],
 
-   GridLoop[
+   loopFunction = Switch[where,
+    Boundary, BoundaryLoop,
+    _, GridLoop];
+
+   loopFunction[
    {CommentedBlock["Assign local copies of grid functions",
                    Map[AssignVariable[localName[#], GridName[#]] &, 
                        gfsInRHS]],
diff --git a/Tools/CodeGen/CodeGen.m b/Tools/CodeGen/CodeGen.m
index 7c8a49d..22c53ff 100644
--- a/Tools/CodeGen/CodeGen.m
+++ b/Tools/CodeGen/CodeGen.m
@@ -100,7 +100,7 @@ CommaNewlineSeparated::usage = ""; (* This should not really be in CodeGen *)
 CommaSeparated::usage = "";
 ReplacePowers::usage = "";
 CFormHideStrings::usage = "";
-
+BoundaryLoop::usage = "";
 
 NameRoot::usage = "";
 PartitionVarList::usage = "";
@@ -385,7 +385,7 @@ DeclareGridLoopVariables[] :=
      Map[DeclareVariables[#, "CCTK_INT"] &, 
          {{"i", "j", "k"}, {"istart", "jstart", "kstart"}, 
           {"iend", "jend", "kend"},
-          {"index_offset_x", "index_offset_y", "index_offset_z"}}],
+          {"index_offset_x", "index_offset_y", "index_offset_z", "dir", "face"}}],
 
      If[SOURCELANGUAGE == "C", DeclareVariable["index", "CCTK_INT"], "\n"]
   }];
@@ -467,6 +467,51 @@ GridLoop[block_] :=
        }
         ]]]];
 
+switchOptions[{value_, block_}] :=
+{
+  "case ", value, ":\n", block, "break;\n"
+}
+
+SwitchStatement[var_, pairs__] :=
+{
+  "switch(", var, ")\n",
+  CBlock[{Map[switchOptions, {pairs}]}]
+}
+
+
+
+BoundaryLoop[block_] :=
+  CommentedBlock["Loop over all faces",
+   loopOverInteger["dir", "0", "3",
+     loopOverInteger["face", "0", "2",
+     {
+       DeclareArray["bmin", 3, "int"],
+       DeclareArray["bmax", 3, "int"],
+
+       AssignVariable[arrayElement["bmin", 0], 1],
+       AssignVariable[arrayElement["bmin", 1], 1],
+       AssignVariable[arrayElement["bmin", 2], 1],
+       AssignVariable[arrayElement["bmax", 0], arrayElement["cctk_lsh", 0]],
+       AssignVariable[arrayElement["bmax", 1], arrayElement["cctk_lsh", 1]],
+       AssignVariable[arrayElement["bmax", 2], arrayElement["cctk_lsh", 2]],
+       SwitchStatement["face", 
+        {0,  AssignVariable[arrayElement["bmax", "dir"], "boundary_width"]},
+        {1,  AssignVariable[arrayElement["bmin", "dir"], {arrayElement["cctk_lsh", "dir"], 
+                                                          " - boundary_width" }]}],
+
+       loopOverInteger["k", arrayElement["bmin",2], arrayElement["bmax",2],
+         loopOverInteger["j", arrayElement["bmin",1], arrayElement["bmax",1],
+           loopOverInteger["i", arrayElement["bmin",0], arrayElement["bmax",0],
+
+       { If[SOURCELANGUAGE == "C",  
+            AssignVariable["index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"],
+            ""],
+	 block
+       }
+      
+      ]]]
+     }]]];
+
 conditional[condition_, block_] :=
  {"if (", condition, ")\n",
   CBlock[block]};