From d8e9df66d2478ac7ad0b1ed3b24fc0aa048bffa8 Mon Sep 17 00:00:00 2001 From: Erik Schnetter Date: Mon, 21 May 2012 11:37:35 -0400 Subject: Update WaveToy example, generate OpenCL code as well --- m/WaveToy.m | 186 +++++++++++++++++++++++++++++++++++++++--------------------- 1 file changed, 122 insertions(+), 64 deletions(-) diff --git a/m/WaveToy.m b/m/WaveToy.m index de2c73e..7163269 100644 --- a/m/WaveToy.m +++ b/m/WaveToy.m @@ -9,35 +9,18 @@ SetSourceLanguage["C"]; derivOrder = 4; derivatives = -{ - PDstandardNth[i_] -> StandardCenteredDifferenceOperator[1,derivOrder/2,i], - PDstandardNth[i_, i_] -> StandardCenteredDifferenceOperator[2,derivOrder/2,i], - PDstandardNth[i_, j_] -> StandardCenteredDifferenceOperator[1,derivOrder/2,i] - StandardCenteredDifferenceOperator[1,derivOrder/2,j] -(* - PDstandardNth[i_, i_, i_] -> - StandardCenteredDifferenceOperator[3,derivOrder/2,i], - PDstandardNth[i_, i_, j_] -> - StandardCenteredDifferenceOperator[2,derivOrder/2,i] - StandardCenteredDifferenceOperator[1,derivOrder/2,j], - PDstandardNth[i_, j_, i_] -> - StandardCenteredDifferenceOperator[2,derivOrder/2,i] - StandardCenteredDifferenceOperator[1,derivOrder/2,j], - PDstandardNth[j_, i_, i_] -> - StandardCenteredDifferenceOperator[2,derivOrder/2,i] - StandardCenteredDifferenceOperator[1,derivOrder/2,j], - PDstandardNth[i_, j_, k_] -> - StandardCenteredDifferenceOperator[1,derivOrder/2,i] - StandardCenteredDifferenceOperator[1,derivOrder/2,j] - StandardCenteredDifferenceOperator[1,derivOrder/2,k] -*) -}; + {PDstandardNth[i_] -> + StandardCenteredDifferenceOperator[1,derivOrder/2,i], + PDstandardNth[i_, i_] -> + StandardCenteredDifferenceOperator[2,derivOrder/2,i], + PDstandardNth[i_, j_] -> + StandardCenteredDifferenceOperator[1,derivOrder/2,i] + StandardCenteredDifferenceOperator[1,derivOrder/2,j]}; FD = PDstandardNth; ResetJacobians; DefineJacobian[PD, FD, KD, Zero3]; -(* timelevels *) evolutionTimelevels = 2; KD = KroneckerDelta; @@ -47,19 +30,23 @@ KD = KroneckerDelta; (******************************************************************************) (* Register the tensor quantities with the TensorTools package *) -Map [DefineTensor, {u, rho, v, w}]; +Map[DefineTensor, {u, rho, eps}]; (******************************************************************************) (* Groups *) (******************************************************************************) evolvedGroups = - {SetGroupName [CreateGroupFromTensor [u ], "WT_u" ], - SetGroupName [CreateGroupFromTensor [rho], "WT_rho"]}; -evaluatedGroups = {}; + {SetGroupName[CreateGroupFromTensor[u ], "WT_u" ], + SetGroupName[CreateGroupFromTensor[rho], "WT_rho"]}; +evaluatedGroups = + {SetGroupName[CreateGroupFromTensor[eps], "WT_eps"]}; -declaredGroups = Join [evolvedGroups, evaluatedGroups]; -declaredGroupNames = Map [First, declaredGroups]; +declaredGroups = Join[evolvedGroups, evaluatedGroups]; +declaredGroupNames = Map[First, declaredGroups]; + +extraGroups = + {{"grid::coordinates", {x, y, z, r}}}; groups = declaredGroups; @@ -67,49 +54,120 @@ groups = declaredGroups; (* Initial data *) (******************************************************************************) -initialCalc = -{ - Name -> "WT_Gaussian", - Schedule -> {"AT initial"}, - (* Where -> Boundary, *) - (* Where -> Interior, *) - Equations -> - { - u -> 0, - rho -> 0 - } -}; +initialCalcGaussian = + {Name -> "WT_Gaussian", + Schedule -> {"AT initial"}, + ConditionalOnKeyword -> {"initial_data", "Gaussian"}, + Equations -> + {u -> amplitude Exp[-(1/2) (r/width)^2], + rho -> 0}}; + +initialCalcStanding = + {Name -> "WT_Standing", + Schedule -> {"AT initial"}, + ConditionalOnKeyword -> {"initial_data", "Standing"}, + Shorthands -> {k, omega}, + Equations -> + {k -> Pi / width, + omega -> Sqrt[3 k^2], + u -> amplitude Cos[k x] Cos[k y] Cos[k z] Cos[omega t], + rho -> amplitude Cos[k x] Cos[k y] Cos[k z] Sin[omega t] (-omega)}}; (******************************************************************************) (* Evolution equations *) (******************************************************************************) evolCalc = -{ - Name -> "WT_RHS", - Schedule -> {"IN MoL_CalcRHS", "AT analysis"}, - Where -> Interior, - Equations -> - { - dot[u] -> rho, - dot[rho] -> KD[ua,ub] PD[u,la,lb] - } -}; + {Name -> "WT_RHS", + Schedule -> {"IN MoL_CalcRHS" (*"AT analysis"*)}, + Where -> Interior, + Equations -> + {dot[u] -> rho, + dot[rho] -> KD[ua,ub] PD[u,la,lb]}}; + +(******************************************************************************) +(* Boundary conditions *) +(******************************************************************************) + +boundaryCalc = + {Name -> "WT_Dirichlet", + Schedule -> {"IN MoL_CalcRHS", "AT analysis"}, + Where -> Boundary, + Equations -> + {dot[u] -> 0, + dot[rho] -> 0}}; + +(******************************************************************************) +(* Energy *) +(******************************************************************************) + +energyCalc = + {Name -> "WT_Energy", + Schedule -> {"AT analysis"}, + Where -> Interior, + Equations -> + {eps -> 1/2 (rho^2 + KD[ua,ub] PD[u,la] PD[u,lb])}}; + +energyBoundaryCalc = + {Name -> "WT_EnergyBoundary", + Schedule -> {"AT analysis"}, + Where -> Boundary, + Equations -> + {eps -> 0}}; + +(******************************************************************************) +(* Parameters *) +(******************************************************************************) + +keywordParameters = + {{Name -> "initial_data", + AllowedValues -> {"Gaussian", "Standing"}, + Default -> "Gaussian"}}; + +realParameters = + {{Name -> amplitude, + Description -> "Amplitude of initial Gaussian", + Default -> 1}, + {Name -> width, + Description -> "Width of initial Gaussian", + Default -> 1}}; (******************************************************************************) (* Construct the thorns *) (******************************************************************************) calculations = -{ - initialCalc, - evolCalc -}; - -CreateKrancThornTT [groups, ".", "ML_WaveToy", - Calculations -> calculations, - DeclaredGroups -> declaredGroupNames, - PartialDerivatives -> derivatives, - UseLoopControl -> True, - EvolutionTimelevels -> evolutionTimelevels -]; + {initialCalcGaussian, initialCalcStanding, + evolCalc, boundaryCalc, + energyCalc, energyBoundaryCalc}; + +CreateKrancThornTT[ + groups, ".", "ML_WaveToy", + Calculations -> calculations, + DeclaredGroups -> declaredGroupNames, + PartialDerivatives -> derivatives, + UseLoopControl -> True, + UseVectors -> True, + EvolutionTimelevels -> evolutionTimelevels, + KeywordParameters -> keywordParameters, + RealParameters -> realParameters]; + +renameCalc[calc_, oldprefix_, newprefix_] := Module[ + {name, newname, newcalc}, + name = lookup[calc, Name]; + newname = StringReplace[name, + RegularExpression["^" <> oldprefix] -> newprefix]; + newcalc = mapReplace[calc, Name, newname]; + newcalc]; + +CreateKrancThornTT[ + groups, ".", "ML_WaveToy_CL", + Calculations -> Map[renameCalc[#, "WT_", "WT_CL_"]&, calculations], + DeclaredGroups -> declaredGroupNames, + PartialDerivatives -> derivatives, + UseLoopControl -> True, + UseOpenCL -> True, + UseVectors -> True, + EvolutionTimelevels -> evolutionTimelevels, + KeywordParameters -> keywordParameters, + RealParameters -> realParameters]; -- cgit v1.2.3