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
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
|
# This is an until-end-of-line comment
THORN BSSN
DERIVATIVE
PDstandardNth[i_] -> StandardCenteredDifferenceOperator[1,derivOrder/2,i]
DERIVATIVE
PDstandardNth[i_,i_] -> StandardCenteredDifferenceOperator[2,derivOrder/2,i]
DERIVATIVE
PDstandardNth[i_,j_] -> StandardCenteredDifferenceOperator[1,derivOrder/2,i] *
StandardCenteredDifferenceOperator[1,derivOrder/2,j]
DERIVATIVE
PDdissipationNth[i_] ->
spacing[i]^(derivOrder+1) / 2^(derivOrder+2) *
StandardCenteredDifferenceOperator[derivOrder+2,derivOrder/2+1,i]
JACOBIAN {PD, FD, J, dJ}
TENSOR normal, tangentA, tangentB, dir
TENSOR xx, rr, th, ph
TENSOR J, dJ
TENSOR admg, admK, admalpha, admdtalpha, admbeta, admdtbeta, H, M
TENSOR g, detg, gu, G, R, trR, Km, trK, cdphi, cdphi2
TENSOR phi, gt, At, Xt, Xtn, alpha, A, beta, B, Atm, Atu, trA, Ats, trAts
TENSOR dottrK, dotXt
TENSOR cXt, cS, cA
TENSOR e4phi, em4phi, ddetg, detgt, gtu, ddetgt, dgtu, ddgtu, Gtl, Gtlu, Gt
TENSOR Rt, Rphi, gK
TENSOR T00, T0, T, rho, S
TENSOR x, y, z, r
TENSOR epsdiss
SYMMETRIC admg[la,lb], admK[la,lb]
SYMMETRIC g[la,lb], K[la,lb], R[la,lb], cdphi2[la,lb]
SYMMETRIC gt[la,lb], At[la,lb], Ats[la,lb], Rt[la,lb], Rphi[la,lb], T[la,lb]
SYMMETRIC {dJ[ua,lb,lc], lb, lc}
SYMMETRIC {G[ua,lb,lc], lb, lc}
SYMMETRIC {Gtl[la,lb,lc], lb, lc}
SYMMETRIC {Gt[ua,lb,lc], lb, lc}
SYMMETRIC {gK[la,lb,lc], la, lb}
SYMMETRIC gu[ua,ub], gtu[ua,ub], Atu[ua,ub]
SYMMETRIC {dgtu[ua,ub,lc], ua, ub}
SYMMETRIC {ddgtu[ua,ub,lc,ld], ua, ub}
SYMMETRIC {ddgtu[ua,ub,lc,ld], lc, ld}
CONNECTION {CD, PD, G}
CONNECTION {CDt, PD, Gt}
GROUP {phi , "log_confac"}
GROUP {gt[la,lb], "metric" }
GROUP {Xt[ua ], "Gamma" }
EXTRA_GROUP {"Grid::coordinates", {x, y, z, r}}
EXTRA_GROUP {"ADMBase::metric", {gxx, gxy, gxz, gyy, gyz, gzz}}
EXTRA_GROUP {"ADMBase::curv", {kxx, kxy, kxz, kyy, kyz, kzz}}
DEFINE pi = N[Pi,40]
DEFINE
detgExpr = Det [MatrixOfComponents [g [la,lb]]]
DEFINE
ddetgExpr[la_] =
Sum [D[Det[MatrixOfComponents[g[la, lb]]], X] PD[X, la],
{X, Union[Flatten[MatrixOfComponents[g[la, lb]]]]}]
DEFINE
detgtExpr = Det [MatrixOfComponents [gt[la,lb]]]
DEFINE
ddetgtExpr[la_] =
Sum [D[Det[MatrixOfComponents[gt[la, lb]]], X] PD[X, la],
{X, Union[Flatten[MatrixOfComponents[gt[la, lb]]]]}]
CALCULATION "Minkowski"
Schedule: {"IN ADMBase_InitialData"}
ConditionalOnKeyword: {"my_initial_data", "Minkowski"}
BEGIN EQUATIONS
phi -> IfThen[conformalMethod, 1, 0]
phi -> conformalMethod ? 1 : 0
gt[la,lb] -> KD[la,lb]
trK -> 0
At[la,lb] -> 0
Xt[ua] -> 0
alpha -> 1
A -> 0
beta[ua] -> 0
B[ua] -> 0
END EQUATIONS
END CALCULATION
CALCULATION "convertFromADMBase"
Schedule: {"AT initial AFTER ADMBase_PostInitial"}
ConditionalOnKeyword: {"my_initial_data", "ADMBase"}
SHORTHAND g[la,lb], detg, gu[ua,ub], em4phi
BEGIN EQUATIONS
g[la,lb] -> admg[la,lb]
detg -> detgExpr
gu[ua,ub] -> 1/detg detgExpr MatrixInverse [g[ua,ub]]
phi -> IfThen [conformalMethod, detg^(-1/6), Log[detg]/12]
em4phi -> IfThen [conformalMethod, phi^2, Exp[-4 phi]]
gt[la,lb] -> em4phi g[la,lb]
trK -> gu[ua,ub] admK[la,lb]
At[la,lb] -> em4phi (admK[la,lb] - (1/3) g[la,lb] trK)
alpha -> admalpha
beta[ua] -> admbeta[ua]
END EQUATIONS
END CALCULATION
CALCULATION convertFromADMBaseGammaCalc
Name: BSSN <> "_convertFromADMBaseGamma"
Schedule: {"AT initial AFTER " <> BSSN <> "_convertFromADMBase"}
ConditionalOnKeyword: {"my_initial_data", "ADMBase"}
# Do not synchronise right after this routine; instead, synchronise
# after extrapolating
Where: Interior
# Synchronise after this routine, so that the refinement boundaries
# are set correctly before extrapolating. (We will need to
# synchronise again after extrapolating because extrapolation does
# not fill ghost zones, but this is irrelevant here.)
SHORTHAND dir[ua]
SHORTHAND detgt, gtu[ua,ub], Gt[ua,lb,lc], theta
BEGIN EQUATIONS
dir[ua] -> Sign[beta[ua]]
detgt -> 1 (* detgtExpr *)
gtu[ua,ub] -> 1/detgt detgtExpr MatrixInverse [gt[ua,ub]]
Gt[ua,lb,lc] -> 1/2 gtu[ua,ud]
(PD[gt[lb,ld],lc] + PD[gt[lc,ld],lb] - PD[gt[lb,lc],ld])
Xt[ua] -> gtu[ub,uc] Gt[ua,lb,lc]
# If LapseACoeff=0, then A is not evolved, in the sense that it
# does not influence the time evolution of other variables.
A -> IfThen [LapseACoeff != 0
1 / (- harmonicF alpha^harmonicN)
(+ admdtalpha
- LapseAdvectionCoeff beta[ua] PDua[alpha,la]
- LapseAdvectionCoeff Abs[beta[ua]] PDus[alpha,la])
0]
theta -> thetaExpr
# If ShiftBCoeff=0 or theta ShiftGammaCoeff=0, then B^i is not
# evolved, in the sense that it does not influence the time
# evolution of other variables.
B[ua] -> IfThen [ShiftGammaCoeff ShiftBCoeff != 0
1 / (theta ShiftGammaCoeff)
(+ admdtbeta[ua]
- ShiftAdvectionCoeff beta[ub] PDua[beta[ua],lb]
- ShiftAdvectionCoeff Abs[beta[ub]] PDus[beta[ua],lb])
0]
END EQUATIONS
END CALCULATION
INHERITED_IMPLEMENTATION ADMBase, TmunuBase
KEYWORD_PARAMETER "my_initial_data"
# Visibility: "restricted"
# Description: "ddd"
AllowedValues: {"ADMBase", "Minkowski"}
Default: "ADMBase"
END KEYWORD_PARAMETER
REAL_PARAMETER LapseACoeff
Description: "Whether to evolve A in time"
Default: 0
END REAL_PARAMETER
END THORN
|