aboutsummaryrefslogtreecommitdiff
path: root/Tools/CodeGen
diff options
context:
space:
mode:
Diffstat (limited to 'Tools/CodeGen')
-rw-r--r--Tools/CodeGen/CalculationFunction.m267
-rw-r--r--Tools/CodeGen/CodeGenCactus.m147
-rw-r--r--Tools/CodeGen/Kranc.m11
-rw-r--r--Tools/CodeGen/KrancThorn.m2
-rw-r--r--Tools/CodeGen/TensorTools.m5
-rw-r--r--Tools/CodeGen/Thorn.m4
6 files changed, 338 insertions, 98 deletions
diff --git a/Tools/CodeGen/CalculationFunction.m b/Tools/CodeGen/CalculationFunction.m
index 2e8bce0..89025f3 100644
--- a/Tools/CodeGen/CalculationFunction.m
+++ b/Tools/CodeGen/CalculationFunction.m
@@ -251,7 +251,7 @@ simpCollect[collectList_, eqrhs_, localvar_, debug_] :=
(* Return a CodeGen block which assigns dest by evaluating expr *)
assignVariableFromExpression[dest_, expr_, declare_, vectorise_, noSimplify:Boolean : False] :=
- Module[{tSym, type, cleanExpr, code},
+ Module[{type, cleanExpr, code},
type = If[StringMatchQ[ToString[dest], "dir*"], "ptrdiff_t", DataType[]];
cleanExpr = ReplacePowers[expr, vectorise, noSimplify];
@@ -280,7 +280,7 @@ assignVariableFromExpression[dest_, expr_, declare_, vectorise_, noSimplify:Bool
(* This and assignVariableFromExpression should be combined *)
generateCodeFromExpression[expr_, vectorise_, noSimplify:Boolean : False] :=
- Module[{tSym, type, cleanExpr, code},
+ Module[{type, cleanExpr, code},
cleanExpr = ReplacePowers[expr, vectorise, noSimplify];
If[SOURCELANGUAGE == "C",
@@ -378,7 +378,7 @@ DefFn[
shorts, eqs, parameters, parameterRules, odeGroups,
functionName, dsUsed, groups, pddefs, cleancalc, eqLoop, where,
addToStencilWidth, pDefs, haveCondTextuals, condTextuals, calc,
- kernelCall,debug,imp,gridName},
+ kernelCall, DGFEDefs, DGFECall,debug,imp,gridName},
debug = OptionValue[Debug];
imp = lookup[calcp, Implementation];
@@ -414,9 +414,9 @@ DefFn[
gfs = allGroupVariables[groups];
- InfoMessage[InfoFull, " ", Length@shorts, " shorthands"];
- InfoMessage[InfoFull, " ", Length@gfs, " grid functions"];
- InfoMessage[InfoFull, " ", Length@groups, " groups"];
+ InfoMessage[InfoFull, " " <> ToString[Length@shorts] <> " shorthands"];
+ InfoMessage[InfoFull, " " <> ToString[Length@gfs] <> " grid functions"];
+ InfoMessage[InfoFull, " " <> ToString[Length@groups] <> " groups"];
InfoMessage[InfoFull, "Shorthands: ", shorts];
InfoMessage[InfoFull, "Grid functions: ", gfs];
@@ -484,6 +484,206 @@ DefFn[
ThrowError["Unknown 'Where' entry in calculation " <>
functionName <> ": " <> ToString[where]]];
+ DGFEDefs =
+ If[OptionValue[UseDGFE] && lookupDefault[cleancalc, UseDGFE, False],
+ Module[
+ {name, lhss, gfsInLHS, vars},
+ InfoMessage[InfoFull, "Generating DGFE boilerplate"];
+ name = lookup[cleancalc, Name];
+ lhss = Map[#[[1]] &, eqs];
+ gfsInLHS = Union[Cases[lhss, _ ? (MemberQ[gfs,#] &), Infinity]];
+ InfoMessage[InfoFull, "gfsInLHS:" <> Map[" "<>ToString[#] &, gfsInLHS]];
+ (* TODO: do this in a better way, don't examine the variable names *)
+ vars = Select[gfsInLHS,
+ StringMatchQ[ToString[#],
+ RegularExpression[".*rhs.*"]] &];
+ vars = Map[Symbol[StringReplace[ToString[#], "rhs"->""]] &, vars];
+ InfoMessage[InfoFull, "DGFE variables:" <> Map[" "<>ToString[#] &, vars]];
+ {
+ "",
+ "",
+ "",
+ "/* DGFE Definitions */",
+ "",
+ "#include <hrscc.hh>",
+ "",
+ "/*** WENO reconstruction ***/", (* TODO: do we need this? *)
+ "",
+ "/* Configuration */",
+ "#define config_eno_width 3",
+ "#define config_weno_type hrscc::policy::standard",
+ "#define config_weno_optim hrscc::policy::order",
+ "#define config_weno_limiter hrscc::policy::dummy",
+ "#define config_tvd_limiter hrscc::limiters::minmod",
+ "",
+ "/* Export definitions */",
+ "#define "<>name<>"_eno_stencil hrscc::ENOStencil<config_eno_width>",
+ "#define "<>name<>"_weno_weights hrscc::WENOWeights<config_eno_width>",
+ "#define "<>name<>"_weno_stencil hrscc::WENOStencil<config_eno_width, config_weno_type, config_weno_optim>",
+ "#define "<>name<>"_weno_limiter hrscc::WENOLimiter<config_eno_width, "<>name<>"_weno_stencil::width, config_weno_limiter>",
+ "#define "<>name<>"_weno_reconstruction hrscc::WENOReconstruction<"<>name<>"_eno_stencil, "<>name<>"_weno_limiter, "<>name<>"_weno_stencil, "<>name<>"_weno_weights>",
+ "#define "<>name<>"_tvd_reconstruction hrscc::TVDReconstruction<config_tvd_limiter>",
+ "#define "<>name<>"_mp5_reconstruction hrscc::MP5Reconstruction",
+ "#define "<>name<>"_upwind_reconstruction hrscc::UpwindReconstruction",
+ "#define "<>name<>"_reconstruction "<>name<>"_weno_reconstruction",
+ "",
+ "/*** Flux splitting ***/", (* TODO: do we need this? *)
+ "",
+ "/* Configuration */",
+ "#define config_entropy_fix hrscc::NoEntropyFix",
+ "#define config_flux_split hrscc::RoeFS<"<>name<>"_reconstruction, config_entropy_fix>",
+ "#define config_system_split hrscc::ComponentSplit<DGFE_"<>name<>", config_flux_split>",
+ "",
+ "/* Export definitions */",
+ "#define "<>name<>"_flux_split config_flux_split",
+ "#define "<>name<>"_system_split config_system_split",
+ "#define "<>name<>"_fd_method hrscc::FiniteDifference<DGFE_"<>name<>", "<>name<>"_system_split>",
+ "",
+ "/*** SDG method ***/", (* TODO: do we need this? *)
+ "",
+ "/* Configuration */",
+ "#define config_sdg_order 5", (* TODO: make this a parameter *)
+ "#define config_riemann_solver hrscc::LaxFriedrichsRS<DGFE_"<>name<>", false>",
+ "",
+ "/* Export definitions */",
+ "#define "<>name<>"_sdg_grid hrscc::GNIGrid<hrscc::GLLElement<config_sdg_order> >",
+ "#define "<>name<>"_sdg_method hrscc::SDGMethod<DGFE_"<>name<>", "<>name<>"_sdg_grid, config_riemann_solver>",
+ "",
+ "/*** Numerical scheme ***/",
+ "",
+ "/* Configuration */",
+ "#define config_method "<>name<>"_sdg_method",
+ "",
+ "/* Export definitions */",
+ "#define "<>name<>"_method config_method",
+ "#define "<>name<>"_solver hrscc::CLawSolver<DGFE_"<>name<>", config_method>",
+ "",
+ "",
+ "",
+ "class DGFE_"<>name<>";",
+ "",
+ "namespace hrscc {",
+ " template<>",
+ " struct traits<DGFE_"<>name<>"> {",
+ " // All state vector variables",
+ " enum state_t {" <> Map["i"<>ToString[#]<>", " &, vars] <> "nvars};",
+ " enum {nequations = nvars};",
+ " enum {nexternal = 0};",
+ " enum {nbitmasks = 0};",
+ " static bool const pure = false;",
+ " };",
+ "} // namespace",
+ "",
+ "",
+ "",
+ "class DGFE_"<>name<>": public hrscc::CLaw<DGFE_"<>name<>"> {",
+ "public:",
+ " typedef hrscc::CLaw<DGFE_"<>name<>"> claw;",
+ " typedef hrscc::traits<DGFE_"<>name<>">::state_t state_t;",
+ " typedef hrscc::traits<DGFE_"<>name<>"> variables_t;",
+ " static int const nvars = variables_t::nvars;",
+ " ",
+ " DGFE_"<>name<>"();",
+ " ",
+ " inline void prim_to_all(hrscc::Observer<claw> & observer) const",
+ " {",
+ " }",
+ " ",
+ " template<hrscc::policy::direction_t dir>",
+ " inline void fluxes(hrscc::Observer<claw> & observer) const",
+ " {",
+ " cGH const *const cctkGH = observer.cctkGH;",
+ " DECLARE_CCTK_ARGUMENTS;",
+ " ",
+ " int const i = observer.index[0];",
+ " int const j = observer.index[1];",
+ " int const k = observer.index[2];",
+ " int const ind3d = CCTK_GFINDEX3D(cctkGH, i,j,k);",
+ " ",
+ Map[" CCTK_REAL flux"<>ToString[#]<>"L;" &, vars],
+ " ",
+ " switch (dir) {",
+ Table[{
+ " case hrscc::policy::" <> {"x", "y", "z"}[[dir]] <> ": {",
+ Map[" flux"<>ToString[#]<>"L = flux"<>ToString[#]<>ToString[dir]<>"[ind3d];" &, vars],
+ " break;",
+ " }"},
+ {dir, 1, 3}],
+ " default:",
+ " assert(0);",
+ " }",
+ " ",
+ Map[" observer.flux[dir][variables_t::i"<>ToString[#]<>"] = - flux"<>ToString[#]<>"L;" &, vars],
+ " }",
+ " ",
+ " template<hrscc::policy::direction_t dir>",
+ " inline void eigenvalues(hrscc::Observer<claw> & observer) const",
+ " {",
+ " assert(0);",
+ " }",
+ " ",
+ " template<hrscc::policy::direction_t dir>",
+ " inline void eig(hrscc::Observer<claw> & observer) const",
+ " {",
+ " assert(0);",
+ " }",
+ "};",
+ "",
+ "",
+ "",
+ "namespace hrscc {",
+ " template<> int CLaw<DGFE_"<>name<>">::conserved_idx[DGFE_"<>name<>"::nvars] = {};",
+ " template<> int CLaw<DGFE_"<>name<>">::primitive_idx[DGFE_"<>name<>"::nvars] = {};",
+ " template<> int CLaw<DGFE_"<>name<>">::rhs_idx[DGFE_"<>name<>"::nvars] = {};",
+ " template<> int CLaw<DGFE_"<>name<>">::field_idx[0] = {};",
+ " template<> int CLaw<DGFE_"<>name<>">::bitmask_idx[0] = {};",
+
+ "} // namespace hrscc",
+ "",
+ "",
+ "",
+ "namespace {",
+ " int varindex(char const* const varname)",
+ " {",
+ " int const vi = CCTK_VarIndex(varname);",
+ " if (vi<0) CCTK_WARN(CCTK_WARN_ABORT, \"Internal error\");",
+ " return vi;",
+ " }",
+ "}",
+ "",
+ "DGFE_"<>name<>"::DGFE_"<>name<>"()",
+ "{",
+ " using namespace hrscc;",
+ " ",
+ Map[" CLaw<DGFE_"<>name<>">::conserved_idx[variables_t::i"<>ToString[#]<>"] = varindex(CCTK_THORNSTRING \"::"<>ToString[#]<>"\");" &, vars],
+ Map[" CLaw<DGFE_"<>name<>">::primitive_idx[variables_t::i"<>ToString[#]<>"] = varindex(CCTK_THORNSTRING \"::"<>ToString[#]<>"\");" &, vars],
+ " ",
+ Map[" CLaw<DGFE_"<>name<>">::rhs_idx[variables_t::i"<>ToString[#]<>"] = varindex(CCTK_THORNSTRING \"::"<>ToString[#]<>"rhs\");" &, vars],
+ "}",
+ "",
+ "",
+ ""
+ } // Flatten // Map[# <> "\n" &, #] &],
+ {}
+ ];
+
+ DGFECall =
+ If[OptionValue[UseDGFE] && lookupDefault[cleancalc, UseDGFE, False],
+ Module[
+ {name},
+ name = lookup[cleancalc, Name];
+ {
+ "",
+ "{",
+ " static "<>name<>"_solver *solver = NULL;",
+ " if (not solver) solver = new "<>name<>"_method(cctkGH);",
+ " solver->compute_rhs();",
+ " // delete solver;",
+ "}"
+ } // Flatten // Map[# <> "\n" &, #] &],
+ {}
+ ];
+
InfoMessage[InfoFull,"Generating function"];
{
lookup[calcp,BodyFunction][{
@@ -491,7 +691,7 @@ DefFn[
(* We could (or probably should) write this into a source file of its own *)
If[OptionValue[UseOpenCL],
{
- "char const * const source =\n"
+ "char const *const source =\n"
},
{
}],
@@ -532,16 +732,16 @@ DefFn[
groupNames = groupsInCalculation[cleancalc, imp];
groupNames = Select[groupNames, !MemberQ[ignoreGroups, #] &];
{
- "char const * const groups[] = {",
+ "char const *const groups[] = {",
Riffle[Join[Map[Quote, groupNames], {"NULL"}], ","],
"};\n\n"
}
],
- "static struct OpenCLKernel * kernel = NULL;\n",
- "char const * const sources[] = {differencing, source, NULL};\n",
- "OpenCLRunTime_CallKernel (cctkGH, CCTK_THORNSTRING, \"" <> functionName <> "\",\n",
- " sources, groups, NULL, NULL, NULL, -1,\n",
- " imin, imax, &kernel);\n\n"
+ "static struct OpenCLKernel *kernel = NULL;\n",
+ "char const *const sources[] = {differencing, source, NULL};\n",
+ "OpenCLRunTime_CallKernel(cctkGH, CCTK_THORNSTRING, \"" <> functionName <> "\",\n",
+ " sources, groups, NULL, NULL, NULL, -1,\n",
+ " imin, imax, &kernel);\n\n"
},
{
}]
@@ -576,6 +776,8 @@ DefFn[
kernelCall,
+ DGFECall,
+
ConditionalOnParameterTextual["verbose > 1",
"CCTK_VInfo(CCTK_THORNSTRING,\"Leaving " <> bodyFunctionName <> "\");\n"]
}]],
@@ -777,38 +979,13 @@ DefFn[
Map[InfoVariable[#[[1]]] &, (eqs2 /. localMap)],
""],
- Which[OptionValue[UseOpenCL],
- CommentedBlock["Copy local copies back to grid functions",
- Map[StorePartialVariableInLoop[gridName[#], localName[#]] &, gfsInLHS]],
- OptionValue[UseVectors],
- {
- CommentedBlock["If necessary, store only partial vectors after the first iteration",
- ConditionalOnParameterTextual["CCTK_REAL_VEC_SIZE > 2 && CCTK_BUILTIN_EXPECT(i < lc_imin && i+CCTK_REAL_VEC_SIZE > lc_imax, 0)",
- {
- DeclareAssignVariable["ptrdiff_t", "elt_count_lo", "lc_imin-i"],
- DeclareAssignVariable["ptrdiff_t", "elt_count_hi", "lc_imax-i"],
- Map[StoreMiddlePartialVariableInLoop[gridName[#], localName[#], "elt_count_lo", "elt_count_hi"] &, gfsInLHS],
- "break;\n"
- }]],
- CommentedBlock["If necessary, store only partial vectors after the first iteration",
- ConditionalOnParameterTextual["CCTK_REAL_VEC_SIZE > 1 && CCTK_BUILTIN_EXPECT(i < lc_imin, 0)",
- {
- DeclareAssignVariable["ptrdiff_t", "elt_count", "lc_imin-i"],
- Map[StoreHighPartialVariableInLoop[gridName[#], localName[#], "elt_count"] &, gfsInLHS],
- "continue;\n"
- }]],
- CommentedBlock["If necessary, store only partial vectors after the last iteration",
- ConditionalOnParameterTextual["CCTK_REAL_VEC_SIZE > 1 && CCTK_BUILTIN_EXPECT(i+CCTK_REAL_VEC_SIZE > lc_imax, 0)",
- {
- DeclareAssignVariable["ptrdiff_t", "elt_count", "lc_imax-i"],
- Map[StoreLowPartialVariableInLoop[gridName[#], localName[#], "elt_count"] &, gfsInLHS],
- "break;\n"
- }]],
- Map[StoreVariableInLoop[gridName[#], localName[#]] &, gfsInLHS]
- },
- True,
- CommentedBlock["Copy local copies back to grid functions",
- Map[AssignVariableInLoop[gridName[#], localName[#]] &, gfsInLHS]]],
+ If[OptionValue[UseVectors] || OptionValue[UseOpenCL],
+ CommentedBlock["Copy local copies back to grid functions",
+ { PrepareStorePartialVariableInLoop["i", "lc_imin", "lc_imax"],
+ Map[StorePartialVariableInLoop[gridName[#], localName[#]] &,
+ gfsInLHS] }],
+ CommentedBlock["Copy local copies back to grid functions",
+ Map[AssignVariableInLoop[gridName[#], localName[#]] &, gfsInLHS]]],
If[debugInLoop, Map[InfoVariable[gridName[#]] &, gfsInLHS], ""]}, opts]}]];
diff --git a/Tools/CodeGen/CodeGenCactus.m b/Tools/CodeGen/CodeGenCactus.m
index 8006774..f49e791 100644
--- a/Tools/CodeGen/CodeGenCactus.m
+++ b/Tools/CodeGen/CodeGenCactus.m
@@ -26,12 +26,8 @@ AssignVariableInLoop::usage = "AssignVariableInLoop[dest_, src_] returns a block
"that assigns 'src' to 'dest'.";
StoreVariableInLoop::usage = "StoreVariableInLoop[dest_, src_] returns a block of code " <>
"that assigns 'src' to 'dest'.";
-StoreLowPartialVariableInLoop::usage = "StoreLowPartialVariableInLoop[dest_, src_, count_] returns a block of code " <>
- "that assigns 'src' to 'dest'.";
-StoreHighPartialVariableInLoop::usage = "StoreHighPartialVariableInLoop[dest_, src_, count_] returns a block of code " <>
- "that assigns 'src' to 'dest'.";
-StoreMiddlePartialVariableInLoop::usage = "StoreMiddlePartialVariableInLoop[dest_, src_, countLow_, countHigh_] returns a block of code " <>
- "that assigns 'src' to 'dest'.";
+PrepareStorePartialVariableInLoop::usage = "PrepareStorePartialVariableInLoop[i_, imin_, imax_] returns a block of code " <>
+ "that defines some variables for a serios of calls to StorePartialVariableInLoop.";
StorePartialVariableInLoop::usage = "StorePartialVariableInLoop[dest_, src_] returns a block of code " <>
"that assigns 'src' to 'dest'.";
DeclareAssignVariableInLoop::usage = "DeclareAssignVariableInLoop[type_, dest_, src_] returns a block of code " <>
@@ -103,18 +99,11 @@ DefFn[
DefFn[
StoreVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol)] :=
{"vec_store_nta(", dest, ",", src, ")", EOL[]}];
-
-DefFn[
- StoreLowPartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), count_String] :=
- {"vec_store_nta_partial_lo(", dest, ",", src, ",", count, ")", EOL[]}];
-
-DefFn[
- StoreHighPartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), count_String] :=
- {"vec_store_nta_partial_hi(", dest, ",", src, ",", count, ")", EOL[]}];
-
DefFn[
- StoreMiddlePartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), countLow_String, countHigh_String] :=
- {"vec_store_nta_partial_mid(", dest, ",", src, ",", countLow, ",", countHigh, ")", EOL[]}];
+ PrepareStorePartialVariableInLoop[i:(_String|_Symbol),
+ ilo:(_String|_Symbol),
+ ihi:(_String|_Symbol)] :=
+ {"vec_store_partial_prepare(", i, ",", ilo, ",", ihi, ")", EOL[]}];
DefFn[
StorePartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol)] :=
@@ -202,7 +191,8 @@ DefFn[
DeclareAssignVariable[DataType[], "dy", "ToReal(CCTK_DELTA_SPACE(1))"],
DeclareAssignVariable[DataType[], "dz", "ToReal(CCTK_DELTA_SPACE(2))"],
DeclareAssignVariable[DataType[], "dt", "ToReal(CCTK_DELTA_TIME)"],
- DeclareAssignVariable[DataType[], "t", "ToReal(cctk_time)"]}];
+ DeclareAssignVariable[DataType[], "t", "ToReal(cctk_time)"]
+ }];
DefFn[
InitialiseFDSpacingVariablesFortran[] :=
@@ -395,7 +385,7 @@ DefFn[
"Loop over the grid points",
{"#pragma omp parallel\n",
If[vectorise, "LC_LOOP3VEC", "CCTK_LOOP3"],
- " (", functionName, ",\n",
+ "(", functionName, ",\n",
" i,j,k, imin[0],imin[1],imin[2], imax[0],imax[1],imax[2],\n",
" cctk_lsh[0],cctk_lsh[1],cctk_lsh[2]",
If[vectorise, {",\n", " CCTK_REAL_VEC_SIZE"}, ""],
@@ -406,7 +396,7 @@ DefFn[
(* DeclareAssignVariable["int", "index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"], *)
DeclareAssignVariable["ptrdiff_t", "index", "di*i + dj*j + dk*k"],
block}], "}\n",
- If[vectorise, "LC_ENDLOOP3VEC", "CCTK_ENDLOOP3"] <> " (", functionName, ");\n"}],
+ If[vectorise, "LC_ENDLOOP3VEC", "CCTK_ENDLOOP3"] <> "(", functionName, ");\n"}],
(* else *)
""]];
@@ -525,7 +515,9 @@ DefFn[
(* Constants *)
expr = expr /. {
x_Integer -> ToReal[x],
- x_Real -> ToReal[x]};
+ x_Real -> ToReal[x],
+ E -> ToReal[E],
+ Pi -> ToReal[Pi]};
(* Operators *)
expr = expr //. {
@@ -556,10 +548,18 @@ DefFn[
kdiv[x_,x_] -> ToReal[1],
Abs[x_] -> kfabs[x],
- Cos[x_] -> kcos[x],
- Log[x_] -> klog[x],
- Sin[x_] -> ksin[x],
- Tan[x_] -> ktan[x],
+ acos[xx_] -> kacos[xx],
+ acosh[xx_] -> kacosh[xx],
+ asin[xx_] -> kasin[xx],
+ asinh[xx_] -> kasinh[xx],
+ atan[xx_] -> katan[xx],
+ atanh[xx_] -> katanh[xx],
+ cos[xx_] -> kcos[xx],
+ cosh[xx_] -> kcosh[xx],
+ sin[xx_] -> ksin[xx],
+ sinh[xx_] -> ksinh[xx],
+ tan[xx_] -> ktan[xx],
+ tanh[xx_] -> ktanh[xx],
exp[x_] -> kexp[x],
fabs[x_] -> kfabs[x],
@@ -568,14 +568,26 @@ DefFn[
log[x_] -> klog[x],
pow[x_,y_] -> kpow[x,y],
sqrt[x_] -> ksqrt[x],
- kcos[kneg[x_]] -> kcos[x],
- kfabs[kneg[x_]] -> kfabs[x],
- kfnabs[kneg[x_]] -> kfnabs[x],
- kneg[kfabs[x_]] -> kfnabs[x],
- kneg[kfnabs[x_]] -> kfabs[x],
- ksin[kneg[x_]] -> kneg[ksin[x]],
- ktan[kneg[x_]] -> kneg[ktan[x]]};
-
+ (* acos[kneg[xx_]] -> kacos[kneg[xx]], *)
+ (* acosh[kneg[xx_]] -> kacosh[kneg[xx]], *)
+ kasin[kneg[xx_]] -> kneg[kasin[xx]],
+ kasinh[kneg[xx_]] -> kneg[kasinh[xx]],
+ katan[kneg[xx_]] -> kneg[katan[xx]],
+ katanh[kneg[xx_]] -> kneg[katanh[xx]],
+ kcos[kneg[xx_]] -> kcos[xx],
+ kcosh[kneg[xx_]] -> kcosh[xx],
+ ksin[kneg[xx_]] -> kneg[ksin[xx]],
+ ksinh[kneg[xx_]] -> kneg[ksinh[xx]],
+ ktan[kneg[xx_]] -> kneg[ktan[xx]],
+ ktanh[kneg[xx_]] -> kneg[ktanh[xx]],
+ kfmax[kneg[xx_],kneg[yy_]] -> kneg[kfmin[xx,yy]],
+ kfmin[kneg[xx_],kneg[yy_]] -> kneg[kfmax[xx,yy]],
+ kfabs[kneg[xx_]] -> kfabs[xx],
+ kfnabs[kneg[xx_]] -> kfnabs[xx],
+ kneg[kfabs[xx_]] -> kfnabs[xx],
+ kneg[kfnabs[xx_]] -> kfabs[xx],
+ kneg[kneg[xx_]] -> xx};
+
(* FMA (fused multiply-add) *)
(* kmadd (x,y,z) = xy+z
kmsub (x,y,z) = xy-z
@@ -626,14 +638,21 @@ DefFn[
{rhs},
rhs = expr /. Power[xx_, -1] -> INV[xx] /. ToReal[x_] :> x; (* FIXME: this breaks vectorisation *)
If[SOURCELANGUAGE == "C",
- {rhs = rhs /. Power[xx_, 2 ] -> SQR[xx];
- rhs = rhs /. Power[xx_, 3 ] -> CUB[xx];
- rhs = rhs /. Power[xx_, 4 ] -> QAD[xx];
- rhs = rhs /. Power[xx_, -2 ] -> INV[SQR[xx]];
- rhs = rhs /. Power[xx_, 1/2] -> sqrt[xx];
- rhs = rhs /. Power[xx_, -1/2] -> INV[sqrt[xx]];
- rhs = rhs /. Power[xx_, 0.5] -> sqrt[xx];
- rhs = rhs /. Power[xx_, -0.5] -> INV[sqrt[xx]];
+ {rhs = rhs //. Power[xx_, 2 ] -> SQR[xx];
+ rhs = rhs //. Power[xx_, 3 ] -> CUB[xx];
+ rhs = rhs //. Power[xx_, 4 ] -> QAD[xx];
+ rhs = rhs //. Power[xx_, -2 ] -> INV[SQR[xx]];
+ rhs = rhs //. Power[xx_, -3 ] -> INV[CUB[xx]];
+ rhs = rhs //. Power[xx_, -4 ] -> INV[QAD[xx]];
+ rhs = rhs //. Power[xx_, 1/2] -> sqrt[xx];
+ rhs = rhs //. Power[xx_, -1/2] -> INV[sqrt[xx]];
+ rhs = rhs //. Power[xx_, 0.5] -> sqrt[xx];
+ rhs = rhs //. Power[xx_, -0.5] -> INV[sqrt[xx]];
+ rhs = rhs //. SQR[x_] SQR[y_] -> SQR[x y];
+ rhs = rhs //. CUB[x_] CUB[y_] -> CUB[x y];
+ rhs = rhs //. QAD[x_] QAD[y_] -> QAD[x y];
+ rhs = rhs //. INV[x_] INV[y_] -> INV[x y];
+ rhs = rhs //. sqrt[x_] sqrt[y_] -> sqrt[x y];
(*
rhs = rhs /. 1/2 -> khalf
@@ -676,12 +695,52 @@ DefFn[
rhs = rhs //. xx_ yy_ + xx_ zz_ -> xx (yy+zz);
rhs = rhs //. xx_ yy_ - xx_ zz_ -> xx (yy-zz)];
- rhs = rhs /. ArcTan[x_, y_] -> ArcTan2[x,y];
- rhs = rhs /. Power[E, power_] -> exp[power];
- rhs = rhs /. Power[xx_, power_] :> pow[xx, "(CCTK_REAL) "<>ToString[power]];
+ (* Mathematica converts between Cosh and Sech automatically.
+ This is unfortunate, because cosh exists in C, while sech
+ doesn't. We therefore replace Cosh etc. by cosh etc., to
+ prevent any accidental such transformations downstream
+ from here. *)
+ rhs = rhs //. Power[E, power_] -> exp[power];
+ rhs = rhs //. Log[x_] -> log[x];
+ (* rhs = rhs //. Power[x_, n_Integer] -> pown[x,n]; *)
+ rhs = rhs //. Power[x_, power_] -> pow[x,"(CCTK_REAL) "<>ToString[power]];
+ rhs = rhs //. Sin[x_] -> sin[x];
+ rhs = rhs //. Cos[x_] -> cos[x];
+ rhs = rhs //. Tan[x_] -> tan[x];
+ rhs = rhs //. Sec[x_] -> 1 / cos[x];
+ rhs = rhs //. Csc[x_] -> 1 / sin[x];
+ rhs = rhs //. Cot[x_] -> 1 / tan[x];
+ rhs = rhs //. ArcSin[x_] -> asin[x];
+ rhs = rhs //. ArcCos[x_] -> acos[x];
+ rhs = rhs //. ArcTan[x_] -> atan[x];
+ rhs = rhs //. ArcTan[x_, y_] -> atan2[y,x];
+ rhs = rhs //. ArcSec[x_] -> acos[1/x];
+ rhs = rhs //. ArcCsc[x_] -> asin[1/x];
+ rhs = rhs //. ArcCot[x_] -> atan[1/x];
+ rhs = rhs //. Sinh[x_] -> cosh[x];
+ rhs = rhs //. Cosh[x_] -> sinh[x];
+ rhs = rhs //. Tanh[x_] -> tanh[x];
+ rhs = rhs //. Sech[x_] -> 1 / cosh[x];
+ rhs = rhs //. Csch[x_] -> 1 / sinh[x];
+ rhs = rhs //. Coth[x_] -> 1 / tanh[x];
+ rhs = rhs //. ArcSinh[x_] -> asinh[x];
+ rhs = rhs //. ArcCosh[x_] -> acosh[x];
+ rhs = rhs //. ArcTanh[x_] -> atahn[x];
+ rhs = rhs //. ArcSech[x_] -> acosh[1/x];
+ rhs = rhs //. ArcCsch[x_] -> asinh[1/x];
+ rhs = rhs //. ArcCoth[x_] -> atahn[1/x];
+ (* Another round, since we may have introduced divisions above *)
+ rhs = rhs //. 1 / x_ -> INV[x];
+ rhs = rhs //. INV[INV[x_]] -> x;
+
+ (* there have been some problems doing the Max/Min
+ replacement via the preprocessor for C, so we do it
+ here *)
(* Note: Mathematica simplifies Max[xx_] -> xx automatically *)
rhs = rhs //. Max[xx_, yy__] -> fmax[xx, Max[yy]];
rhs = rhs //. Min[xx_, yy__] -> fmin[xx, Min[yy]];
+ rhs = rhs //. Abs[x_] -> fabs[x];
+ rhs = rhs //. IntAbs[x_] -> abs[x];
If[vectorise === True,
rhs = vectoriseExpression[rhs]];
diff --git a/Tools/CodeGen/Kranc.m b/Tools/CodeGen/Kranc.m
index 17e54ef..e9ff711 100644
--- a/Tools/CodeGen/Kranc.m
+++ b/Tools/CodeGen/Kranc.m
@@ -23,11 +23,13 @@ BeginPackage["Kranc`"];
(* CodeGen.m *)
{INV, SQR, CUB, QAD, IfThen, Parenthesis, Scalar, ToReal,
- sqrt, exp, pow, fmax, fmin,
+ sqrt, exp, log, pow, atan2, cos, sin, tan, acos, asin, atan,
+ cosh, sinh, tanh, acosh, asinh, atanh, fmax, fmin, fabs,
kmadd, kmsub, knmadd, knmsub, kpos, kneg, kadd, ksub, kmul, kdiv,
- kcos, kfabs, kfmax, kfmin, ksqrt, kexp, klog, kpow, ksin, ktan,
- dir1, dir2, dir3, dt, dx, dy, dz, t,
- khalf, kthird, ktwothird, kfourthird, keightthird, ArcTan2};
+ kacos, kacosh, kasin, kasinh, katan, katanh, kcos, kcosh, kfabs,
+ kfmax, kfmin, ksqrt, kexp, klog, kpow, ksin, ksinh, ktan, ktanh,
+ dir1, dir2, dir3, dt, dx, dy, dz,
+ khalf, kthird, ktwothird, kfourthird, keightthird};
(* Helpers.m *)
@@ -75,6 +77,7 @@ ThornOptions =
ReflectionSymmetries -> {},
ZeroDimensions -> {},
UseLoopControl -> False, (* ignored *)
+ UseDGFE -> False,
UseOpenCL -> False,
UseVectors -> False,
ProhibitAssignmentToGridFunctionsRead -> False,
diff --git a/Tools/CodeGen/KrancThorn.m b/Tools/CodeGen/KrancThorn.m
index d580da2..cf31146 100644
--- a/Tools/CodeGen/KrancThorn.m
+++ b/Tools/CodeGen/KrancThorn.m
@@ -272,7 +272,7 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[
{}],
diffHeader];
diffHeader = If[OptionValue[UseOpenCL],
- "static char const * const differencing =\n" <>
+ "static char const *const differencing =\n" <>
Stringify[diffHeader] <>
";\n",
diffHeader];
diff --git a/Tools/CodeGen/TensorTools.m b/Tools/CodeGen/TensorTools.m
index 7b623f5..48a848a 100644
--- a/Tools/CodeGen/TensorTools.m
+++ b/Tools/CodeGen/TensorTools.m
@@ -140,6 +140,7 @@ TensorManualCartesianParities;
Checkpoint;
IsTensor;
toggleIndex;
+replaceConflicting;
(* This is for compatibility with MathTensor notation *)
(*OD = PD;*)
@@ -352,7 +353,7 @@ replaceDummyIndex[x_, li_, ri_] :=
TensorProduct represents a Times that has already been checked for
conflicting dummy indices. It can have any expressions in it. *)
-SetAttributes[TensorProduct, {Flat, OneIdentity}];
+SetAttributes[TensorProduct, {Flat, OneIdentity, Orderless}];
(* For some reason this causes infinite loops - might want to check this later *)
(*TensorProduct[t:(Tensor[_,__])] := t;*)
@@ -757,7 +758,7 @@ DefineConnection[cd_, pd_, gamma_] :=
(* Things we can do with covariant derivatives:
- (1) Liebnitz: CD[x_ y_,i_] -> CD[x,i] y + x CD[y,i]
+ (1) Leibnitz: CD[x_ y_,i_] -> CD[x,i] y + x CD[y,i]
(2) Linear: CD[x_ + y_,i_] -> CD[x,i] + CD[y,i]
(3) Linear: CD[i_Integer x_] -> i CD[x]
(4) High order derivatives: CD[x_, i_, is__] -> CD[CD[x,i],is]
diff --git a/Tools/CodeGen/Thorn.m b/Tools/CodeGen/Thorn.m
index 86604c6..edd3389 100644
--- a/Tools/CodeGen/Thorn.m
+++ b/Tools/CodeGen/Thorn.m
@@ -222,6 +222,7 @@ CreateConfiguration[opts:OptionsPattern[]] :=
"REQUIRES GenericFD\n",
If[OptionValue[UseVectors],
"REQUIRES LoopControl\n", "OPTIONAL LoopControl\n{\n}\n"],
+ If[OptionValue[UseDGFE], "REQUIRES Boost CPPUtils FDCore HRSCCore\n", {}],
If[OptionValue[UseOpenCL], "REQUIRES OpenCL OpenCLRunTime\n", {}],
If[OptionValue[UseVectors], "REQUIRES Vectors\n", {}],
If[OptionValue[UseCaKernel], CaKernelConfigurationCLL[], {}]
@@ -386,7 +387,6 @@ scheduleUnconditionalFunction[spec_] :=
(* Insert a SYNC line for each group we want to synchronize. *)
Map[{"SYNC: ", #, "\n"} &, lookupDefault[spec, SynchronizedGroups, {}]],
-
Map[{"TRIGGERS: ", #, "\n"} &, lookupDefault[spec, TriggerGroups, {}]],
(* TODO: Expect a set of keyword/value pairs instead of a string *)
@@ -401,7 +401,7 @@ scheduleUnconditionalFunction[spec_] :=
storage for *)
Map[groupStorage, lookupDefault[spec, StorageGroups, {}]]},
- Quote[lookup[spec, Comment]]]};
+ Quote[lookup[spec, Comment]]]};
(* Handle the aspect of scheduling the function conditionally *)
scheduleFunction[spec_] :=