diff options
| -rw-r--r-- | Auxiliary/Cactus/KrancNumericalTools/GenericFD/param.ccl | 12 | ||||
| -rw-r--r-- | Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.c | 36 | ||||
| -rw-r--r-- | Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.h | 2 | ||||
| -rw-r--r-- | Doc/KrancDoc.tex | 97 | ||||
| -rw-r--r-- | Tools/CodeGen/CalculationFunction.m | 14 | ||||
| -rw-r--r-- | Tools/CodeGen/Jacobian.m | 124 | ||||
| -rw-r--r-- | Tools/CodeGen/Kranc.m | 3 | ||||
| -rw-r--r-- | Tools/CodeGen/KrancThorn.m | 2 | ||||
| -rw-r--r-- | Tools/CodeGen/Param.m | 7 |
9 files changed, 288 insertions, 9 deletions
diff --git a/Auxiliary/Cactus/KrancNumericalTools/GenericFD/param.ccl b/Auxiliary/Cactus/KrancNumericalTools/GenericFD/param.ccl index 91075c9..85de1f9 100644 --- a/Auxiliary/Cactus/KrancNumericalTools/GenericFD/param.ccl +++ b/Auxiliary/Cactus/KrancNumericalTools/GenericFD/param.ccl @@ -30,3 +30,15 @@ CCTK_INT boundary_width "width of boundary (fix later to use Cactus boundary ca { -1:* :: "Any integer" } 1 + +restricted: +CCTK_STRING jacobian_group "Name of group containing Jacobian" STEERABLE=RECOVER +{ + "" :: "String of the form <implementation>::<groupname>" +} "" + +restricted: +CCTK_STRING jacobian_derivative_group "Name of group containing Jacobian derivative" STEERABLE=RECOVER +{ + "" :: "String of the form <implementation>::<groupname>" +} "" diff --git a/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.c b/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.c index f23f8dc..1663231 100644 --- a/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.c +++ b/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.c @@ -508,3 +508,39 @@ void GenericFD_AssertGroupStorage(cGH const * restrict const cctkGH, const char } } } + +/* Return a list of pointers to the members of a named group */ +void GenericFD_GroupDataPointers(cGH const * restrict const cctkGH, const char *group_name, + int nvars, CCTK_REAL **ptrs) +{ + int group_index, status; + cGroup group_info; + + group_index = CCTK_GroupIndex(group_name); + if (group_index < 0) + CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING, + "Error return %d trying to get group index for group \'%s\'", + group_index, + group_name); + + status = CCTK_GroupData(group_index, &group_info); + if (status < 0) + CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING, + "Error return %d trying to get info for group \'%s\'", + status, + group_name); + + if (group_info.numvars != nvars) + { + CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING, + "Group \'%s\' has %d variables but %d were expected", + group_name, group_info.numvars, nvars); + } + + int v1 = CCTK_FirstVarIndex(group_name); + + for (int v = 0; v < nvars; v++) + { + ptrs[v] = (CCTK_REAL *) CCTK_VarDataPtrI(cctkGH, 0 /* timelevel */, v1+v); + } +} diff --git a/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.h b/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.h index b7a4239..5974e63 100644 --- a/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.h +++ b/Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/GenericFD.h @@ -152,6 +152,8 @@ void GenericFD_LoopOverBoundary(cGH const * restrict cctkGH, Kranc_Calculation c void GenericFD_LoopOverBoundaryWithGhosts(cGH const * restrict cctkGH, Kranc_Calculation calc); void GenericFD_LoopOverInterior(cGH const * restrict cctkGH, Kranc_Calculation calc); +void GenericFD_GroupDataPointers(cGH const * restrict const cctkGH, const char *group_name, + int nvars, CCTK_REAL **ptrs); #ifdef __cplusplus diff --git a/Doc/KrancDoc.tex b/Doc/KrancDoc.tex index e357de5..131f5c2 100644 --- a/Doc/KrancDoc.tex +++ b/Doc/KrancDoc.tex @@ -3,6 +3,7 @@ \usepackage{tabularx} \usepackage{graphicx} \usepackage{alltt} +\usepackage{hyperref} \addtolength{\oddsidemargin}{-0.25in} \addtolength{\textwidth}{1in} @@ -243,6 +244,7 @@ by the Kranc system, such as a name for the calculation. %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \chapter{Using Kranc} +\label{chp:usingkranc} %% \section{Types of arguments} @@ -925,6 +927,101 @@ same arguments, but they can be tensorial in nature. \end{tabularx} \end{center} +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +\chapter{Additional Features} + +In addition to the basic functionality described in Chapter +\ref{chp:usingkranc}, Kranc provides a number of additional features. +Each of these features can be used independently, and are typically +enabled by an option of the form \verb|Use* -> True| in the call to +\verb|CreateThorn|. + +\section{Jacobians} + +Kranc allows the user to relate finite differencing operators to +partial derivatives via an arbitrary user-defined Jacobian +transformation provided in a grid function. This feature is enabled +by setting \verb|UseJacobian -> True| in the options to \verb|CreateThorn| or +\verb|CreateThornTT|. + +Kranc does not provide the Jacobian grid function; it might be +provided by an external infrastructure (for example for multi-block +schemes), or could be provided easily in the user's thorn, or another +Kranc-generated thorn. Wherever the Jacobian is provided, it must +adhere to the following conventions. There should be one Cactus group +for the components of the Jacobian matrix and another for the +components of its first spatial derivative (this is necessary for +systems containing second spatial derivatives). The components should +be real-valued grid functions declared in a similar manner to the +following: + +\begin{verbatim} +CCTK_REAL jac type=GF timelevels=1 +{ + J11, J12, J13, J21, J22, J23, J31, J32, J33 +} "Jacobian of the coordinate transformation" + +CCTK_REALd djac type=GF timelevels=1 +{ + dJ111, dJ112, dJ113, dJ122, dJ123, dJ133, + dJ211, dJ212, dJ213, dJ222, dJ223, dJ233, + dJ311, dJ312, dJ313, dJ322, dJ323, dJ333, +} "Derivative of the Jacobian" +\end{verbatim} + +The names of the groups and variables are not important, but the order +of the variables within the groups is critical. + +The GenericFD thorn provides two parameters, \verb|jacobian_group| and +\verb|jacobian_derivative_group| which should be set by the user in their +parameter file to the names of the Jacobian and Jacobian derivative +groups. With the above Jacobian definition, provided by a thorn with +implementation \verb|MyCoordTransform| the user would set + +\begin{verbatim} +GenericFD::jacobian = "MyCoordTransform::jac" +GenericFD::jacobian_derivative = "MyCoordTransform::djac" +\end{verbatim} + +The partial derivatives, associated with certain finite difference +operators, specified in the user's calculation, will then be +multiplied by the Jacobian. If the user specifies the following, + +\begin{center} +\begin{minipage}{0.8 \textwidth} +\begin{verbatim} +derivs = { + PDstandard2nd[i_] -> DZero[i]} + +... + +dot[v] -> PDstandard2nd[v,1] +\end{verbatim} +\end{minipage} +\end{center} + +the code that will actually be generated will be + +\begin{center} +\begin{minipage}{0.8 \textwidth} +\begin{verbatim} +dot[v] -> J11 PDstandard2nd[v,1] + J21 PDstandard2nd[v,2] + + J31 PDstandard2nd[v,3] +\end{verbatim} +\end{minipage} +\end{center} + +Note: + +\begin{itemize} +\item The Jacobian multiplication introduces an additional performance + cost to the simulation, so it should not be enabled unless + necessary. +\item It is currently not possible to generate a thorn which applies + the Jacobians optionally based on a run-time parameter or depending + on the current Carpet map. This is planned for the future. +\end{itemize} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% diff --git a/Tools/CodeGen/CalculationFunction.m b/Tools/CodeGen/CalculationFunction.m index 99d4343..1ca9513 100644 --- a/Tools/CodeGen/CalculationFunction.m +++ b/Tools/CodeGen/CalculationFunction.m @@ -20,7 +20,7 @@ BeginPackage["CalculationFunction`", {"CodeGen`", "MapLookup`", "KrancGroups`", "Differencing`", "Errors`", - "Helpers`", "Kranc`", "Optimize`"}]; + "Helpers`", "Kranc`", "Optimize`", "Jacobian`"}]; CreateCalculationFunction::usage = ""; VerifyCalculation::usage = ""; @@ -344,11 +344,13 @@ pdCanonicalOrdering[name_[inds___] -> x_] := Options[CreateCalculationFunction] = ThornOptions; -CreateCalculationFunction[calc_, debug_, imp_, opts:OptionsPattern[]] := +CreateCalculationFunction[calcp_, debug_, imp_, opts:OptionsPattern[]] := Module[{gfs, allSymbols, knownSymbols, shorts, eqs, parameters, parameterRules, functionName, dsUsed, groups, pddefs, cleancalc, eqLoop, where, - addToStencilWidth, pDefs, haveCondTextuals, condTextuals}, + addToStencilWidth, pDefs, haveCondTextuals, condTextuals, calc}, + + calc = If[OptionValue[UseJacobian], InsertJacobian[calcp, opts], calcp]; cleancalc = removeUnusedShorthands[calc]; If[OptionValue[CSE], @@ -359,6 +361,7 @@ CreateCalculationFunction[calc_, debug_, imp_, opts:OptionsPattern[]] := eqs = lookup[cleancalc, Equations]; parameters = lookupDefault[cleancalc, Parameters, {}]; groups = lookup[cleancalc, Groups]; + If[OptionValue[UseJacobian], groups = Join[groups, JacobianGroups[]]]; pddefs = lookupDefault[cleancalc, PartialDerivatives, {}]; where = lookupDefault[cleancalc, Where, Everywhere]; addToStencilWidth = lookupDefault[cleancalc, AddToStencilWidth, 0]; @@ -420,7 +423,8 @@ CreateCalculationFunction[calc_, debug_, imp_, opts:OptionsPattern[]] := allSymbols = calculationSymbols[cleancalc]; knownSymbols = Join[lookupDefault[cleancalc, AllowedSymbols, {}], gfs, shorts, parameters, {dx,dy,dz,idx,idy,idz,t, Pi, E, Symbol["i"], Symbol["j"], Symbol["k"], normal1, normal2, - normal3, tangentA1, tangentA2, tangentA3, tangentB1, tangentB2, tangentB3}]; + normal3, tangentA1, tangentA2, tangentA3, tangentB1, tangentB2, tangentB3}, + If[OptionValue[UseJacobian], JacobianSymbols[], {}]]; unknownSymbols = Complement[allSymbols, knownSymbols]; @@ -452,6 +456,8 @@ CreateCalculationFunction[calc_, debug_, imp_, opts:OptionsPattern[]] := InitialiseFDVariables[OptionValue[UseVectors]], definePreDefinitions[pDefs], + If[OptionValue[UseJacobian], CreateJacobianVariables[], {}], + If[Cases[{pddefs}, SBPDerivative[_], Infinity] != {}, CommentedBlock["Compute Summation By Parts derivatives", IncludeFile["sbp_calc_coeffs.h"]], diff --git a/Tools/CodeGen/Jacobian.m b/Tools/CodeGen/Jacobian.m new file mode 100644 index 0000000..59a57f9 --- /dev/null +++ b/Tools/CodeGen/Jacobian.m @@ -0,0 +1,124 @@ + +(* Copyright 2011 Ian Hinder + + This file is part of Kranc. + + Kranc is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + Kranc is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with Kranc; if not, write to the Free Software + Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +*) + +BeginPackage["Jacobian`", {"Errors`", "Helpers`", "Kranc`", "Differencing`", "MapLookup`", "CodeGen`"}]; + +JacobianQ; +InsertJacobian; +CreateJacobianVariables; +JacobianGenericFDParameters; +JacobianSymbols; +JacobianGroups; + +Begin["`Private`"]; + +Options[JacobianQ] = ThornOptions; +JacobianQ[opts:OptionsPattern[]] := + Length[OptionValue[Jacobian]] > 0; + +(* Assign a shorthand containing the Jacobian multiplied by the passed + 1st derivative *) +jacobianShorthand[d:(deriv_[var_, i_])] := + Module[{}, + derivToJacDeriv[d] -> + Sum[Symbol["J"<>ToString[j]<>ToString[i]] deriv[var, j], {j, 1 3}] + ]; + +(* Assign a shorthand containing the Jacobian multiplied by the passed + 2nd derivative *) +jacobianShorthand[d:(deriv_[var_, i_,j_])] := + Module[{ip,jp}, + {ip,jp} = Sort[{i,j}]; (* dJ is symmetric in the last two indices *) + derivToJacDeriv[d] -> + Sum[Symbol["dJ"<>ToString[a]<>ToString[ip]<>ToString[jp]] deriv[var, a], {a, 1 3}] + + Sum[Symbol["J"<>ToString[a]<>ToString[i]] Symbol["J"<>ToString[b]<>ToString[j]] deriv[var, a, b], {a, 1 3}, {b, 1, 3}] + ]; + +(* Convert a 1st derivative to a Jacobian-multiplied derivative *) +derivToJacDeriv[deriv_[var_, i_]] := + Symbol["Global`Jac"<>ToString[deriv]<>ToString[i]<>ToString[var]]; + +(* Convert a 2nd derivative to a Jacobian-multiplied derivative *) +derivToJacDeriv[deriv_[var_, i_, j_]] := + Symbol["Global`Jac"<>ToString[deriv]<>ToString[i]<>ToString[j]<>ToString[var]]; + +(* Given a calculation containing partial derivatives, return a + version of the calculation with all the partial derivatives multiplied + by the Jacobian *) +Options[InsertJacobian] = ThornOptions; +InsertJacobian[calc_List, opts:OptionsPattern[]] := + Module[{pdDefs, derivs, newShortDefs, newShorts, combinedShorts, combinedEqs, combinedCalc, eqs, newEqs}, + pdDefs = OptionValue[PartialDerivatives]; + derivs = GridFunctionDerivativesInExpression[pdDefs, lookup[calc, Equations]]; + If[Length[derivs] === 0, Return[calc]]; + newShortDefs = Map[jacobianShorthand, derivs]; + newShorts = Map[First, newShortDefs]; + combinedShorts = Join[lookupDefault[calc, Shorthands, {}], newShorts]; + eqs = lookup[calc, Equations]; + newEqs = eqs /. (x_?(MemberQ[derivs, #] &) :> derivToJacDeriv[x]); + combinedEqs = Join[newShortDefs, newEqs]; + combinedCalc = mapReplace[mapReplace[mapEnsureKey[calc, Shorthands, {}], Shorthands, combinedShorts], Equations, combinedEqs]; + combinedCalc]; + +(* Define local pointers to the members of the Jacobian and Jacobian + derivatives groups *) +CreateJacobianVariables[] := +CommentedBlock["Jacobian variable pointers", + {"if (strlen(jacobian_group) == 0 || strlen(jacobian_derivative_group) == 0)\n", + "{\n", + " CCTK_WARN (1, \"GenericFD::jacobian_group and GenericFD::jacobian_derivative_group must both be set to valid group names\");\n", + "}\n\n", + "CCTK_REAL *jacobian_ptrs[9];\n", + "GenericFD_GroupDataPointers(cctkGH, jacobian_group,\n", + " 9, jacobian_ptrs);\n", + "\n", + Table[{"CCTK_REAL *J",i,j," = jacobian_ptrs[",(i-1)*3+j-1,"];\n"},{i,1,3},{j,1,3}], + "\n", + "CCTK_REAL *jacobian_derivative_ptrs[18];\n", + "GenericFD_GroupDataPointers(cctkGH, jacobian_derivative_group,\n", + " 18, jacobian_derivative_ptrs);\n", + "\n", + Module[{syms = Flatten[Table[{"dJ",i,j,k},{i,1,3},{j,1,3},{k,j,3}],2]}, + MapIndexed[{"CCTK_REAL *", #1, " = jacobian_derivative_ptrs[", #2-1, "];\n"} &, syms]]}]; + +(* List of symbols which should be allowed in a calculation *) +JacobianSymbols[] := + Map[Symbol, Join[Flatten[Table[FlattenBlock[{"dJ",i,j,k}],{i,1,3},{j,1,3},{k,j,3}],2], + Flatten[Table[FlattenBlock[{"J",i,j}],{i,1,3},{j,1,3}],1]]]; + +(* Parameters to inherit from GenericFD *) +JacobianGenericFDParameters[] := + {{Name -> "jacobian_group", Type -> "CCTK_STRING"}, + {Name -> "jacobian_derivative_group", Type -> "CCTK_STRING"}}; + +(* The symbols which are used for the Jacobian variables in the + generated source code. These do not have to coincide with the + actual variable names, as the variable pointers are read using + CCTK_VarDataPtr. *) +JacobianGroups[] := + {{"unknown::unknown", {Global`J11, Global`J12, Global`J13, Global`J21, Global`J22, Global`J23, Global`J31, Global`J32, Global`J33}}, + {"unknown::unknown", {Global`dJ111, Global`dJ112, Global`dJ113, Global`dJ122, Global`dJ123, Global`dJ133, + Global`dJ211, Global`dJ212, Global`dJ213, Global`dJ222, Global`dJ223, Global`dJ233, + Global`dJ311, Global`dJ312, Global`dJ313, Global`dJ322, Global`dJ323, Global`dJ333}}}; + + +End[]; + +EndPackage[]; diff --git a/Tools/CodeGen/Kranc.m b/Tools/CodeGen/Kranc.m index bd2c9c7..1077305 100644 --- a/Tools/CodeGen/Kranc.m +++ b/Tools/CodeGen/Kranc.m @@ -73,7 +73,8 @@ ThornOptions = UseVectors -> False, ProhibitAssignmentToGridFunctionsRead -> False, IncludeFiles -> {}, - CSE -> False}; + CSE -> False, + UseJacobian -> False}; (* Thorn.m *) diff --git a/Tools/CodeGen/KrancThorn.m b/Tools/CodeGen/KrancThorn.m index 6463345..bb1b115 100644 --- a/Tools/CodeGen/KrancThorn.m +++ b/Tools/CodeGen/KrancThorn.m @@ -179,7 +179,7 @@ CreateKrancThorn[groupsOrig_, parentDirectory_, thornName_, opts:OptionsPattern[ inheritedRealParams, inheritedIntParams, inheritedKeywordParams, extendedRealParams, extendedIntParams, extendedKeywordParams, evolutionTimelevels, defaultEvolutionTimelevels, - calcs]; + calcs, opts]; (* Construct the schedule file *) InfoMessage[Terse, "Creating schedule file"]; diff --git a/Tools/CodeGen/Param.m b/Tools/CodeGen/Param.m index 90cbb7d..1b4ad37 100644 --- a/Tools/CodeGen/Param.m +++ b/Tools/CodeGen/Param.m @@ -18,7 +18,7 @@ Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA *) -BeginPackage["Param`", {"Thorn`", "Errors`", "Helpers`", "MapLookup`", "KrancGroups`", "Kranc`"}]; +BeginPackage["Param`", {"Thorn`", "Errors`", "Helpers`", "MapLookup`", "KrancGroups`", "Kranc`", "Jacobian`"}]; CreateKrancParam; MakeFullParamDefs; @@ -118,12 +118,13 @@ extendParameters[imp_, reals_, ints_, keywords_] := Return[{Name -> imp, ExtendedParameters -> Join[realStructs, intStructs, keywordStructs]}], Return[{}]]]; +Options[CreateKrancParam] = ThornOptions; CreateKrancParam[evolvedGroups_, nonevolvedGroups_, groups_, thornName_, reals_, ints_, keywords_, inheritedReals_, inheritedInts_, inheritedKeywords_, extendedReals_, extendedInts_, extendedKeywords_, evolutionTimelevels_, defaultEvolutionTimelevels_, - calcs_] := + calcs_, opts:OptionsPattern[]] := Module[{nEvolved, evolvedMoLParam, evolvedGFs, (*constrainedMoLParam,*) genericfdStruct, realStructs, intStructs, allInherited, allExtended, implementationNames, molImplementation, @@ -194,7 +195,7 @@ CreateKrancParam[evolvedGroups_, nonevolvedGroups_, groups_, thornName_, { Name -> "GenericFD", UsedParameters -> - {} + If[OptionValue[UseJacobian], JacobianGenericFDParameters[], {}] }; realStructs = Map[krancParamStruct[#, "CCTK_REAL", False] &, reals]; |