aboutsummaryrefslogtreecommitdiff
path: root/Tools/CodeGen/CalculationFunction.m
diff options
context:
space:
mode:
Diffstat (limited to 'Tools/CodeGen/CalculationFunction.m')
-rw-r--r--Tools/CodeGen/CalculationFunction.m270
1 files changed, 224 insertions, 46 deletions
diff --git a/Tools/CodeGen/CalculationFunction.m b/Tools/CodeGen/CalculationFunction.m
index bdbf040..a2cda03 100644
--- a/Tools/CodeGen/CalculationFunction.m
+++ b/Tools/CodeGen/CalculationFunction.m
@@ -275,7 +275,7 @@ assignVariableFromExpression[dest_, expr_, declare_, vectorise_, noSimplify:Bool
code = StringReplace[code, "normal2" -> "normal[1]"];
code = StringReplace[code, "normal3" -> "normal[2]"];
code = StringReplace[code, "BesselJ"-> "gsl_sf_bessel_Jn"];
- code = StringReplace[code, ToString@tSym -> "cctk_time"];
+ code = StringReplace[code, ToString@tSym -> "ToReal(cctk_time)"];
code = StringReplace[code, "\"" -> ""];
{code}];
@@ -303,7 +303,7 @@ generateCodeFromExpression[expr_, vectorise_, noSimplify:Boolean : False] :=
code = StringReplace[code, "normal2" -> "normal[1]"];
code = StringReplace[code, "normal3" -> "normal[2]"];
code = StringReplace[code, "BesselJ"-> "gsl_sf_bessel_Jn"];
- code = StringReplace[code, ToString@tSym -> "cctk_time"];
+ code = StringReplace[code, ToString@tSym -> "ToReal(cctk_time)"];
code = StringReplace[code, "\"" -> ""];
{code}];
@@ -382,7 +382,7 @@ DefFn[
shorts, eqs, parameters, parameterRules, odeGroups,
functionName, dsUsed, groups, pddefs, cleancalc, eqLoop, where,
addToStencilWidth, pDefs, haveCondTextuals, condTextuals, calc,
- kernelCall},
+ kernelCall, DGFEDefs, DGFECall},
functionName = ToString@lookup[calcp, Name];
bodyFunctionName = functionName <> "_Body";
@@ -416,9 +416,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];
@@ -486,9 +486,210 @@ 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"];
{
- DefineFunction[bodyFunctionName, "static void", "cGH const * restrict const cctkGH, int const dir, int const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], int const imin[3], int const imax[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]",
+ DGFEDefs,
+ DefineFunction[bodyFunctionName, "static void", "cGH const *restrict const cctkGH, int const dir, int const face, CCTK_REAL const normal[3], CCTK_REAL const tangentA[3], CCTK_REAL const tangentB[3], int const imin[3], int const imax[3], int const n_subblock_gfs, CCTK_REAL * restrict const subblock_gfs[]",
{
"DECLARE_CCTK_ARGUMENTS;\n",
"DECLARE_CCTK_PARAMETERS;\n\n",
@@ -497,7 +698,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"
},
{
}],
@@ -539,16 +740,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"
},
{
}]
@@ -582,6 +783,8 @@ DefFn[
kernelCall,
+ DGFECall,
+
ConditionalOnParameterTextual["verbose > 1",
"CCTK_VInfo(CCTK_THORNSTRING,\"Leaving " <> bodyFunctionName <> "\");\n"]
}]]
@@ -778,38 +981,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]}]];
generated by cgit v1.2.3 (git 2.39.1) at 2024-06-27 12:21:57 +0000