1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
|
Get["KrancThorn`"];
SetEnhancedTimes[False];
(**************************************************************************************)
(* Tensors *)
(**************************************************************************************)
(* Register the tensor quantities with the TensorTools package *)
Map[DefineTensor, {Den, S, En, rho, v, p}];
(**************************************************************************************)
(* Groups *)
(**************************************************************************************)
evolvedGroups = Map[CreateGroupFromTensor, {Den, S[uj], En}];
nonevolvedGroups = Map[CreateGroupFromTensor,
{
rho, v[uj], p
}];
declaredGroups = Join[evolvedGroups, nonevolvedGroups];
declaredGroupNames = Map[First, declaredGroups];
groups = declaredGroups;
(**************************************************************************************)
(* Initial data *)
(**************************************************************************************)
initialShockCalc =
{
Name -> "eulerauto_initial_shock",
Schedule -> {"at CCTK_INITIAL as eulerauto_initial"},
ConditionalOnKeyword -> {"initial_data", "shock"},
Equations ->
{
rho -> rhoR0 StepFunction[x-0.5] + rhoL0 (1-StepFunction[x-0.5]),
v[1] -> vR0 StepFunction[x-0.5] + vL0 (1-StepFunction[x-0.5]),
v[2] -> 0,
v[3] -> 0,
p -> pR0 StepFunction[x-0.5] + pL0 (1-StepFunction[x-0.5])
}
};
(**************************************************************************************)
(* Evolution equations *)
(**************************************************************************************)
(* Euler's equation is dot[u] + PD[F[ui],li] = 0
with
u = {D, S, E}
and
DF[ui] = rho v[ui]
SF[ui,uj] = rho v[ui] v[uj] + p Euc[ui,uj]
EnF[ui] = v[ui](En + p)
*)
eulerCons =
{
Name -> "eulerauto_cons_calc",
Primitives -> {rho, v[ui], p},
Equations ->
{
flux[Den,ui] -> rho v[ui],
flux[S[uj],ui] -> rho v[ui] v[uj] + p Euc[ui,uj],
flux[En,ui] -> v[ui](En + p)
},
ConservedEquations ->
{
Den -> rho,
S[ui] -> rho v[ui],
En -> p/(gamma-1) + 1/2 rho v[ui] v[uj] Euc[li,lj]
},
PrimitiveEquations ->
{
rho -> Den,
v[ui] -> S[ui] / Den,
p -> (gamma-1)(En - 1/2 Euc[li,lj] S[ui] S[uj]/Den)
}
}
(**************************************************************************************)
(* Parameters *)
(**************************************************************************************)
realParameters = {sigma, v0, amp, rhoR0, rhoL0, vR0, vL0, pR0, pL0, gamma};
keywordParameters = {
{
Name -> "initial_data",
Default -> "shock",
AllowedValues -> {"shock"}
}
};
(**************************************************************************************)
(* Construct the thorn *)
(**************************************************************************************)
calculations =
{
initialShockCalc
};
consCalculations = {eulerCons};
CreateKrancThornTT[groups, ".", "EulerAuto",
Calculations -> calculations,
ConservationCalculations -> consCalculations,
DeclaredGroups -> declaredGroupNames,
RealParameters -> realParameters,
KeywordParameters -> keywordParameters];
|
