diff options
author | Ian Hinder <ian.hinder@aei.mpg.de> | 2011-10-07 00:58:45 +0200 |
---|---|---|
committer | Ian Hinder <ian.hinder@aei.mpg.de> | 2011-10-07 00:58:45 +0200 |
commit | 4c20e5c17b6c343bd8d8c4285591d2c67afad0bf (patch) | |
tree | 12324736c9b7a8e07102afe04bac6b395e86a669 /Tools/CodeGen | |
parent | ad5d11b6f6ee7d69ce3947c265c9e1e429df739a (diff) |
Split CodeGen.m into CodeGen.m, CodeGenC.m and CodeGenCactus.m
- CodeGenC.m still contains Fortran variants. These might be moved into CodeGenFortran.m.
- Several other modules essentially need to pull in all three packages due to newly broken abstraction barriers. This can be fixed incrementally.
- Indent and format all three files consistently
Diffstat (limited to 'Tools/CodeGen')
-rw-r--r-- | Tools/CodeGen/CalculationBoundaries.m | 2 | ||||
-rw-r--r-- | Tools/CodeGen/CalculationFunction.m | 2 | ||||
-rw-r--r-- | Tools/CodeGen/CodeGen.m | 1062 | ||||
-rw-r--r-- | Tools/CodeGen/CodeGenC.m | 243 | ||||
-rw-r--r-- | Tools/CodeGen/CodeGenCactus.m | 727 | ||||
-rw-r--r-- | Tools/CodeGen/Differencing.m | 4 | ||||
-rw-r--r-- | Tools/CodeGen/Jacobian.m | 3 | ||||
-rw-r--r-- | Tools/CodeGen/KrancThorn.m | 2 | ||||
-rw-r--r-- | Tools/CodeGen/Schedule.m | 6 | ||||
-rw-r--r-- | Tools/CodeGen/Thorn.m | 44 |
10 files changed, 1064 insertions, 1031 deletions
diff --git a/Tools/CodeGen/CalculationBoundaries.m b/Tools/CodeGen/CalculationBoundaries.m index 3d331ab..3a66aba 100644 --- a/Tools/CodeGen/CalculationBoundaries.m +++ b/Tools/CodeGen/CalculationBoundaries.m @@ -19,7 +19,7 @@ *) BeginPackage["CalculationBoundaries`", {"Errors`", "Helpers`", - "Kranc`", "MapLookup`", "KrancGroups`", "CodeGen`"}]; + "Kranc`", "MapLookup`", "KrancGroups`", "CodeGen`", "CodeGenCactus`", "CodeGenC`"}]; CalculationBoundariesFunction; diff --git a/Tools/CodeGen/CalculationFunction.m b/Tools/CodeGen/CalculationFunction.m index d1c02bb..7bf9fec 100644 --- a/Tools/CodeGen/CalculationFunction.m +++ b/Tools/CodeGen/CalculationFunction.m @@ -19,7 +19,7 @@ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA *) -BeginPackage["CalculationFunction`", {"CodeGen`", +BeginPackage["CalculationFunction`", {"CodeGenCactus`", "CodeGenC`", "CodeGen`", "MapLookup`", "KrancGroups`", "Differencing`", "Errors`", "Helpers`", "Kranc`", "Optimize`", "Jacobian`"}]; diff --git a/Tools/CodeGen/CodeGen.m b/Tools/CodeGen/CodeGen.m index 09e9863..d478d40 100644 --- a/Tools/CodeGen/CodeGen.m +++ b/Tools/CodeGen/CodeGen.m @@ -22,13 +22,6 @@ BeginPackage["CodeGen`", {"Errors`", "Kranc`"}]; -SOURCELANGUAGE::usage = "global variable == \"C\" or \"Fortran\" determines language - for code generation"; -SOURCESUFFIX::usage = "file suffix for source files"; - -EOL::usage = "the end of line termination string"; -SetSourceLanguage::usage = "set the source language to \"C\" or \"Fortran\""; - FlattenBlock::usage = "FlattenBlock[block] converts 'block' to a string."; SeparatedBlock::usage = "SeparatedBlock[block] returns a version of 'block' with " <> "a newline before it."; @@ -36,1059 +29,128 @@ GenerateFile::usage = "GenerateFile[name, block] writes 'block' to a file of the "specified 'name'."; AddToFile::usage = "AddToFile[name, block] appends 'block' to a file of the " <> "specified 'name'."; -IncludeFile::usage = "IncludeFile[name] returns a block of code" <> - "that includes a header file (i.e '#include \"name\"')."; -IncludeSystemFile::usage = "IncludeFile[name] returns a block of code" <> - "that includes a system header file (i.e '#include <name>')."; -DeclareVariable::usage = "DeclareVariable[name, type] returns a block of code " <> - "that declares a variable of given name and type. 'name' and 'type' should be " <> - "strings."; -DeclareVariableNoInit::usage = "DeclareVariableNoInit[name, type] returns a block of code " <> - "that declares a variable of given name and type without initialising it. 'name' and 'type' should be " <> - "strings."; -DeclareVariables::usage = "DeclareVariables[names, type] returns a block of code " <> - "that declares a list of variables of given name and type. 'names' should be a list" <> - " of strings and 'type' should be a string string."; -DeclarePointer::usage = "DeclarePointer[name, type] returns a block of code " <> - "that declares a pointer of given name and type. 'name' and 'type' should be " <> - "strings."; -DeclarePointers::usage = "DeclarePointers[names, type] returns a block of code " <> - "that declares a list of pointers of given name and type. 'names' should be a list" <> - " of strings and 'type' should be a string string."; -DefineVariable::usage = "DefineVariable[name, type, value] returns a block of " <> - "code that declares and initialised a variable 'name' of type 'type' to value 'value'."; -AssignVariable::usage = "AssignVariable[dest_, src_] returns a block of code " <> - "that assigns 'src' to 'dest'."; -DeclareAssignVariable::usage = "DeclareAssignVariable[type_, dest_, src_] returns a block of code " <> - "that declares and sets a constant variable of given name and type."; -AssignVariableInLoop::usage = "AssignVariableInLoop[dest_, src_] returns a block of code " <> - "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'."; -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 " <> - "that assigns 'src' to 'dest'."; -MaybeAssignVariableInLoop::usage = "MaybeAssignVariableInLoop[dest_, src_, cond_] returns a block of code " <> - "that assigns 'src' to 'dest'."; -DeclareMaybeAssignVariableInLoop::usage = "DeclareMaybeAssignVariableInLoop[type_, dest_, src_, cond_] returns a block of code " <> - "that assigns 'src' to 'dest'."; -TestForNaN::usage = "TestForNaN[expr_] returns a block of code " <> - "that tests 'expr' for nan."; -CommentedBlock::usage = "CommentedBlock[comment, block] returns a block consisting " <> - "of 'comment' followed by 'block'."; -DefineCCTKFunction::usage = "DefineCCTKFunction[name, type, block] returns a block " <> - "of code that defines a CCTK function of name 'name' returning type 'type' with " <> - "body 'block'."; -DefineCCTKSubroutine::usage = "DefineCCTKSubroutine[name, block] returns a block " <> - "of code that defines a CCTK Fortran subroutine of name 'name' with body 'block'."; -GridName::usage = "GridName[variable] returns the name needed to access variable " <> - "assuming it is a grid variable when inside a grid loop."; -DeclareGridLoopVariables::usage = "DeclareGridLoopVariables[] returns a block " <> - "that defines the variables needed during a grid loop."; -InitialiseGridLoopVariables::usage = "InitialiseGridLoopVariables[] returns a block " <> - "that initialises variables needed by a grid loop."; -InitialiseGridLoopVariablesWithStencil::usage = "InitialiseGridLoopVariables[] returns a block " <> - "that initialises variables needed by a grid loop using the evolution stencils."; -ConditionalOnParameter::usage = "ConditionalOnParameter[name, value, block] returns " <> - "a block that introduces a conditional expression whereby 'block' is only executed " <> - "if the Cactus parameter 'name' has value 'value'."; -(* -GridLoop::usage = "GridLoop[block] returns a block that is looped over for every " <> - "grid point. Must have previously set up the grid loop variables (see " <> - "DeclareGridLoopVariables and InitialiseGridLoopVariables."; -*) -GridLoop::usage = "GridLoop[block] returns a block that is looped over for every " <> - "grid point. Must have previously set up the grid loop variables (see " <> - "InitialiseGridLoopVariables."; - -DeclareArray::usage = ""; - -DefineFunction::usage = ""; -ConditionalOnParameterTextual::usage = ""; - -SpaceSeparated::usage = ""; (* This should not really be in CodeGen *) -NewlineSeparated::usage = ""; (* This should not really be in CodeGen *) - -CBlock::usage = ""; -SuffixedCBlock::usage = ""; +SpaceSeparated::usage = ""; +NewlineSeparated::usage = ""; InfoVariable::usage = ""; - -DeclareFDVariables::usage = ""; - -InitialiseFDVariables::usage = ""; - -CommaNewlineSeparated::usage = ""; (* This should not really be in CodeGen *) +CommaNewlineSeparated::usage = ""; CommaSeparated::usage = ""; -ReplacePowers::usage = ""; -CFormHideStrings::usage = ""; -BoundaryLoop::usage = ""; -BoundaryWithGhostsLoop::usage = ""; -GenericGridLoop::usage = ""; - Stringify::usage = ""; - -NameRoot::usage = ""; -PartitionVarList::usage = ""; Quote::usage = "Quote[x] returns x surrounded by quotes"; -DataType::usage = "DataType[] returns a string for the grid function data type (e.g. CCTK_REAL)"; -SetDataType::usage = "SetDataType[type] sets a string for the grid function data type (e.g. CCTK_REAL)"; -Conditional; -SwitchStatement; +IndentBlock::usage = ""; +CheckBlock::usage = ""; CodeGenBlock := _String | _?AtomQ | List[(_?(MatchQ[#, CodeGenBlock] &)) ...]; - Boolean = (True | False); Begin["`Private`"]; -SOURCELANGUAGE = "C"; -SOURCESUFFIX = ".cc"; - -setSourceSuffix[lang_] := -If[ (lang == "C"), - SOURCESUFFIX = ".cc"; -, - SOURCESUFFIX = ".F90"; -]; - - -SetSourceLanguage[lang_] := -If[ (lang == "C" || lang == "Fortran"), - SOURCELANGUAGE = lang; - setSourceSuffix[lang]; - InfoMessage[Terse, "User set source language to " <> lang], - - SOURCELANGUAGE = "C"; - setSourceSuffix[".cc"]; - InfoMessage[Terse, "Setting Source Language to C"]; -]; - -SetDataType[type_String] := - dataType = type; - -DataType[] := - If[dataType === Symbol["datatype"], - Throw["DataType: Have not set a data type"], - dataType]; - (* Code generation utilities; not specific to any language *) -FlattenBlock[x_String] := x; -FlattenBlock[l_List] := StringJoin@@Map[FlattenBlock, l]; -FlattenBlock[a_?AtomQ] := ToString[a]; +CheckBlock[s_String] := s; +CheckBlock[a_?AtomQ] := a; +CheckBlock[l_List] := Map[CheckBlock, l]; +ErrorDefinition[CheckBlock]; + +FlattenBlock[b_] := + Module[ + {flattenBlock}, + flattenBlock[x_String] := x; + flattenBlock[l_List] := StringJoin@@Map[FlattenBlock, l]; + flattenBlock[a_?AtomQ] := ToString[a]; + + CheckBlock[b]; + flattenBlock[b]]; + ErrorDefinition[FlattenBlock]; -indentBlock[block:CodeGenBlock] := +IndentBlock[block:CodeGenBlock] := StringDrop[" " <> StringReplace[FlattenBlock[block], {"\n" -> "\n "}],-2]; -ErrorDefinition[indentBlock]; +ErrorDefinition[IndentBlock]; SeparatedBlock[block:CodeGenBlock] := {"\n", block}; ErrorDefinition[SeparatedBlock]; -GenerateFile[filename_String, contents:CodeGenBlock] := - Module[{fp = OpenWrite[filename]}, +GenerateFile[filename_String, contents_] := + Module[ + {fp = OpenWrite[filename]}, + CheckBlock[contents]; WriteString[fp, FlattenBlock[contents]]; Close[fp]]; ErrorDefinition[GenerateFile]; AddToFile[filename_String, contents:CodeGenBlock] := - Module[{fp = OpenAppend[filename]}, - WriteString[fp, FlattenBlock[contents]]; - Close[fp]]; + Module[ + {fp = OpenAppend[filename]}, + WriteString[fp, FlattenBlock[contents]]; + Close[fp]]; ErrorDefinition[AddToFile]; -CommaNewlineSeparated[l_List] := Riffle[l, ",\n"]; +CommaNewlineSeparated[l_List] := + Riffle[l, ",\n"]; ErrorDefinition[CommaNewlineSeparated]; -SpaceSeparated[l_List] := +SpaceSeparated[l_List] := Riffle[l, " "]; ErrorDefinition[SpaceSeparated]; -CommaSeparated[l_List] := +CommaSeparated[l_List] := Riffle[l, ", "]; ErrorDefinition[CommaSeparated]; -NewlineSeparated[l_List] := +NewlineSeparated[l_List] := Riffle[l, "\n"]; ErrorDefinition[NewlineSeparated]; CommaInitSeparated[l_List] := Riffle[Map[{#," = INITVALUE"} &, l], ", "]; -(* Riffle[l, " = INITVALUE, "];*) ErrorDefinition[CommaInitSeparated]; - - (* Turn a section of code into a string: 1. quote all quotes (replace all quotes with backslash-quote) 2. break the string into lines to make it readable (replace all newlines with quote-newline-quote) 3. surround the result with quotes *) -Stringify[x:CodeGenBlock] := "\"" <> - StringReplace[StringReplace[FlattenBlock[x], - "\"" -> "\\\""], "\n" -> "\\n\"\n\""] <> - "\"\n"; +Stringify[x:CodeGenBlock] := + "\"" <> StringReplace[StringReplace[FlattenBlock[x], "\"" -> "\\\""], + "\n" -> "\\n\"\n\""] <> "\"\n"; ErrorDefinition[Stringify]; +PartitionVarList[list_List] := + Module[ + {partition, split}, + partition[locallist_] := + Module[ + {cutoff}, + cutoff = 6; + If[Length@locallist > cutoff, Partition[locallist, cutoff, cutoff, {1,1}, {}], + {locallist}]]; -NameRoot[name_Symbol] := Module[{dropNumberRule, root}, - - dropNumberRule = {"1" -> "", "2" -> "", "3" -> "", "4" -> "", "5" -> "", - "6" -> "", "7" -> "", "8" -> "", "9" -> "", "0" -> "", "rhs" -> ""}; - - root = StringReplace[ToString@name, dropNumberRule] - ]; -ErrorDefinition[NameRoot]; - -PartitionVarList[list_List]:= Module[{partition, split}, - -partition[locallist_] := Module[{cutoff}, - cutoff = 6; - If[Length@locallist > cutoff, Partition[locallist, cutoff, cutoff, {1,1}, {}], {locallist}] -]; - -split = Split[list, NameRoot[#1] == NameRoot[#2] &]; -split = Flatten[Map[partition, split], 1]; + split = Split[list, NameRoot[#1] == NameRoot[#2] &]; + split = Flatten[Map[partition, split], 1]; -split -]; + split]; ErrorDefinition[PartitionVarList]; - -(* Code generation for generic C and C-preprocessed Fortran *) - -EOL[dummy___] := If[SOURCELANGUAGE == "C" || SOURCELANGUAGE == "C++", ";\n", "\n"]; - -IncludeFile[filename_String] := - {"#include \"", filename, "\"\n"}; -ErrorDefinition[IncludeFile]; - -IncludeSystemFile[filename_String] := - {"#include <", filename, ">\n"}; -ErrorDefinition[IncludeSystemFile]; - -DeclareVariable[name:(_String|_Symbol), type_String] := -If[SOURCELANGUAGE == "C", - {type, " ", name, " = INITVALUE" <> EOL[]}, - {type, " :: ", name, EOL[]} (* no value init here to avoid implicit SAVE attribute *) - ]; -ErrorDefinition[DeclareVariable]; - -DeclareVariableNoInit[name:(_String|_Symbol), type_String] := -If[SOURCELANGUAGE == "C", - {type, " ", name, EOL[]}, - {type, " :: ", name, EOL[]} (* no value init here to avoid implicit SAVE attribute *) - ]; -ErrorDefinition[DeclareVariableNoInit]; - -DeclareVariables[names_?ListQ, type_String] := -If[SOURCELANGUAGE == "C", - {type, " ", CommaSeparated@names, EOL[]}, - {type, " :: ", CommaSeparated@names, EOL[]} (* no value init avoids implicit SAVE attribute *) - ]; -ErrorDefinition[DeclareVariables]; - -DeclarePointer[name:(_String|_Symbol), type_String] := -If[SOURCELANGUAGE == "C", - {type, " *", name, EOL[]}, - {type, ", target :: ", name, EOL[]} - ]; -ErrorDefinition[DeclarePointer]; - -DeclarePointers[names_?ListQ, type_String] := -If[SOURCELANGUAGE == "C", - {type, " *", CommaInitSeparated@names, EOL[]}, - {type, ", target :: ", CommaSeparated@names, EOL[]} - ]; -ErrorDefinition[DeclarePointers]; - -DeclareArray[name:(_String|_Symbol), dim_Integer, type_String] := - If[SOURCELANGUAGE == "C", - DeclareArrayC[name, dim, type], - DeclareArrayFortran[name, dim, type]]; -ErrorDefinition[DeclareArray]; - -DeclareArrayC[name:(_String|_Symbol), dim_Integer, type_String] := - {type, " ", name, "[", dim, "];","\n"}; -ErrorDefinition[DeclareArrayC]; - -DeclareArrayFortran[name:(_String|_Symbol), dim_Integer, type_String] := - {type, " :: ", name, "(", dim, ")","\n"}; -ErrorDefinition[DeclareArrayFortran]; - -DefineVariable[name:(_String|_Symbol), type_String, value:CodeGenBlock] := - {type, " ", name, " = ", value, EOL[]}; -ErrorDefinition[DefineVariable]; - -AssignVariable[dest:(_String|_Symbol), src:CodeGenBlock] := - {dest, " = ", src, EOL[]}; -ErrorDefinition[AssignVariable]; - -DeclareAssignVariable[type_String, dest:(_String|_Symbol), src:CodeGenBlock] := - {type, " const ", dest, " = ", src, EOL[]}; -ErrorDefinition[DeclareAssignVariable]; - -AssignVariableInLoop[dest:(_String|_Symbol), src:CodeGenBlock, vectorise_:False] := - Module[{loader}, - loader[x_] := If[vectorise, {"vec_load(", x, ")"}, x]; - {dest, " = ", loader[src], EOL[]}]; -(* - {dest, " = ", src, EOL[], - TestForNaN[dest]}; -*) -ErrorDefinition[AssignVariableInLoop]; - -StoreVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol)] := - {"vec_store_nta(", dest, ",", src, ")", EOL[]}; -ErrorDefinition[StoreVariableInLoop]; - -StoreLowPartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), count_String] := - {"vec_store_nta_partial_lo(", dest, ",", src, ",", count, ")", EOL[]}; -ErrorDefinition[StoreLowPartialVariableInLoop]; - -StoreHighPartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), count_String] := - {"vec_store_nta_partial_hi(", dest, ",", src, ",", count, ")", EOL[]}; -ErrorDefinition[StoreHighPartialVariableInLoop]; - -StoreMiddlePartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), countLow_String, countHigh_String] := - {"vec_store_nta_partial_mid(", dest, ",", src, ",", countLow, ",", countHigh, ")", EOL[]}; -ErrorDefinition[StoreMiddlePartialVariableInLoop]; - -StorePartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol)] := - {"vec_store_nta_partial(", dest, ",", src, ")", EOL[]}; -ErrorDefinition[StorePartialVariableInLoop]; - -DeclareAssignVariableInLoop[type_String, dest:(_String|_Symbol), src:(_String|_Symbol)] := - {type, " const ", dest, " = vec_load(", src, ")", EOL[]}; -ErrorDefinition[DeclareAssignVariableInLoop]; - -MaybeAssignVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), cond:Boolean] := - If [cond, - {dest, " = useMatter ? vec_load(", src, ") : ToReal(0.0)", EOL[]}, - {dest, " = vec_load(", src, ")", EOL[]}]; -ErrorDefinition[MaybeAssignVariableInLoop]; - -DeclareMaybeAssignVariableInLoop[type_String, dest:(_String|_Symbol), src:(_String|_Symbol), mmaCond:Boolean, codeCond:CodeGenBlock, - vectorise:Boolean:False] := - Module[{loader}, - loader[x_] := If[vectorise, {"vec_load(", x, ")"}, x]; - If [mmaCond, - {type, " ", dest, " = (", codeCond, ") ? ", loader[src], " : ToReal(0.0)", EOL[]}, - {type, " ", dest, " = ", loader[src], EOL[]}]]; -ErrorDefinition[DeclareMaybeAssignVariableInLoop]; - -TestForNaN[expr:CodeGenBlock] := - {"if (isnan(", expr, ")) {\n", - " CCTK_VInfo(CCTK_THORNSTRING, \"NaN found\");\n", - " CCTK_VInfo(CCTK_THORNSTRING, \"ipos: %d %d %d\", i, j, k);\n", - " CCTK_VInfo(CCTK_THORNSTRING, \"lbnd: %d %d %d\", cctk_lbnd[0], cctk_lbnd[1], cctk_lbnd[2]);\n", - " CCTK_VInfo(CCTK_THORNSTRING, \"lsh: %d %d %d\", cctk_lsh[0], cctk_lsh[1], cctk_lsh[2]);\n", - " CCTK_VInfo(CCTK_THORNSTRING, \"LSSH: %d %d %d\", CCTK_LSSH(0,0), CCTK_LSSH(0,1), CCTK_LSSH(0,2));\n", - " CCTK_VInfo(CCTK_THORNSTRING, \"", expr, ": %.17g\", (double)", expr, ");\n", - "}\n"}; -ErrorDefinition[TestForNaN]; - -(* comments are always done C-style because they are killed by cpp anyway *) -insertComment[text:CodeGenBlock] := {"/* ", text, " */\n"}; -ErrorDefinition[insertComment]; - -CBlock[block:CodeGenBlock] := - {"{\n", - indentBlock[block], - "}\n"}; -ErrorDefinition[CBlock]; - -SuffixedCBlock[block:CodeGenBlock, suffix_] := - {"{\n", - indentBlock[block], - "} ", suffix, "\n"}; -ErrorDefinition[SuffixedCBlock]; - - -loopOverInteger[name_String, start_String, endplusone_String, block_CodeGenBlock] := -If[SOURCELANGUAGE == "C" || SOURCELANGUAGE == "C++", - - {"for (", name, " = ", start, "; ", name, " < ", endplusone, "; ", name, "++)\n", - "{\n", - indentBlock[block], - "}\n"}, - - {"Do ", name, " = ", start, ", ", endplusone, "\n", - "\n", - indentBlock[block], - "End Do\n"} -]; -ErrorDefinition[loopOverInteger]; - - -CommentedBlock[comment:CodeGenBlock, block:CodeGenBlock] := - SeparatedBlock[{insertComment[comment], - block}]; -ErrorDefinition[CommentedBlock]; - -(* FUNCTIONS *) - -defineFunctionC[name_String, type_String, args:CodeGenBlock, contents:CodeGenBlock] := - SeparatedBlock[ - {type, " ", name, "(", args, ")\n", - CBlock[contents]}]; -ErrorDefinition[defineFunctionC]; - -defineFunctionF[name_String, args:CodeGenBlock, contents:CodeGenBlock] := - SeparatedBlock[ - {"FUNCTION", " ", name, "(", args, ")\n", - indentBlock[contents]}]; -ErrorDefinition[defineFunctionF]; - -DefineFunction[name_String, type_String, args:CodeGenBlock, contents:CodeGenBlock] := - If[SOURCELANGUAGE == "C", - defineFunctionC[name, type, args, contents], - defineFunctionF[name, args, contents]]; -ErrorDefinition[DefineFunction]; - -(* SUBROUTINES *) - -defineSubroutine[name_String, args:CodeGenBlock, contents:CodeGenBlock] := - If[SOURCELANGUAGE == "C", - defineSubroutineC[name, args, contents], - defineSubroutineF[name, args, contents]]; -ErrorDefinition[defineSubroutine]; - -defineSubroutineC[name_String, args:CodeGenBlock, contents:CodeGenBlock] := - SeparatedBlock[ - {"extern \"C\" void ", name, "(", args, ")", "\n", - CBlock[contents]}]; -ErrorDefinition[defineSubroutineC]; - -defineSubroutineF[name_String, args:CodeGenBlock, contents:CodeGenBlock] := - SeparatedBlock[ - {"subroutine ", name, "(", args, ")", "\n", - "\nimplicit none\n\n", - contents, - "end subroutine\n"}]; -ErrorDefinition[defineSubroutineF]; - - - - -(********* Code generation for Cactus C or Fortran code **********) - - - - -(* This is a Cactus-callable function *) -DefineCCTKFunction[name_String, type_String, contents:CodeGenBlock] := - DefineFunction[name, "extern \"C\" " <> type, "CCTK_ARGUMENTS", - { - "DECLARE_CCTK_ARGUMENTS;\n", - "DECLARE_CCTK_PARAMETERS;\n\n", - contents - }]; -ErrorDefinition[DefineCCTKFunction]; - -(* This is a Cactus-callable subroutine *) -DefineCCTKSubroutine[name_String, contents:CodeGenBlock] := - defineSubroutine[ - name, "CCTK_ARGUMENTS", - { - "DECLARE_CCTK_ARGUMENTS;\n", - "DECLARE_CCTK_PARAMETERS;\n\n", - contents - }]; -ErrorDefinition[DefineCCTKSubroutine]; - -DeclareFDVariables[] := -(* - CommentedBlock["Declare finite differencing variables", - {Map[DeclareVariables[#, "CCTK_REAL"] &, {{"dx", "dy", "dz"}, - {"dxi", "dyi", "dzi"}, - {khalf,kthird,ktwothird,kfourthird,keightthird}, - {"hdxi", "hdyi", "hdzi"}}], - "\n"}, - {Map[DeclareVariables[#, "ptrdiff_t"] &, {{"di", "dj", "dk"}}], - "\n"}]; -*) - CommentedBlock["Declare finite differencing variables", {}]; - -InitialiseFDSpacingVariablesC[] := - { - (* DeclareAssignVariable["ptrdiff_t", "di", "CCTK_GFINDEX3D(cctkGH,1,0,0) - CCTK_GFINDEX3D(cctkGH,0,0,0)"], *) - DeclareAssignVariable["ptrdiff_t", "di", "1"], - DeclareAssignVariable["ptrdiff_t", "dj", "CCTK_GFINDEX3D(cctkGH,0,1,0) - CCTK_GFINDEX3D(cctkGH,0,0,0)"], - DeclareAssignVariable["ptrdiff_t", "dk", "CCTK_GFINDEX3D(cctkGH,0,0,1) - CCTK_GFINDEX3D(cctkGH,0,0,0)"], - DeclareAssignVariable["ptrdiff_t", "cdi", "sizeof(CCTK_REAL) * di"], - DeclareAssignVariable["ptrdiff_t", "cdj", "sizeof(CCTK_REAL) * dj"], - DeclareAssignVariable["ptrdiff_t", "cdk", "sizeof(CCTK_REAL) * dk"], - DeclareAssignVariable[DataType[], "dx", "ToReal(CCTK_DELTA_SPACE(0))"], - 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)"] - }; - -InitialiseFDSpacingVariablesFortran[] := - { - AssignVariable["dt", "CCTK_DELTA_TIME"], - AssignVariable["dx", "CCTK_DELTA_SPACE(1)"], - AssignVariable["dy", "CCTK_DELTA_SPACE(2)"], - AssignVariable["dz", "CCTK_DELTA_SPACE(3)"] - } - -InitialiseFDVariables[vectorise:Boolean] := - CommentedBlock["Initialise finite differencing variables", - { If[SOURCELANGUAGE == "Fortran", - InitialiseFDSpacingVariablesFortran[], - InitialiseFDSpacingVariablesC[]], - - DeclareAssignVariable[DataType[], "dxi", "INV(dx)"], - DeclareAssignVariable[DataType[], "dyi", "INV(dy)"], - DeclareAssignVariable[DataType[], "dzi", "INV(dz)"], - If[vectorise, - {DeclareAssignVariable[DataType[], "khalf", "ToReal(0.5)"], - DeclareAssignVariable[DataType[], "kthird", "ToReal(1.0/3.0)"], - DeclareAssignVariable[DataType[], "ktwothird", "ToReal(2.0/3.0)"], - DeclareAssignVariable[DataType[], "kfourthird", "ToReal(4.0/3.0)"], - DeclareAssignVariable[DataType[], "keightthird", "ToReal(8.0/3.0)"], - DeclareAssignVariable[DataType[], "hdxi", "kmul(ToReal(0.5), dxi)"], - DeclareAssignVariable[DataType[], "hdyi", "kmul(ToReal(0.5), dyi)"], - DeclareAssignVariable[DataType[], "hdzi", "kmul(ToReal(0.5), dzi)"]}, - {DeclareAssignVariable[DataType[], "khalf", "0.5"], - DeclareAssignVariable[DataType[], "kthird", "1/3.0"], - DeclareAssignVariable[DataType[], "ktwothird", "2.0/3.0"], - DeclareAssignVariable[DataType[], "kfourthird", "4.0/3.0"], - DeclareAssignVariable[DataType[], "keightthird", "8.0/3.0"], - DeclareAssignVariable[DataType[], "hdxi", "0.5 * dxi"], - DeclareAssignVariable[DataType[], "hdyi", "0.5 * dyi"], - DeclareAssignVariable[DataType[], "hdzi", "0.5 * dzi"]}]}]; -ErrorDefinition[InitialiseFDVariables]; - -GridName[x_Symbol] := If[SOURCELANGUAGE == "C", - ToString[x] <> "[index]", - ToString[x] <> "(i,j,k)" - ]; -ErrorDefinition[GridName]; - -DeclareGridLoopVariables[] := - SeparatedBlock[ - {insertComment["Declare the variables used for looping over grid points"], - Map[DeclareVariables[#, "CCTK_INT"] &, - {{"i", "j", "k"}(*, {"istart", "jstart", "kstart"}, - {"iend", "jend", "kend"}, - {"index_offset_x", "index_offset_y", "index_offset_z", "dir", "face"} *)}] (*, - Map[DeclareArray[#, 6, "CCTK_INT"] &, {"is_symbnd", "is_physbnd", "is_ipbnd"}], - Map[DeclareArray[#, 3, "CCTK_INT"] &, {"imin", "imax", "bmin", "bmax"}] *), - - If[SOURCELANGUAGE == "C", DeclareVariable["index", "// CCTK_INT"], "\n"] - }]; - -(* Access an element of an array; syntax is different between C and - Fortran. Always give this function a C-style array index. *) -arrayElement[var_String, i_String] := - If[SOURCELANGUAGE == "C", - {var, "[", arrayIndex[i], "]"}, - {var, "(", arrayIndex[i], ")"}]; -ErrorDefinition[arrayElement]; - -(* Given a C-style variable index, return the corresponding index for - the language currently in use. The idea is that the caller does not - need to know what language is being used. *) -arrayIndex[i:(_Integer|_String|_Symbol)] := - If[SOURCELANGUAGE == "C", - i, - If[NumberQ[i], i+1, {i, " + 1"}]]; -ErrorDefinition[arrayIndex]; - -max[]:= If[SOURCELANGUAGE == "C", "IMAX", "max"]; - -InitialiseGridLoopVariables[derivativesUsedSwitch:Boolean, addToStencilWidth_Integer] := - CommentedBlock["Set up variables used in the grid loop for the physical grid points", - - If[ (derivativesUsedSwitch), - { - AssignVariable["index_offset_x", max[] <>"(stencil_width, stencil_width_x) + " <> ToString[addToStencilWidth]], - AssignVariable["index_offset_y", max[] <>"(stencil_width, stencil_width_y) + " <> ToString[addToStencilWidth]], - AssignVariable["index_offset_z", max[] <>"(stencil_width, stencil_width_z) + " <> ToString[addToStencilWidth]], - - "\n", - AssignVariable["istart", arrayIndex["index_offset_x"]], - AssignVariable["jstart", arrayIndex["index_offset_y"]], - AssignVariable["kstart", arrayIndex["index_offset_z"]], - - "\n", - AssignVariable["iend", {"CCTK_LSSH(0,0) - index_offset_x"}], - AssignVariable["jend", {"CCTK_LSSH(0,1) - index_offset_y"}], - AssignVariable["kend", {"CCTK_LSSH(0,2) - index_offset_z"}] - }, - - { - AssignVariable["istart", arrayIndex[0]], - AssignVariable["jstart", arrayIndex[0]], - AssignVariable["kstart", arrayIndex[0]], - - "\n", - AssignVariable["iend", "CCTK_LSSH(0,0)"], - AssignVariable["jend", "CCTK_LSSH(0,1)"], - AssignVariable["kend", "CCTK_LSSH(0,2)"] - }] -]; -ErrorDefinition[InitialiseGridLoopVariables]; - - -ConditionalOnParameter[name_String, value_String, block:CodeGenBlock] := - SeparatedBlock[ - {"if (CCTK_EQUALS(", name, ", \"", value, "\"))\n", - "{\n", - indentBlock[block], - "}\n"}]; -ErrorDefinition[ConditionalOnParameter]; - -ConditionalOnParameterTextual[text:CodeGenBlock, block:CodeGenBlock] := - SeparatedBlock[ - {"if (", text, ")\n", - "{\n", - indentBlock[block], - "}\n"}]; -ErrorDefinition[ConditionalOnParameterTextual]; - -(* -GridLoop[block_] := - CommentedBlock["Loop over the grid points", - loopOverInteger["k", "kstart", "kend", - loopOverInteger["j", "jstart", "jend", - loopOverInteger["i", "istart", "iend", - - { If[SOURCELANGUAGE == "C", - AssignVariable["index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"], - ""], - block - } - ]]]]; - -*) - -(* -GridLoop[block_] := - If[SOURCELANGUAGE == "C", - CommentedBlock["Loop over the grid points", - { - "#pragma omp parallel\n", - "LC_LOOP3 (unnamed,\n", - " i,j,k, istart,jstart,kstart, iend,jend,kend,\n", - " cctk_lsh[0],cctk_lsh[1],cctk_lsh[2])\n", - "{\n", - indentBlock[ - { - DeclareVariable["index", "int"], - AssignVariable["index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"], - block - } - ], - "}\n", - "LC_ENDLOOP3 (unnamed);\n" - } - ], - CommentedBlock["Loop over the grid points", - { - "#pragma omp parallel\n", - "LC_LOOP3 (unnamed,\n", - " i,j,k, istart,jstart,kstart, iend,jend,kend,\n", - " cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))\n", - indentBlock[block], - "LC_ENDLOOP3 (unnamed)\n" - } - ] - ]; -*) - -Options[GenericGridLoop] = ThornOptions; - -GenericGridLoop[functionName_String, block:CodeGenBlock, opts:OptionsPattern[]] := - If[OptionValue[UseLoopControl], - GenericGridLoopUsingLoopControl[functionName, block, OptionValue[UseVectors]], - GenericGridLoopTraditional[block]]; -ErrorDefinition[GenericGridLoop]; - -GenericGridLoopTraditional[block:CodeGenBlock] := - CommentedBlock["Loop over the grid points", - loopOverInteger["k", "imin[2]", "imax[2]", - loopOverInteger["j", "imin[1]", "imax[1]", - loopOverInteger["i", "imin[0]", "imax[0]", - - { If[SOURCELANGUAGE == "C", - { - DeclareAssignVariable["int", "index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"] - }, - ""], - block - } - ]]]]; -ErrorDefinition[GenericGridLoopTraditional]; - -GenericGridLoopUsingLoopControl[functionName_String, block:CodeGenBlock, vectorise:Boolean] := - If[SOURCELANGUAGE == "C", - CommentedBlock["Loop over the grid points", - { - "#pragma omp parallel\n", - If[vectorise, "LC_LOOP3VEC", "LC_LOOP3"] <> " (", 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"},""] <> ")\n", - "{\n", - indentBlock[ - { - (* DeclareVariable["index", "// int"], *) - (* 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", "LC_ENDLOOP3"] <> " (", functionName, ");\n" - } - ], - "" - ]; -ErrorDefinition[GenericGridLoopUsingLoopControl]; - -switchOptions[{value:(_String|_Symbol|_?NumberQ), block:CodeGenBlock}] := -{ - "case ", value, ":\n", indentBlock[{block,"break;\n"}] -} -ErrorDefinition[switchOptions]; - -SwitchStatement[var:(_String|_Symbol), pairs__] := -{ - "switch(", var, ")\n", - CBlock[{Riffle[Map[switchOptions, {pairs}],"\n"]}] -} -ErrorDefinition[SwitchStatement]; - - - -BoundaryLoop[block:CodeGenBlock] := -{ - "\nGenericFD_GetBoundaryInfo(cctkGH, cctk_lsh, cctk_lssh, cctk_bbox, cctk_nghostzones, imin, imax, is_symbnd, is_physbnd, is_ipbnd);\n", - - CommentedBlock["Start by looping over the whole grid, minus the NON-PHYSICAL boundary points, which are set by synchronization. ", { - AssignVariable[arrayElement["bmin", 0], "is_physbnd[0*2+0] ? 0 : imin[0]"], - AssignVariable[arrayElement["bmin", 1], "is_physbnd[1*2+0] ? 0 : imin[1]"], - AssignVariable[arrayElement["bmin", 2], "is_physbnd[2*2+0] ? 0 : imin[2]"], - AssignVariable[arrayElement["bmax", 0], "is_physbnd[0*2+1] ? CCTK_LSSH(0,0) : imax[0]"], - AssignVariable[arrayElement["bmax", 1], "is_physbnd[1*2+1] ? CCTK_LSSH(0,1) : imax[1]"], - AssignVariable[arrayElement["bmax", 2], "is_physbnd[2*2+1] ? CCTK_LSSH(0,2) : imax[2]"]}], - - CommentedBlock["Loop over all faces", - loopOverInteger["dir", "0", "3", - loopOverInteger["face", "0", "2", - { - CommentedBlock["Now restrict to only the boundary points on the current face", - SwitchStatement["face", - {0, {AssignVariable[arrayElement["bmax", "dir"], {arrayElement["imin", "dir"], ""}], - AssignVariable[arrayElement["bmin", "dir"], {0, ""}]}}, - {1, {AssignVariable[arrayElement["bmin", "dir"], {arrayElement["imax", "dir"], "" }], - AssignVariable[arrayElement["bmax", "dir"], {"CCTK_LSSH(0,dir)", ""}]}}]], - conditional[arrayElement["is_physbnd", "dir * 2 + face"], - 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 - } - - ]]] - ]} - ]]]}; -ErrorDefinition[BoundaryLoop]; - -BoundaryWithGhostsLoop[block:CodeGenBlock] := -{ - "\nGenericFD_GetBoundaryInfo(cctkGH, cctk_lsh, cctk_lssh, cctk_bbox, cctk_nghostzones, imin, imax, is_symbnd, is_physbnd, is_ipbnd);\n", - - CommentedBlock["Start by looping over the whole grid, including the NON-PHYSICAL boundary points. ", { - AssignVariable[arrayElement["bmin", 0], "0"], - AssignVariable[arrayElement["bmin", 1], "0"], - AssignVariable[arrayElement["bmin", 2], "0"], - AssignVariable[arrayElement["bmax", 0], "CCTK_LSSH(0,0)"], - AssignVariable[arrayElement["bmax", 1], "CCTK_LSSH(0,1)"], - AssignVariable[arrayElement["bmax", 2], "CCTK_LSSH(0,2)"]}], - - CommentedBlock["Loop over all faces", - loopOverInteger["dir", "0", "3", - loopOverInteger["face", "0", "2", - { - CommentedBlock["Now restrict to only the boundary points on the current face", - SwitchStatement["face", - {0, {AssignVariable[arrayElement["bmax", "dir"], {arrayElement["imin", "dir"], ""}], - AssignVariable[arrayElement["bmin", "dir"], {0, ""}]}}, - {1, {AssignVariable[arrayElement["bmin", "dir"], {arrayElement["imax", "dir"], "" }], - AssignVariable[arrayElement["bmax", "dir"], {"CCTK_LSSH(0,dir)]", ""}]}}]], - conditional[arrayElement["is_physbnd", "dir * 2 + face"], - 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 - } - - ]]] - ]} - ]]]}; -ErrorDefinition[BoundaryWithGhostsLoop]; - -Conditional[condition:CodeGenBlock, block:CodeGenBlock] := - {"if (", condition, ")\n", - CBlock[block]}; -ErrorDefinition[Conditional]; - -Conditional[condition:CodeGenBlock, block1:CodeGenBlock, block2:CodeGenBlock] := - {"if (", condition, ")\n", - CBlock[block1], "else\n", CBlock[block2]}; -ErrorDefinition[Conditional]; - -onceInGridLoop[block:CodeGenBlock] := - Conditional["i == 5 && j == 5 && k == 5", - block]; -ErrorDefinition[onceInGridLoop]; - -InfoVariable[name_String] := - onceInGridLoop[ - { "char buffer[255];\n", - "sprintf(buffer,\"" , name , " == %f\", " , name , ");\n", - "CCTK_INFO(buffer);\n"}]; -ErrorDefinition[InfoVariable]; - -(* Code generation for Cactus .par files *) - -activeThorns[list:{_String...}] := - {"ActiveThorns = \"", SpaceSeparated[list], "\"\n"}; -ErrorDefinition[activeThorns]; - -setParameter[thorn_String, par_String, value_] := - {thorn, " = ", If[NumberQ[value], ToString[value], "\"" <> value <> "\""], "\n"}; -ErrorDefinition[setParameter]; - -setParametersForThorn[thorn_String, map_List] := - Map[setParameter[thorn, #[[1]], #[[2]]] &, map]; -ErrorDefinition[setParametersForThorn]; - -insertFile[name_String] := - Module[{istream_, contents_}, +insertFile[name_String] := + Module[ + {istream_, contents_}, istream = OpenRead[name]; contents = ReadList[istream, String]; Close[istream]; contents]; ErrorDefinition[insertFile]; -vectoriseExpression[exprp_] := - Module[{isNotMinusOneQ, isNotTimesMinusOneQ, fmaRules, isNotKneg, arithRules, undoRules, expr}, - expr = exprp; - - (* Constants *) - (* expr = expr /. xx_Integer/; xx!=-1 :> ToReal[xx]; *) - expr = expr /. xx_Integer -> ToReal[xx]; - expr = expr /. xx_Real -> ToReal[xx]; - removeToRealRules = { - - ToReal[xx_] -> ToReal[- xx], - ToReal[xx_] + ToReal[yy_] -> ToReal[xx + yy], - ToReal[xx_] - ToReal[yy_] -> ToReal[xx - yy], - ToReal[xx_] * ToReal[yy_] -> ToReal[xx * yy], - ToReal[xx_] == ToReal[yy_] -> ToReal[xx == yy], - ToReal[xx_] != ToReal[yy_] -> ToReal[xx != yy], - pow[xx_, ToReal[power_]] -> pow[xx, power] - }; - expr = expr //. removeToRealRules; - - (* Replace all operators and functions *) - (* kneg, kadd, ksub, kmul, kdiv *) - isNotKneg[n_] := ! MatchQ[n,kneg[_]]; - arithRules = { - - xx_ -> kneg[xx], - xx_ * yy_ -> kmul[xx,yy], - xx_ / yy_ -> kdiv[xx,yy], - xx_ + yy_ -> kadd[xx,yy], - xx_ - yy_ -> ksub[xx,yy], - kmul[-1,xx_] -> kneg[xx], - kmul[xx_,-1] -> kneg[xx], - kmul[ToReal[-1],xx_] -> kneg[xx], - kmul[xx_,ToReal[-1]] -> kneg[xx], - ToReal[- xx_] -> kneg[ToReal[xx]], - (* kmul[xx_,INV[yy_]] -> kdiv[xx,yy], *) - (* kmul[INV[xx_],yy_] -> kdiv[yy,xx], *) - kdiv[xx_,kdiv[yy_,zz_]] -> kdiv[kmul[xx,zz],yy], - kdiv[kdiv[xx_,yy_],zz_] -> kdiv[xx,kmul[yy,zz]], - kmul[kneg[xx_],yy_] -> kneg[kmul[xx,yy]], - kmul[xx_,kneg[yy_]] -> kneg[kmul[xx,yy]], - kdiv[kneg[xx_],yy_] -> kneg[kdiv[xx,yy]], - kdiv[xx_,kneg[yy_]] -> kneg[kdiv[xx,yy]], - kadd[kneg[xx_],yy_] -> ksub[yy,xx], - ksub[kneg[xx_],yy_] -> kneg[kadd[xx,yy]], - kadd[xx_,kneg[yy_]] -> ksub[xx,yy], - ksub[xx_,kneg[yy_]] -> kadd[xx,yy], - kneg[ksub[xx_,yy_]] -> ksub[yy,xx], - Abs[xx_] -> kfabs[xx], - Cos[xx_] -> kcos[xx], - Log[xx_] -> klog[xx], - Sin[xx_] -> ksin[xx], - Tan[xx_] -> ktan[xx], - exp[xx_] -> kexp[xx], - fabs[xx_] -> kfabs[xx], - fmax[xx_,yy_] -> kfmax[xx,yy], - fmin[xx_,yy_] -> kfmin[xx,yy], - log[xx_] -> klog[xx], - pow[xx_,yy_] -> kpow[xx,yy], - sqrt[xx_] -> ksqrt[xx], - kcos[kneg[xx_]] -> kcos[xx], - kfabs[kneg[xx_]] -> kfabs[xx], - kfnabs[kneg[xx_]] -> kfnabs[xx], - kneg[kfabs[xx_]] -> kfnabs[xx], - kneg[kfnabs[xx_]] -> kfabs[xx], - kneg[kneg[xx_]] -> xx, - ksin[kneg[xx_]] -> kneg[ksin[xx]], - ktan[kneg[xx_]] -> kneg[ktan[xx]] - }; - expr = expr //. arithRules; - - (* Undo some transformations *) - undoRules = { - IfThen[_, aa_, aa_] -> aa, - IfThen[xx_?IntegerQ, aa_, bb_] /; xx!=0 :> aa, - IfThen[xx_?IntegerQ, aa_, bb_] /; xx==0 :> bb, - IfThen[kmul[xx_,yy_], aa_, bb_] -> IfThen[xx*yy, aa, bb], - IfThen[kmul[xx_,yy_] != zz_, aa_, bb_] -> IfThen[xx*yy!=zz, aa, bb], - IfThen[ToReal[xx_], aa_, bb_] -> IfThen[xx, aa, bb], - Scalar[kmul[xx_,yy_]] -> Scalar[xx*yy], - Scalar[kmul[xx_,yy_] != zz_] -> Scalar[xx*yy!=zz], - Scalar[ToReal[xx_]] -> Scalar[xx], - Scalar[xx_ != ToReal[yy_]] -> Scalar[xx != yy], - ToReal[kneg[xx_]] -> ToReal[-xx], - ToReal[kadd[xx_,yy_]] -> ToReal[xx+yy], - ToReal[ksub[xx_,yy_]] -> ToReal[xx-yy], - ToReal[kmul[xx_,yy_]] -> ToReal[xx*yy], - ToReal[xx_*kadd[yy_,zz_]] -> ToReal[xx*(yy+zz)], - kpow[xx_, kneg[power_]] -> kpow[xx, -power] - }; - expr = expr //. undoRules; - - (* FMA (fused multiply-add) instructions *) - (* kmadd (x,y,z) = xy+z - kmsub (x,y,z) = xy-z - knmadd(x,y,z) = -(xy+z) - knmsub(x,y,z) = -(xy-z) *) - fmaRules = { - kadd[kmul[xx_,yy_],zz_] -> kmadd[xx,yy,zz], - kadd[zz_,kmul[xx_,yy_]] -> kmadd[xx,yy,zz], - ksub[kmul[xx_,yy_],zz_] -> kmsub[xx,yy,zz], - ksub[zz_,kmul[xx_,yy_]] -> knmsub[xx,yy,zz], - kneg[kmadd [xx_,yy_,zz_]] -> knmadd[xx,yy,zz], - kneg[kmsub [xx_,yy_,zz_]] -> knmsub[xx,yy,zz], - kneg[knmadd[xx_,yy_,zz_]] -> kmadd [xx,yy,zz], - kneg[knmsub[xx_,yy_,zz_]] -> kmsub [xx,yy,zz] - (* we could match this and similar patterns - kmul[xx_, kadd[yy_, ToReal[+1]]] -> kmadd[xx, yy, xx], - kmul[xx_, kadd[yy_, ToReal[-1]]] -> kmsub[xx, yy, xx], - *) - }; - expr = expr //. fmaRules; - - Return[expr]]; -ErrorDefinition[vectoriseExpression]; - -(* Take an expression x and replace occurrences of Powers with the C -macros SQR, CUB, QAD *) -ReplacePowers[expr_, vectorise:Boolean] := - Module[{rhs}, - rhs = expr /. Power[xx_, -1] -> INV[xx]; - - If[SOURCELANGUAGE == "C", - Module[{}, - 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 /. 1/2 -> khalf - rhs = rhs /. -1/2 -> -khalf; - - rhs = rhs /. 1/3 -> kthird; - rhs = rhs /. -1/3 -> -kthird; - - rhs = rhs /. 2/3 -> ktwothird; - rhs = rhs /. -2/3 -> -ktwothird; - - rhs = rhs /. 4/3 -> kfourthird; - rhs = rhs /. -4/3 -> -kfourthird; - - rhs = rhs /. 8/3 -> keightthird; - rhs = rhs /. -8/3 -> -keightthird; - *) - - (* Remove parentheses *) - rhs = rhs //. Parenthesis[xx_] -> xx; - - (* Avoid rational numbers *) - rhs = rhs /. Rational[xx_,yy_] :> N[xx/yy, 30]; - - rhs = rhs //. IfThen[cond1_,xx1_,yy1_] + IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> IfThen[cond1, Simplify[ xx1 + xx2], Simplify[ yy1 + yy2]]; - rhs = rhs //. ff1_ IfThen[cond1_,xx1_,yy1_] + IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> IfThen[cond1, Simplify[ff1 xx1 + xx2], Simplify[ff1 yy1 + yy2]]; - rhs = rhs //. IfThen[cond1_,xx1_,yy1_] + ff2_ IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> IfThen[cond1, Simplify[ xx1 + ff2 xx2], Simplify[ yy1 + ff2 yy2]]; - rhs = rhs //. ff1_ IfThen[cond1_,xx1_,yy1_] + ff2_ IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> IfThen[cond1, Simplify[ff1 xx1 + ff2 xx2], Simplify[ff1 yy1 + ff2 yy2]]; - - (* Is this still a good idea when FMA instructions are used? *) - rhs = rhs //. xx_ yy_ + xx_ zz_ -> xx (yy+zz); - rhs = rhs //. xx_ yy_ - xx_ zz_ -> xx (yy-zz); - - rhs = rhs /. Power[E, power_] -> exp[power]; - - (* 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 /. Power[xx_, power_] -> pow[xx, power]; - - If[vectorise === True, - rhs = vectoriseExpression[rhs]]; - - (* Remove Scalar[] after vectorising *) - rhs = rhs /. Scalar[xx_] -> xx; - ], - - rhs = rhs /. Power[xx_, power_] -> xx^power - ]; -(* Print[rhs//FullForm];*) - rhs - ]; -ErrorDefinition[ReplacePowers]; - -(* Convert an expression to CForm, but remove the quotes from any - strings present *) -CFormHideStrings[x_, opts___] := StringReplace[ToString[CForm[x,opts]], "\"" -> ""]; -ErrorDefinition[CFormHideStrings]; - - +NameRoot[name_Symbol] := + Module[ + {dropNumberRule, root}, + dropNumberRule = {"1" -> "", "2" -> "", "3" -> "", "4" -> "", "5" -> "", + "6" -> "", "7" -> "", "8" -> "", "9" -> "", "0" -> "", "rhs" -> ""}; + root = StringReplace[ToString@name, dropNumberRule]]; +ErrorDefinition[NameRoot]; -Quote[x:CodeGenBlock] := {"\"", x, "\""}; +Quote[x:CodeGenBlock] := + {"\"", x, "\""}; ErrorDefinition[Quote]; End[]; diff --git a/Tools/CodeGen/CodeGenC.m b/Tools/CodeGen/CodeGenC.m new file mode 100644 index 0000000..d19f616 --- /dev/null +++ b/Tools/CodeGen/CodeGenC.m @@ -0,0 +1,243 @@ + +(* $Id$ *) + +(* Copyright 2004 Sascha Husa, Ian Hinder, Christiane Lechner + + 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["CodeGenC`", {"Errors`", "Kranc`", "CodeGen`"}]; + +SOURCELANGUAGE::usage = "global variable == \"C\" or \"Fortran\" determines language + for code generation"; +SOURCESUFFIX::usage = "file suffix for source files"; +SetSourceLanguage::usage = "set the source language to \"C\" or \"Fortran\""; + +EOL::usage = "the end of line termination string"; + +CommentedBlock::usage = "CommentedBlock[comment, block] returns a block consisting " <> + "of 'comment' followed by 'block'."; +CBlock::usage = ""; +SuffixedCBlock::usage = ""; +IncludeFile::usage = "IncludeFile[name] returns a block of code" <> + "that includes a header file (i.e '#include \"name\"')."; +IncludeSystemFile::usage = "IncludeFile[name] returns a block of code" <> + "that includes a system header file (i.e '#include <name>')."; +DeclareVariable::usage = "DeclareVariable[name, type] returns a block of code " <> + "that declares a variable of given name and type. 'name' and 'type' should be " <> + "strings."; +DeclareVariableNoInit::usage = "DeclareVariableNoInit[name, type] returns a block of code " <> + "that declares a variable of given name and type without initialising it. 'name' and 'type' should be " <> + "strings."; +DeclareVariables::usage = "DeclareVariables[names, type] returns a block of code " <> + "that declares a list of variables of given name and type. 'names' should be a list" <> + " of strings and 'type' should be a string string."; +DeclarePointer::usage = "DeclarePointer[name, type] returns a block of code " <> + "that declares a pointer of given name and type. 'name' and 'type' should be " <> + "strings."; +DeclarePointers::usage = "DeclarePointers[names, type] returns a block of code " <> + "that declares a list of pointers of given name and type. 'names' should be a list" <> + " of strings and 'type' should be a string string."; +DefineVariable::usage = "DefineVariable[name, type, value] returns a block of " <> + "code that declares and initialised a variable 'name' of type 'type' to value 'value'."; +AssignVariable::usage = "AssignVariable[dest_, src_] returns a block of code " <> + "that assigns 'src' to 'dest'."; +DeclareAssignVariable::usage = "DeclareAssignVariable[type_, dest_, src_] returns a block of code " <> + "that declares and sets a constant variable of given name and type."; +DeclareArray::usage = ""; +DefineFunction::usage = ""; +DefineSubroutine::usage = ""; +Conditional::usage = ""; +SwitchStatement::usage = ""; +CFormHideStrings::usage = ""; + +Begin["`Private`"]; + +SOURCELANGUAGE = "C"; +SOURCESUFFIX = ".cc"; + +setSourceSuffix[lang_] := + If[lang == "C", SOURCESUFFIX = ".cc", SOURCESUFFIX = ".F90"]; + +SetSourceLanguage[lang_] := + If[lang == "C" || lang == "Fortran", + SOURCELANGUAGE = lang; + setSourceSuffix[lang]; + InfoMessage[Terse, "User set source language to " <> lang], + (* else *) + SOURCELANGUAGE = "C"; + setSourceSuffix[".cc"]; + InfoMessage[Terse, "Setting Source Language to C"]]; + +EOL[dummy___] := + If[SOURCELANGUAGE == "C" || SOURCELANGUAGE == "C++", ";\n", "\n"]; + +IncludeFile[filename_String] := + {"#include \"", filename, "\"\n"}; +ErrorDefinition[IncludeFile]; + +IncludeSystemFile[filename_String] := + {"#include <", filename, ">\n"}; +ErrorDefinition[IncludeSystemFile]; + +DeclareVariable[name:(_String|_Symbol), type_String] := + If[SOURCELANGUAGE == "C", + {type, " ", name, " = INITVALUE" <> EOL[]}, + {type, " :: ", name, EOL[]} (* no value init here to avoid implicit SAVE attribute *)]; +ErrorDefinition[DeclareVariable]; + +DeclareVariableNoInit[name:(_String|_Symbol), type_String] := + If[SOURCELANGUAGE == "C", + {type, " ", name, EOL[]}, + {type, " :: ", name, EOL[]} (* no value init here to avoid implicit SAVE attribute *)]; +ErrorDefinition[DeclareVariableNoInit]; + +DeclareVariables[names_?ListQ, type_String] := + If[SOURCELANGUAGE == "C", + {type, " ", CommaSeparated@names, EOL[]}, + {type, " :: ", CommaSeparated@names, EOL[]} (* no value init avoids implicit SAVE attribute *)]; +ErrorDefinition[DeclareVariables]; + +DeclarePointer[name:(_String|_Symbol), type_String] := + If[SOURCELANGUAGE == "C", + {type, " *", name, EOL[]}, + {type, ", target :: ", name, EOL[]}]; +ErrorDefinition[DeclarePointer]; + +DeclarePointers[names_?ListQ, type_String] := + If[SOURCELANGUAGE == "C", + {type, " *", CommaInitSeparated@names, EOL[]}, + {type, ", target :: ", CommaSeparated@names, EOL[]}]; +ErrorDefinition[DeclarePointers]; + +DeclareArray[name:(_String|_Symbol), dim_Integer, type_String] := + If[SOURCELANGUAGE == "C", + DeclareArrayC[name, dim, type], + DeclareArrayFortran[name, dim, type]]; +ErrorDefinition[DeclareArray]; + +DeclareArrayC[name:(_String|_Symbol), dim_Integer, type_String] := + {type, " ", name, "[", dim, "];","\n"}; +ErrorDefinition[DeclareArrayC]; + +DeclareArrayFortran[name:(_String|_Symbol), dim_Integer, type_String] := + {type, " :: ", name, "(", dim, ")","\n"}; +ErrorDefinition[DeclareArrayFortran]; + +DefineVariable[name:(_String|_Symbol), type_String, value:CodeGenBlock] := + {type, " ", name, " = ", value, EOL[]}; +ErrorDefinition[DefineVariable]; + +AssignVariable[dest:(_String|_Symbol), src:CodeGenBlock] := + {dest, " = ", src, EOL[]}; +ErrorDefinition[AssignVariable]; + +DeclareAssignVariable[type_String, dest:(_String|_Symbol), src:CodeGenBlock] := + {type, " const ", dest, " = ", src, EOL[]}; +ErrorDefinition[DeclareAssignVariable]; + +(* comments are always done C-style because they are killed by cpp anyway *) +insertComment[text:CodeGenBlock] := {"/* ", text, " */\n"}; +ErrorDefinition[insertComment]; + +CBlock[block:CodeGenBlock] := + {"{\n", + IndentBlock[block], + "}\n"}; +ErrorDefinition[CBlock]; + +SuffixedCBlock[block:CodeGenBlock, suffix_] := + {"{\n", + IndentBlock[block], + "} ", suffix, "\n"}; +ErrorDefinition[SuffixedCBlock]; + +CommentedBlock[comment:CodeGenBlock, block:CodeGenBlock] := + SeparatedBlock[{insertComment[comment], + block}]; +ErrorDefinition[CommentedBlock]; + +(* FUNCTIONS *) + +defineFunctionC[name_String, type_String, args:CodeGenBlock, contents:CodeGenBlock] := + SeparatedBlock[ + {type, " ", name, "(", args, ")\n", + CBlock[contents]}]; +ErrorDefinition[defineFunctionC]; + +defineFunctionF[name_String, args:CodeGenBlock, contents:CodeGenBlock] := + SeparatedBlock[ + {"FUNCTION", " ", name, "(", args, ")\n", + IndentBlock[contents]}]; +ErrorDefinition[defineFunctionF]; + +DefineFunction[name_String, type_String, args:CodeGenBlock, contents:CodeGenBlock] := + If[SOURCELANGUAGE == "C", + defineFunctionC[name, type, args, contents], + defineFunctionF[name, args, contents]]; +ErrorDefinition[DefineFunction]; + +(* SUBROUTINES *) + +DefineSubroutine[name_String, args:CodeGenBlock, contents:CodeGenBlock] := + If[SOURCELANGUAGE == "C", + DefineSubroutineC[name, args, contents], + DefineSubroutineF[name, args, contents]]; +ErrorDefinition[DefineSubroutine]; + +DefineSubroutineC[name_String, args:CodeGenBlock, contents:CodeGenBlock] := + SeparatedBlock[ + {"extern \"C\" void ", name, "(", args, ")", "\n", + CBlock[contents]}]; +ErrorDefinition[DefineSubroutineC]; + +DefineSubroutineF[name_String, args:CodeGenBlock, contents:CodeGenBlock] := + SeparatedBlock[ + {"subroutine ", name, "(", args, ")", "\n", + "\nimplicit none\n\n", + contents, + "end subroutine\n"}]; +ErrorDefinition[DefineSubroutineF]; + +switchOption[{value:(_String|_Symbol|_?NumberQ), block:CodeGenBlock}] := + {"case ", value, ":\n", IndentBlock[{block,"break;\n"}]}; (* Outer list unnecessary? *) +ErrorDefinition[switchOptions]; + +SwitchStatement[var:(_String|_Symbol), pairs__] := + {"switch(", var, ")\n", + CBlock[{Riffle[Map[switchOption, {pairs}],"\n"]}]}; +ErrorDefinition[SwitchStatement]; + +Conditional[condition:CodeGenBlock, block:CodeGenBlock] := + {"if (", condition, ")\n", + CBlock[block]}; +ErrorDefinition[Conditional]; + +Conditional[condition:CodeGenBlock, block1:CodeGenBlock, block2:CodeGenBlock] := + {"if (", condition, ")\n", + CBlock[block1], "else\n", CBlock[block2]}; +ErrorDefinition[Conditional]; + +(* Convert an expression to CForm, but remove the quotes from any + strings present *) +CFormHideStrings[x_, opts___] := + StringReplace[ToString[CForm[x,opts]], "\"" -> ""]; +ErrorDefinition[CFormHideStrings]; + +End[]; + +EndPackage[]; diff --git a/Tools/CodeGen/CodeGenCactus.m b/Tools/CodeGen/CodeGenCactus.m new file mode 100644 index 0000000..16f410f --- /dev/null +++ b/Tools/CodeGen/CodeGenCactus.m @@ -0,0 +1,727 @@ + +(* $Id$ *) + +(* Copyright 2004 Sascha Husa, Ian Hinder, Christiane Lechner + + 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["CodeGenCactus`", {"Errors`", "Kranc`", "CodeGenC`", "CodeGen`"}]; + +AssignVariableInLoop::usage = "AssignVariableInLoop[dest_, src_] returns a block of code " <> + "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'."; +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 " <> + "that assigns 'src' to 'dest'."; +MaybeAssignVariableInLoop::usage = "MaybeAssignVariableInLoop[dest_, src_, cond_] returns a block of code " <> + "that assigns 'src' to 'dest'."; +DeclareMaybeAssignVariableInLoop::usage = "DeclareMaybeAssignVariableInLoop[type_, dest_, src_, cond_] returns a block of code " <> + "that assigns 'src' to 'dest'."; +TestForNaN::usage = "TestForNaN[expr_] returns a block of code " <> + "that tests 'expr' for nan."; +DefineCCTKFunction::usage = "DefineCCTKFunction[name, type, block] returns a block " <> + "of code that defines a CCTK function of name 'name' returning type 'type' with " <> + "body 'block'."; +DefineCCTKSubroutine::usage = "DefineCCTKSubroutine[name, block] returns a block " <> + "of code that defines a CCTK Fortran subroutine of name 'name' with body 'block'."; +GridName::usage = "GridName[variable] returns the name needed to access variable " <> + "assuming it is a grid variable when inside a grid loop."; +DeclareGridLoopVariables::usage = "DeclareGridLoopVariables[] returns a block " <> + "that defines the variables needed during a grid loop."; +InitialiseGridLoopVariables::usage = "InitialiseGridLoopVariables[] returns a block " <> + "that initialises variables needed by a grid loop."; +InitialiseGridLoopVariablesWithStencil::usage = "InitialiseGridLoopVariables[] returns a block " <> + "that initialises variables needed by a grid loop using the evolution stencils."; +ConditionalOnParameter::usage = "ConditionalOnParameter[name, value, block] returns " <> + "a block that introduces a conditional expression whereby 'block' is only executed " <> + "if the Cactus parameter 'name' has value 'value'."; +(* +GridLoop::usage = "GridLoop[block] returns a block that is looped over for every " <> + "grid point. Must have previously set up the grid loop variables (see " <> + "DeclareGridLoopVariables and InitialiseGridLoopVariables."; +*) +GridLoop::usage = "GridLoop[block] returns a block that is looped over for every " <> + "grid point. Must have previously set up the grid loop variables (see " <> + "InitialiseGridLoopVariables."; +ConditionalOnParameterTextual::usage = ""; +DeclareFDVariables::usage = ""; +InitialiseFDVariables::usage = ""; +ReplacePowers::usage = ""; +BoundaryLoop::usage = ""; +BoundaryWithGhostsLoop::usage = ""; +GenericGridLoop::usage = ""; +NameRoot::usage = ""; +PartitionVarList::usage = ""; +DataType::usage = "DataType[] returns a string for the grid function data type (e.g. CCTK_REAL)"; +SetDataType::usage = "SetDataType[type] sets a string for the grid function data type (e.g. CCTK_REAL)"; + +Begin["`Private`"]; + +SetDataType[type_String] := + dataType = type; + +DataType[] := + If[dataType === Symbol["datatype"], + Throw["DataType: Have not set a data type"], + dataType]; + +AssignVariableInLoop[dest:(_String|_Symbol), src:CodeGenBlock, vectorise_:False] := + Module[ + {loader}, + loader[x_] := If[vectorise, {"vec_load(", x, ")"}, x]; + {dest, " = ", loader[src], EOL[]}]; +ErrorDefinition[AssignVariableInLoop]; + +StoreVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol)] := + {"vec_store_nta(", dest, ",", src, ")", EOL[]}; +ErrorDefinition[StoreVariableInLoop]; + +StoreLowPartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), count_String] := + {"vec_store_nta_partial_lo(", dest, ",", src, ",", count, ")", EOL[]}; +ErrorDefinition[StoreLowPartialVariableInLoop]; + +StoreHighPartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), count_String] := + {"vec_store_nta_partial_hi(", dest, ",", src, ",", count, ")", EOL[]}; +ErrorDefinition[StoreHighPartialVariableInLoop]; + +StoreMiddlePartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), countLow_String, countHigh_String] := + {"vec_store_nta_partial_mid(", dest, ",", src, ",", countLow, ",", countHigh, ")", EOL[]}; +ErrorDefinition[StoreMiddlePartialVariableInLoop]; + +StorePartialVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol)] := + {"vec_store_nta_partial(", dest, ",", src, ")", EOL[]}; +ErrorDefinition[StorePartialVariableInLoop]; + +DeclareAssignVariableInLoop[type_String, dest:(_String|_Symbol), src:(_String|_Symbol)] := + {type, " const ", dest, " = vec_load(", src, ")", EOL[]}; +ErrorDefinition[DeclareAssignVariableInLoop]; + +MaybeAssignVariableInLoop[dest:(_String|_Symbol), src:(_String|_Symbol), cond:Boolean] := + If[cond, + {dest, " = useMatter ? vec_load(", src, ") : ToReal(0.0)", EOL[]}, + {dest, " = vec_load(", src, ")", EOL[]}]; +ErrorDefinition[MaybeAssignVariableInLoop]; + +DeclareMaybeAssignVariableInLoop[type_String, dest:(_String|_Symbol), src:(_String|_Symbol), + mmaCond:Boolean, codeCond:CodeGenBlock, + vectorise:Boolean:False] := + Module[ + {loader}, + loader[x_] := If[vectorise, {"vec_load(", x, ")"}, x]; + If[mmaCond, + {type, " ", dest, " = (", codeCond, ") ? ", loader[src], " : ToReal(0.0)", EOL[]}, + {type, " ", dest, " = ", loader[src], EOL[]}]]; +ErrorDefinition[DeclareMaybeAssignVariableInLoop]; + +TestForNaN[expr:CodeGenBlock] := + {"if (isnan(", expr, ")) {\n", + " CCTK_VInfo(CCTK_THORNSTRING, \"NaN found\");\n", + " CCTK_VInfo(CCTK_THORNSTRING, \"ipos: %d %d %d\", i, j, k);\n", + " CCTK_VInfo(CCTK_THORNSTRING, \"lbnd: %d %d %d\", cctk_lbnd[0], cctk_lbnd[1], cctk_lbnd[2]);\n", + " CCTK_VInfo(CCTK_THORNSTRING, \"lsh: %d %d %d\", cctk_lsh[0], cctk_lsh[1], cctk_lsh[2]);\n", + " CCTK_VInfo(CCTK_THORNSTRING, \"LSSH: %d %d %d\", CCTK_LSSH(0,0), CCTK_LSSH(0,1), CCTK_LSSH(0,2));\n", + " CCTK_VInfo(CCTK_THORNSTRING, \"", expr, ": %.17g\", (double)", expr, ");\n", + "}\n"}; +ErrorDefinition[TestForNaN]; + +loopOverInteger[name_String, start_String, endplusone_String, block_CodeGenBlock] := + If[SOURCELANGUAGE == "C" || SOURCELANGUAGE == "C++", + {"for (", name, " = ", start, "; ", name, " < ", endplusone, "; ", name, "++)\n", + "{\n", + IndentBlock[block], + "}\n"}, + + {"Do ", name, " = ", start, ", ", endplusone, "\n", + "\n", + IndentBlock[block], + "End Do\n"} + ]; +ErrorDefinition[loopOverInteger]; + +(* This is a Cactus-callable function *) +DefineCCTKFunction[name_String, type_String, contents:CodeGenBlock] := + DefineFunction[ + name, "extern \"C\" " <> type, "CCTK_ARGUMENTS", + { + "DECLARE_CCTK_ARGUMENTS;\n", + "DECLARE_CCTK_PARAMETERS;\n\n", + contents + }]; +ErrorDefinition[DefineCCTKFunction]; + +(* This is a Cactus-callable subroutine *) +DefineCCTKSubroutine[name_String, contents:CodeGenBlock] := + DefineSubroutine[ + name, "CCTK_ARGUMENTS", + { + "DECLARE_CCTK_ARGUMENTS;\n", + "DECLARE_CCTK_PARAMETERS;\n\n", + contents + }]; +ErrorDefinition[DefineCCTKSubroutine]; + +DeclareFDVariables[] := +(* + CommentedBlock["Declare finite differencing variables", + {Map[DeclareVariables[#, "CCTK_REAL"] &, {{"dx", "dy", "dz"}, + {"dxi", "dyi", "dzi"}, + {khalf,kthird,ktwothird,kfourthird,keightthird}, + {"hdxi", "hdyi", "hdzi"}}], + "\n"}, + {Map[DeclareVariables[#, "ptrdiff_t"] &, {{"di", "dj", "dk"}}], + "\n"}]; +*) + CommentedBlock["Declare finite differencing variables", {}]; + +InitialiseFDSpacingVariablesC[] := + { + (* DeclareAssignVariable["ptrdiff_t", "di", "CCTK_GFINDEX3D(cctkGH,1,0,0) - CCTK_GFINDEX3D(cctkGH,0,0,0)"], *) + DeclareAssignVariable["ptrdiff_t", "di", "1"], + DeclareAssignVariable["ptrdiff_t", "dj", "CCTK_GFINDEX3D(cctkGH,0,1,0) - CCTK_GFINDEX3D(cctkGH,0,0,0)"], + DeclareAssignVariable["ptrdiff_t", "dk", "CCTK_GFINDEX3D(cctkGH,0,0,1) - CCTK_GFINDEX3D(cctkGH,0,0,0)"], + DeclareAssignVariable["ptrdiff_t", "cdi", "sizeof(CCTK_REAL) * di"], + DeclareAssignVariable["ptrdiff_t", "cdj", "sizeof(CCTK_REAL) * dj"], + DeclareAssignVariable["ptrdiff_t", "cdk", "sizeof(CCTK_REAL) * dk"], + DeclareAssignVariable[DataType[], "dx", "ToReal(CCTK_DELTA_SPACE(0))"], + 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)"] + }; + +InitialiseFDSpacingVariablesFortran[] := + { + AssignVariable["dt", "CCTK_DELTA_TIME"], + AssignVariable["dx", "CCTK_DELTA_SPACE(1)"], + AssignVariable["dy", "CCTK_DELTA_SPACE(2)"], + AssignVariable["dz", "CCTK_DELTA_SPACE(3)"] + }; + +InitialiseFDVariables[vectorise:Boolean] := + CommentedBlock[ + "Initialise finite differencing variables", + { + If[SOURCELANGUAGE == "Fortran", + InitialiseFDSpacingVariablesFortran[], + InitialiseFDSpacingVariablesC[]], + + DeclareAssignVariable[DataType[], "dxi", "INV(dx)"], + DeclareAssignVariable[DataType[], "dyi", "INV(dy)"], + DeclareAssignVariable[DataType[], "dzi", "INV(dz)"], + If[vectorise, + {DeclareAssignVariable[DataType[], "khalf", "ToReal(0.5)"], + DeclareAssignVariable[DataType[], "kthird", "ToReal(1.0/3.0)"], + DeclareAssignVariable[DataType[], "ktwothird", "ToReal(2.0/3.0)"], + DeclareAssignVariable[DataType[], "kfourthird", "ToReal(4.0/3.0)"], + DeclareAssignVariable[DataType[], "keightthird", "ToReal(8.0/3.0)"], + DeclareAssignVariable[DataType[], "hdxi", "kmul(ToReal(0.5), dxi)"], + DeclareAssignVariable[DataType[], "hdyi", "kmul(ToReal(0.5), dyi)"], + DeclareAssignVariable[DataType[], "hdzi", "kmul(ToReal(0.5), dzi)"]}, + (* else *) + {DeclareAssignVariable[DataType[], "khalf", "0.5"], + DeclareAssignVariable[DataType[], "kthird", "1/3.0"], + DeclareAssignVariable[DataType[], "ktwothird", "2.0/3.0"], + DeclareAssignVariable[DataType[], "kfourthird", "4.0/3.0"], + DeclareAssignVariable[DataType[], "keightthird", "8.0/3.0"], + DeclareAssignVariable[DataType[], "hdxi", "0.5 * dxi"], + DeclareAssignVariable[DataType[], "hdyi", "0.5 * dyi"], + DeclareAssignVariable[DataType[], "hdzi", "0.5 * dzi"]}]}]; +ErrorDefinition[InitialiseFDVariables]; + +GridName[x_Symbol] := + If[SOURCELANGUAGE == "C", + ToString[x] <> "[index]", + ToString[x] <> "(i,j,k)"]; +ErrorDefinition[GridName]; + +DeclareGridLoopVariables[] := + SeparatedBlock[ + {insertComment["Declare the variables used for looping over grid points"], + Map[DeclareVariables[#, "CCTK_INT"] &, + {{"i", "j", "k"} + (*, {"istart", "jstart", "kstart"}, + {"iend", "jend", "kend"}, + {"index_offset_x", "index_offset_y", "index_offset_z", "dir", "face"} *)}] + (*, Map[DeclareArray[#, 6, "CCTK_INT"] &, {"is_symbnd", "is_physbnd", "is_ipbnd"}], + Map[DeclareArray[#, 3, "CCTK_INT"] &, {"imin", "imax", "bmin", "bmax"}] *), + If[SOURCELANGUAGE == "C", DeclareVariable["index", "// CCTK_INT"], "\n"]}]; + +(* Access an element of an array; syntax is different between C and + Fortran. Always give this function a C-style array index. *) +arrayElement[var_String, i_String] := + If[SOURCELANGUAGE == "C", + {var, "[", arrayIndex[i], "]"}, + {var, "(", arrayIndex[i], ")"}]; +ErrorDefinition[arrayElement]; + +(* Given a C-style variable index, return the corresponding index for + the language currently in use. The idea is that the caller does not + need to know what language is being used. *) +arrayIndex[i:(_Integer|_String|_Symbol)] := + If[SOURCELANGUAGE == "C", + i, + If[NumberQ[i], i+1, {i, " + 1"}]]; +ErrorDefinition[arrayIndex]; + +max[]:= If[SOURCELANGUAGE == "C", "IMAX", "max"]; + +InitialiseGridLoopVariables[derivativesUsedSwitch:Boolean, addToStencilWidth_Integer] := + CommentedBlock[ + "Set up variables used in the grid loop for the physical grid points", + + If[derivativesUsedSwitch, + {AssignVariable["index_offset_x", max[] <>"(stencil_width, stencil_width_x) + " <> ToString[addToStencilWidth]], + AssignVariable["index_offset_y", max[] <>"(stencil_width, stencil_width_y) + " <> ToString[addToStencilWidth]], + AssignVariable["index_offset_z", max[] <>"(stencil_width, stencil_width_z) + " <> ToString[addToStencilWidth]], + "\n", + AssignVariable["istart", arrayIndex["index_offset_x"]], + AssignVariable["jstart", arrayIndex["index_offset_y"]], + AssignVariable["kstart", arrayIndex["index_offset_z"]], + "\n", + AssignVariable["iend", {"CCTK_LSSH(0,0) - index_offset_x"}], + AssignVariable["jend", {"CCTK_LSSH(0,1) - index_offset_y"}], + AssignVariable["kend", {"CCTK_LSSH(0,2) - index_offset_z"}]}, + (* else *) + {AssignVariable["istart", arrayIndex[0]], + AssignVariable["jstart", arrayIndex[0]], + AssignVariable["kstart", arrayIndex[0]], + "\n", + AssignVariable["iend", "CCTK_LSSH(0,0)"], + AssignVariable["jend", "CCTK_LSSH(0,1)"], + AssignVariable["kend", "CCTK_LSSH(0,2)"]}]]; +ErrorDefinition[InitialiseGridLoopVariables]; + +ConditionalOnParameter[name_String, value_String, block:CodeGenBlock] := + SeparatedBlock[ + {"if (CCTK_EQUALS(", name, ", \"", value, "\"))\n", + "{\n", + IndentBlock[block], + "}\n"}]; +ErrorDefinition[ConditionalOnParameter]; + +ConditionalOnParameterTextual[text:CodeGenBlock, block:CodeGenBlock] := + SeparatedBlock[ + {"if (", text, ")\n", + "{\n", + IndentBlock[block], + "}\n"}]; +ErrorDefinition[ConditionalOnParameterTextual]; + +(* +GridLoop[block_] := + CommentedBlock["Loop over the grid points", + loopOverInteger["k", "kstart", "kend", + loopOverInteger["j", "jstart", "jend", + loopOverInteger["i", "istart", "iend", + + { If[SOURCELANGUAGE == "C", + AssignVariable["index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"], + ""], + block + } + ]]]]; + +*) + +(* +GridLoop[block_] := + If[SOURCELANGUAGE == "C", + CommentedBlock["Loop over the grid points", + { + "#pragma omp parallel\n", + "LC_LOOP3 (unnamed,\n", + " i,j,k, istart,jstart,kstart, iend,jend,kend,\n", + " cctk_lsh[0],cctk_lsh[1],cctk_lsh[2])\n", + "{\n", + IndentBlock[ + { + DeclareVariable["index", "int"], + AssignVariable["index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"], + block + } + ], + "}\n", + "LC_ENDLOOP3 (unnamed);\n" + } + ], + CommentedBlock["Loop over the grid points", + { + "#pragma omp parallel\n", + "LC_LOOP3 (unnamed,\n", + " i,j,k, istart,jstart,kstart, iend,jend,kend,\n", + " cctk_lsh(1),cctk_lsh(2),cctk_lsh(3))\n", + IndentBlock[block], + "LC_ENDLOOP3 (unnamed)\n" + } + ] + ]; +*) + +Options[GenericGridLoop] = ThornOptions; + +GenericGridLoop[functionName_String, block:CodeGenBlock, opts:OptionsPattern[]] := + If[OptionValue[UseLoopControl], + GenericGridLoopUsingLoopControl[functionName, block, OptionValue[UseVectors]], + GenericGridLoopTraditional[block]]; +ErrorDefinition[GenericGridLoop]; + +GenericGridLoopTraditional[block:CodeGenBlock] := + CommentedBlock[ + "Loop over the grid points", + loopOverInteger[ + "k", "imin[2]", "imax[2]", + loopOverInteger[ + "j", "imin[1]", "imax[1]", + loopOverInteger[ + "i", "imin[0]", "imax[0]", + {If[SOURCELANGUAGE == "C", + DeclareAssignVariable["int", "index", "CCTK_GFINDEX3D(cctkGH,i,j,k)"], + ""], + block}]]]]; +ErrorDefinition[GenericGridLoopTraditional]; + +GenericGridLoopUsingLoopControl[functionName_String, block:CodeGenBlock, vectorise:Boolean] := + If[SOURCELANGUAGE == "C", + CommentedBlock[ + "Loop over the grid points", + {"#pragma omp parallel\n", + If[vectorise, "LC_LOOP3VEC", "LC_LOOP3"], + " (", 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"}, ""], + ")\n", + "{\n", + IndentBlock[ + {(* DeclareVariable["index", "// int"], *) + (* 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", "LC_ENDLOOP3"] <> " (", functionName, ");\n"}], + (* else *) + ""]; +ErrorDefinition[GenericGridLoopUsingLoopControl]; + +BoundaryLoop[block:CodeGenBlock] := + { + "\nGenericFD_GetBoundaryInfo(cctkGH, cctk_lsh, cctk_lssh, cctk_bbox, cctk_nghostzones, imin, imax, is_symbnd, is_physbnd, is_ipbnd);\n", + + CommentedBlock[ + "Start by looping over the whole grid, minus the NON-PHYSICAL boundary points, which are set by synchronization. ", { + AssignVariable[arrayElement["bmin", 0], "is_physbnd[0*2+0] ? 0 : imin[0]"], + AssignVariable[arrayElement["bmin", 1], "is_physbnd[1*2+0] ? 0 : imin[1]"], + AssignVariable[arrayElement["bmin", 2], "is_physbnd[2*2+0] ? 0 : imin[2]"], + AssignVariable[arrayElement["bmax", 0], "is_physbnd[0*2+1] ? CCTK_LSSH(0,0) : imax[0]"], + AssignVariable[arrayElement["bmax", 1], "is_physbnd[1*2+1] ? CCTK_LSSH(0,1) : imax[1]"], + AssignVariable[arrayElement["bmax", 2], "is_physbnd[2*2+1] ? CCTK_LSSH(0,2) : imax[2]"]}], + + CommentedBlock[ + "Loop over all faces", + loopOverInteger[ + "dir", "0", "3", + loopOverInteger[ + "face", "0", "2", + {CommentedBlock[ + "Now restrict to only the boundary points on the current face", + SwitchStatement[ + "face", + {0, {AssignVariable[arrayElement["bmax", "dir"], {arrayElement["imin", "dir"], ""}], + AssignVariable[arrayElement["bmin", "dir"], {0, ""}]}}, + {1, {AssignVariable[arrayElement["bmin", "dir"], {arrayElement["imax", "dir"], "" }], + AssignVariable[arrayElement["bmax", "dir"], {"CCTK_LSSH(0,dir)", ""}]}}]], + conditional[ + arrayElement["is_physbnd", "dir * 2 + face"], + 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}]]]]}]]]}; +ErrorDefinition[BoundaryLoop]; + +BoundaryWithGhostsLoop[block:CodeGenBlock] := +{ + "\nGenericFD_GetBoundaryInfo(cctkGH, cctk_lsh, cctk_lssh, cctk_bbox, cctk_nghostzones, imin, imax, is_symbnd, is_physbnd, is_ipbnd);\n", + + CommentedBlock[ + "Start by looping over the whole grid, including the NON-PHYSICAL boundary points. ", { + AssignVariable[arrayElement["bmin", 0], "0"], + AssignVariable[arrayElement["bmin", 1], "0"], + AssignVariable[arrayElement["bmin", 2], "0"], + AssignVariable[arrayElement["bmax", 0], "CCTK_LSSH(0,0)"], + AssignVariable[arrayElement["bmax", 1], "CCTK_LSSH(0,1)"], + AssignVariable[arrayElement["bmax", 2], "CCTK_LSSH(0,2)"]}], + + CommentedBlock[ + "Loop over all faces", + loopOverInteger[ + "dir", "0", "3", + loopOverInteger[ + "face", "0", "2", + { + CommentedBlock[ + "Now restrict to only the boundary points on the current face", + SwitchStatement[ + "face", + {0, {AssignVariable[arrayElement["bmax", "dir"], {arrayElement["imin", "dir"], ""}], + AssignVariable[arrayElement["bmin", "dir"], {0, ""}]}}, + {1, {AssignVariable[arrayElement["bmin", "dir"], {arrayElement["imax", "dir"], "" }], + AssignVariable[arrayElement["bmax", "dir"], {"CCTK_LSSH(0,dir)]", ""}]}}]], + conditional[arrayElement["is_physbnd", "dir * 2 + face"], + 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}]]]]}]]]}; +ErrorDefinition[BoundaryWithGhostsLoop]; + +onceInGridLoop[block:CodeGenBlock] := + Conditional[ + "i == 5 && j == 5 && k == 5", + block]; +ErrorDefinition[onceInGridLoop]; + +InfoVariable[name_String] := + onceInGridLoop[ + {"char buffer[255];\n", + "sprintf(buffer,\"" , name , " == %f\", " , name , ");\n", + "CCTK_INFO(buffer);\n"}]; +ErrorDefinition[InfoVariable]; + +(* Code generation for Cactus .par files *) + +activeThorns[list:{_String...}] := + {"ActiveThorns = \"", SpaceSeparated[list], "\"\n"}; +ErrorDefinition[activeThorns]; + +setParameter[thorn_String, par_String, value_] := + {thorn, " = ", If[NumberQ[value], ToString[value], "\"" <> value <> "\""], "\n"}; +ErrorDefinition[setParameter]; + +setParametersForThorn[thorn_String, map_List] := + Map[setParameter[thorn, #[[1]], #[[2]]] &, map]; +ErrorDefinition[setParametersForThorn]; + +vectoriseExpression[exprp_] := + Module[ + {isNotMinusOneQ, isNotTimesMinusOneQ, fmaRules, isNotKneg, arithRules, undoRules, expr}, + expr = exprp; + + (* Constants *) + (* expr = expr /. xx_Integer/; xx!=-1 :> ToReal[xx]; *) + expr = expr /. xx_Integer -> ToReal[xx]; + expr = expr /. xx_Real -> ToReal[xx]; + + removeToRealRules = + {-ToReal[xx_] -> ToReal[- xx], + ToReal[xx_] + ToReal[yy_] -> ToReal[xx + yy], + ToReal[xx_] - ToReal[yy_] -> ToReal[xx - yy], + ToReal[xx_] * ToReal[yy_] -> ToReal[xx * yy], + ToReal[xx_] == ToReal[yy_] -> ToReal[xx == yy], + ToReal[xx_] != ToReal[yy_] -> ToReal[xx != yy], + pow[xx_, ToReal[power_]] -> pow[xx, power]}; + + expr = expr //. removeToRealRules; + + (* Replace all operators and functions *) + (* kneg, kadd, ksub, kmul, kdiv *) + isNotKneg[n_] := ! MatchQ[n,kneg[_]]; + + arithRules = + { - xx_ -> kneg[xx], + xx_ * yy_ -> kmul[xx,yy], + xx_ / yy_ -> kdiv[xx,yy], + xx_ + yy_ -> kadd[xx,yy], + xx_ - yy_ -> ksub[xx,yy], + kmul[-1,xx_] -> kneg[xx], + kmul[xx_,-1] -> kneg[xx], + kmul[ToReal[-1],xx_] -> kneg[xx], + kmul[xx_,ToReal[-1]] -> kneg[xx], + ToReal[- xx_] -> kneg[ToReal[xx]], + (* kmul[xx_,INV[yy_]] -> kdiv[xx,yy], *) + (* kmul[INV[xx_],yy_] -> kdiv[yy,xx], *) + kdiv[xx_,kdiv[yy_,zz_]] -> kdiv[kmul[xx,zz],yy], + kdiv[kdiv[xx_,yy_],zz_] -> kdiv[xx,kmul[yy,zz]], + kmul[kneg[xx_],yy_] -> kneg[kmul[xx,yy]], + kmul[xx_,kneg[yy_]] -> kneg[kmul[xx,yy]], + kdiv[kneg[xx_],yy_] -> kneg[kdiv[xx,yy]], + kdiv[xx_,kneg[yy_]] -> kneg[kdiv[xx,yy]], + kadd[kneg[xx_],yy_] -> ksub[yy,xx], + ksub[kneg[xx_],yy_] -> kneg[kadd[xx,yy]], + kadd[xx_,kneg[yy_]] -> ksub[xx,yy], + ksub[xx_,kneg[yy_]] -> kadd[xx,yy], + kneg[ksub[xx_,yy_]] -> ksub[yy,xx], + Abs[xx_] -> kfabs[xx], + Cos[xx_] -> kcos[xx], + Log[xx_] -> klog[xx], + Sin[xx_] -> ksin[xx], + Tan[xx_] -> ktan[xx], + exp[xx_] -> kexp[xx], + fabs[xx_] -> kfabs[xx], + fmax[xx_,yy_] -> kfmax[xx,yy], + fmin[xx_,yy_] -> kfmin[xx,yy], + log[xx_] -> klog[xx], + pow[xx_,yy_] -> kpow[xx,yy], + sqrt[xx_] -> ksqrt[xx], + kcos[kneg[xx_]] -> kcos[xx], + kfabs[kneg[xx_]] -> kfabs[xx], + kfnabs[kneg[xx_]] -> kfnabs[xx], + kneg[kfabs[xx_]] -> kfnabs[xx], + kneg[kfnabs[xx_]] -> kfabs[xx], + kneg[kneg[xx_]] -> xx, + ksin[kneg[xx_]] -> kneg[ksin[xx]], + ktan[kneg[xx_]] -> kneg[ktan[xx]]}; + expr = expr //. arithRules; + + (* Undo some transformations *) + undoRules = + { IfThen[_, aa_, aa_] -> aa, + IfThen[xx_?IntegerQ, aa_, bb_] /; xx!=0 :> aa, + IfThen[xx_?IntegerQ, aa_, bb_] /; xx==0 :> bb, + IfThen[kmul[xx_,yy_], aa_, bb_] -> IfThen[xx*yy, aa, bb], + IfThen[kmul[xx_,yy_] != zz_, aa_, bb_] -> IfThen[xx*yy!=zz, aa, bb], + IfThen[ToReal[xx_], aa_, bb_] -> IfThen[xx, aa, bb], + Scalar[kmul[xx_,yy_]] -> Scalar[xx*yy], + Scalar[kmul[xx_,yy_] != zz_] -> Scalar[xx*yy!=zz], + Scalar[ToReal[xx_]] -> Scalar[xx], + Scalar[xx_ != ToReal[yy_]] -> Scalar[xx != yy], + ToReal[kneg[xx_]] -> ToReal[-xx], + ToReal[kadd[xx_,yy_]] -> ToReal[xx+yy], + ToReal[ksub[xx_,yy_]] -> ToReal[xx-yy], + ToReal[kmul[xx_,yy_]] -> ToReal[xx*yy], + ToReal[xx_*kadd[yy_,zz_]] -> ToReal[xx*(yy+zz)], + kpow[xx_, kneg[power_]] -> kpow[xx, -power]}; + expr = expr //. undoRules; + + (* FMA (fused multiply-add) instructions *) + (* kmadd (x,y,z) = xy+z + kmsub (x,y,z) = xy-z + knmadd(x,y,z) = -(xy+z) + knmsub(x,y,z) = -(xy-z) *) + fmaRules = + { kadd[kmul[xx_,yy_],zz_] -> kmadd[xx,yy,zz], + kadd[zz_,kmul[xx_,yy_]] -> kmadd[xx,yy,zz], + ksub[kmul[xx_,yy_],zz_] -> kmsub[xx,yy,zz], + ksub[zz_,kmul[xx_,yy_]] -> knmsub[xx,yy,zz], + kneg[kmadd [xx_,yy_,zz_]] -> knmadd[xx,yy,zz], + kneg[kmsub [xx_,yy_,zz_]] -> knmsub[xx,yy,zz], + kneg[knmadd[xx_,yy_,zz_]] -> kmadd [xx,yy,zz], + kneg[knmsub[xx_,yy_,zz_]] -> kmsub [xx,yy,zz] + (* we could match this and similar patterns + kmul[xx_, kadd[yy_, ToReal[+1]]] -> kmadd[xx, yy, xx], + kmul[xx_, kadd[yy_, ToReal[-1]]] -> kmsub[xx, yy, xx], + *)}; + expr = expr //. fmaRules; + + Return[expr]]; +ErrorDefinition[vectoriseExpression]; + +(* Take an expression x and replace occurrences of Powers with the C +macros SQR, CUB, QAD *) +ReplacePowers[expr_, vectorise:Boolean] := + Module[ + {rhs}, + rhs = expr /. Power[xx_, -1] -> INV[xx]; + 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 /. 1/2 -> khalf + rhs = rhs /. -1/2 -> -khalf; + + rhs = rhs /. 1/3 -> kthird; + rhs = rhs /. -1/3 -> -kthird; + + rhs = rhs /. 2/3 -> ktwothird; + rhs = rhs /. -2/3 -> -ktwothird; + + rhs = rhs /. 4/3 -> kfourthird; + rhs = rhs /. -4/3 -> -kfourthird; + + rhs = rhs /. 8/3 -> keightthird; + rhs = rhs /. -8/3 -> -keightthird; + *) + + (* Remove parentheses *) + rhs = rhs //. Parenthesis[xx_] -> xx; + + (* Avoid rational numbers *) + rhs = rhs /. Rational[xx_,yy_] :> N[xx/yy, 30]; + + rhs = rhs //. IfThen[cond1_,xx1_,yy1_] + IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> + IfThen[cond1, Simplify[ xx1 + xx2], Simplify[ yy1 + yy2]]; + + rhs = rhs //. ff1_ IfThen[cond1_,xx1_,yy1_] + IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> + IfThen[cond1, Simplify[ff1 xx1 + xx2], Simplify[ff1 yy1 + yy2]]; + + rhs = rhs //. IfThen[cond1_,xx1_,yy1_] + ff2_ IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> + IfThen[cond1, Simplify[ xx1 + ff2 xx2], Simplify[ yy1 + ff2 yy2]]; + + rhs = rhs //. ff1_ IfThen[cond1_,xx1_,yy1_] + ff2_ IfThen[cond2_,xx2_,yy2_] /; cond1==cond2 :> + IfThen[cond1, Simplify[ff1 xx1 + ff2 xx2], Simplify[ff1 yy1 + ff2 yy2]]; + + (* Is this still a good idea when FMA instructions are used? *) + rhs = rhs //. xx_ yy_ + xx_ zz_ -> xx (yy+zz); + rhs = rhs //. xx_ yy_ - xx_ zz_ -> xx (yy-zz); + + rhs = rhs /. Power[E, power_] -> exp[power]; + + (* 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 /. Power[xx_, power_] -> pow[xx, power]; + + If[vectorise === True, + rhs = vectoriseExpression[rhs]]; + + (* Remove Scalar[] after vectorising *) + rhs = rhs /. Scalar[xx_] -> xx}, + (* else *) + rhs = rhs /. Power[xx_, power_] -> xx^power]; + (* Print[rhs//FullForm];*) + rhs]; +ErrorDefinition[ReplacePowers]; + +End[]; + +EndPackage[]; diff --git a/Tools/CodeGen/Differencing.m b/Tools/CodeGen/Differencing.m index dbfed8e..c0ebcf1 100644 --- a/Tools/CodeGen/Differencing.m +++ b/Tools/CodeGen/Differencing.m @@ -125,7 +125,7 @@ point. Should be checked by someone competent! *) -BeginPackage["Differencing`", {"CodeGen`", "Kranc`", "MapLookup`", +BeginPackage["Differencing`", {"CodeGen`", "CodeGenC`", "CodeGenCactus`", "Kranc`", "MapLookup`", (* "LinearAlgebra`MatrixManipulation`", *) "Errors`"}]; CreateDifferencingHeader::usage = ""; @@ -500,7 +500,7 @@ DifferenceGFTerm[op_, i_, j_, k_, vectorise_] := If[Cases[{remaining}, shift[_], Infinity] != {}, ThrowError["Could not parse difference operator", op]]; - If[CodeGen`SOURCELANGUAGE == "C", + If[CodeGenC`SOURCELANGUAGE == "C", If[vectorise, remaining "KRANC_GFOFFSET3D(u," <> diff --git a/Tools/CodeGen/Jacobian.m b/Tools/CodeGen/Jacobian.m index 20f24df..a7085e5 100644 --- a/Tools/CodeGen/Jacobian.m +++ b/Tools/CodeGen/Jacobian.m @@ -18,7 +18,8 @@ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA *) -BeginPackage["Jacobian`", {"Errors`", "Helpers`", "Kranc`", "Differencing`", "MapLookup`", "CodeGen`", "KrancGroups`"}]; +BeginPackage["Jacobian`", {"Errors`", "Helpers`", "Kranc`", "Differencing`", "MapLookup`", + "CodeGen`", "CodeGenC`", "KrancGroups`"}]; JacobianQ; InsertJacobian; diff --git a/Tools/CodeGen/KrancThorn.m b/Tools/CodeGen/KrancThorn.m index 3d1d708..26ca6ef 100644 --- a/Tools/CodeGen/KrancThorn.m +++ b/Tools/CodeGen/KrancThorn.m @@ -223,7 +223,7 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[ (* Add the predefinitions into the calcs *) calcs = Map[Join[#, {PreDefinitions -> pDefs}] &, calcs]; - ext = CodeGen`SOURCESUFFIX; + ext = CodeGenC`SOURCESUFFIX; (* Construct a source file for each calculation *) allParams = Join[Map[ParamName, realParamDefs], diff --git a/Tools/CodeGen/Schedule.m b/Tools/CodeGen/Schedule.m index 66dec5a..f596ad1 100644 --- a/Tools/CodeGen/Schedule.m +++ b/Tools/CodeGen/Schedule.m @@ -109,7 +109,7 @@ scheduleCalc[calc_, groups_] := SynchronizedGroups -> If[StringMatchQ[#, "*MoL_CalcRHS*", IgnoreCase -> True] || StringMatchQ[#, "*MoL_RHSBoundaries*", IgnoreCase -> True], {}, groupsToSync], - Language -> CodeGen`SOURCELANGUAGE, + Language -> CodeGenC`SOURCELANGUAGE, Comment -> lookup[calc, Name] }, If[triggered, {TriggerGroups -> lookup[calc, TriggerGroups]}, @@ -147,7 +147,7 @@ scheduleCalc[calc_, groups_] := fnSched = { Name -> lookup[calc, Name], SchedulePoint -> "in " <> groupName, - Language -> CodeGen`SOURCELANGUAGE, + Language -> CodeGenC`SOURCELANGUAGE, Comment -> lookup[calc, Name] }; @@ -165,7 +165,7 @@ scheduleCalc[calc_, groups_] := SchedulePoint -> "in " <> bcGroupName, SynchronizedGroups -> groupsToSync, Options -> "level", - Language -> CodeGen`SOURCELANGUAGE, + Language -> CodeGenC`SOURCELANGUAGE, Comment -> lookup[calc, Name] <> "_SelectBCs" }; diff --git a/Tools/CodeGen/Thorn.m b/Tools/CodeGen/Thorn.m index 20c41bc..afffb47 100644 --- a/Tools/CodeGen/Thorn.m +++ b/Tools/CodeGen/Thorn.m @@ -23,7 +23,7 @@ (* This package provides a set of functions to create the various parts of a Cactus thorn and assemble them. *) -BeginPackage["Thorn`", "CodeGen`", "CalculationFunction`", +BeginPackage["Thorn`", "CodeGen`", "CodeGenC`", "CodeGenCactus`", "CalculationFunction`", "CalculationBoundaries`", "MapLookup`", "KrancGroups`", "Helpers`", "Errors`", "Kranc`"]; @@ -502,11 +502,11 @@ CreateSetterSource[calcs_, debug_, include_, imp_, If[!MatchQ[include, _List], Throw["CreateSetterSource: Include should be a list but is in fact " <> ToString[include]]]; - {whoWhen[CodeGen`SOURCELANGUAGE], + {whoWhen[CodeGenC`SOURCELANGUAGE], - "#define KRANC_" <> ToUpperCase[CodeGen`SOURCELANGUAGE] <> "\n\n", + "#define KRANC_" <> ToUpperCase[CodeGenC`SOURCELANGUAGE] <> "\n\n", - If[CodeGen`SOURCELANGUAGE == "C", + If[CodeGenC`SOURCELANGUAGE == "C", {IncludeSystemFile["assert.h"], IncludeSystemFile["math.h"], IncludeSystemFile["stdio.h"], @@ -593,8 +593,8 @@ CreateSymmetriesRegistrationSource[thornName_, implementationName_, GFs_, reflec Print["Registering Symmetries for: ", GFs]; ]; - lang = CodeGen`SOURCELANGUAGE; - CodeGen`SOURCELANGUAGE = "C"; + lang = CodeGenC`SOURCELANGUAGE; + CodeGenC`SOURCELANGUAGE = "C"; spec = Table[{FullName -> implementationName <> "::" <> ToString@GFs[[j]], Sym -> If[reflectionSymmetries === False, @@ -617,7 +617,7 @@ CreateSymmetriesRegistrationSource[thornName_, implementationName_, GFs_, reflec ] }; - CodeGen`SOURCELANGUAGE = lang; + CodeGenC`SOURCELANGUAGE = lang; tmp ]; @@ -644,8 +644,8 @@ CreateMoLRegistrationSource[spec_, debug_] := Print[" Primitive Gridfunctions: ", lookup[spec, PrimitiveGFs] ]; ]; - lang = CodeGen`SOURCELANGUAGE; - CodeGen`SOURCELANGUAGE= "C"; + lang = CodeGenC`SOURCELANGUAGE; + CodeGenC`SOURCELANGUAGE= "C"; tmp = {whoWhen["C"], @@ -676,7 +676,7 @@ CreateMoLRegistrationSource[spec_, debug_] := lookup[spec, PrimitiveGFs]]], *) "return;\n"}]}; - CodeGen`SOURCELANGUAGE = lang; + CodeGenC`SOURCELANGUAGE = lang; tmp ]; @@ -909,8 +909,8 @@ CreateMoLBoundariesSource[spec_] := "\n}\n"}]; - lang = CodeGen`SOURCELANGUAGE; - CodeGen`SOURCELANGUAGE = "C"; + lang = CodeGenC`SOURCELANGUAGE; + CodeGenC`SOURCELANGUAGE = "C"; tmp = {whoWhen["C"], @@ -960,7 +960,7 @@ CreateMoLBoundariesSource[spec_] := "*/\n\n" }; - CodeGen`SOURCELANGUAGE = lang; + CodeGenC`SOURCELANGUAGE = lang; tmp ]; @@ -972,8 +972,8 @@ CreateMoLExcisionSource[spec_] := Print["Applying excision to GFs: ", gfs]; - currentlang = CodeGen`SOURCELANGUAGE; - CodeGen`SOURCELANGUAGE = "Fortran"; + currentlang = CodeGenC`SOURCELANGUAGE; + CodeGenC`SOURCELANGUAGE = "Fortran"; excisionExtrap[gf_] := " call ExcisionExtrapolate(ierr, " <> ToString@gf <> ", " <> ToString@gf @@ -1057,7 +1057,7 @@ CreateMoLExcisionSource[spec_] := "return\n"}] }; -CodeGen`SOURCELANGUAGE = currentlang; +CodeGenC`SOURCELANGUAGE = currentlang; body ]; @@ -1244,8 +1244,8 @@ CreateMPCharSource[spec_, debug_] := groups = Map[unqualifiedGroupName, lookup[spec, Groups]]; - lang = CodeGen`SOURCELANGUAGE; - CodeGen`SOURCELANGUAGE = "C"; + lang = CodeGenC`SOURCELANGUAGE; + CodeGenC`SOURCELANGUAGE = "C"; tmp = {whoWhen["C"], @@ -1307,7 +1307,7 @@ charInfoFunction["C2P", spec, debug], "\n\n", charInfoFunction["WHATEVER", spec, debug] }; -CodeGen`SOURCELANGUAGE = lang; +CodeGenC`SOURCELANGUAGE = lang; tmp ]; @@ -1334,8 +1334,8 @@ CreatePrecompMacros[functions_] := CreateStartupFile[thornName_, bannerText_] := Module[{tmp, lang}, - lang = CodeGen`SOURCELANGUAGE; - CodeGen`SOURCELANGUAGE = "C"; + lang = CodeGenC`SOURCELANGUAGE; + CodeGenC`SOURCELANGUAGE = "C"; tmp = {whoWhen["C"], @@ -1345,7 +1345,7 @@ CreateStartupFile[thornName_, bannerText_] := "CCTK_RegisterBanner(banner);\n", "return 0;\n"}]}; - CodeGen`SOURCELANGUAGE = lang; + CodeGenC`SOURCELANGUAGE = lang; tmp ]; |