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
|
/* 1920x1080 - HD screen */
linel: 230;
/* 3840x1600 - Ultrawide screen */
linel: 470;
/* 1366x768 - Poor man's thinkpad */
linel: 150;
derivabbrev: true;
/** **/
load("ctensor");
dim: 2;
cframe_flag: false;
ct_coords: [θ, ϕ];
lg : zeromatrix(2, 2);
/* dθ dθ */
lg[1, 1] : 1;
/* dϕ dϕ */
lg[2, 2] : sin(θ)^2;
/* Load metric and compute inverse */
cmetric();
/* mcs[i, j, k]: Christoffel symbols of the second kind (one upper index - k) */
christof(true);
ricci(true);
/** **/
load("itensor");
imetric(g);
idim(2);
metricconvert: true;
/** **/
depends([c, M, N, P], [u, θ]);
Cdown: zeromatrix(2, 2);
Cdown[1, 1] : 2*c;
Cdown[2, 2] : -2*c*sin(θ)^2;
Ndown: zeromatrix(2, 2);
Ndown[1, 1] : 2*diff(c, u);
Ndown[2, 2] : -2*diff(c, u)*sin(θ)^2;
f(m, n, T) := sum(ug[m, a] * sum(ug[n, b] * T[a, b], b, 1, 2), a, 1, 2);
Nup: zeromatrix(2, 2);
for i: 1 thru 2 do for j: 1 thru 2 do Nup[i, j]: f(i, j, Ndown);
Cup: zeromatrix(2, 2);
for i: 1 thru 2 do for j: 1 thru 2 do Cup[i, j]: f(i, j, Cdown);
/** Mass-loss formula **/
mass_aspect_eqn: Mu = (1/4)*covdiff(Nu([],[A,B]), B, A)
-(1/8)*Nd([A,B], [])*Nu([], [A,B]);
Mu: block([Cu,Cd,Nu,Nd,expr],
expr: ic_convert(mass_aspect_eqn),
Cu: Cup, Cd: Cdown, Nu: Nup, Nd: Ndown,
expand(ev(expr)));
/** Momentum-loss formula **/
momentum_aspect_eqn: Pu([A], []) = covdiff(M([],[]), A)
+ (1/8)*covdiff(Cu([],[B,C])*Nd([C,B],[]), A)
+ (1/4)*covdiff(covdiff(Cu([], [B, C]), B, A)
- g([],[C,E])*g([],[B,D])*covdiff(Cd([A,B],[]), D, E), C)
+ (1/4)*covdiff(Nu([], [B,C])*Cd([A,C],[]) - Cu([],[B,C])*Nd([A,C],[]), B)
- (1/4)*Nu([], [B,C])*covdiff(Cd([B,C],[]), A);
block([Cu,Cd,Nu,Nd,expr],
expr: ic_convert(eqn),
Cu: Cup, Cd: Cdown, Nu: Nup, Nd: Ndown,
ev(expr),
Pu[1]),expand;
sachs_coeff: -(2/3)*g([],[A,B])*(P([B], []) + Cd([B,C],[])*Ubu([],[C]) + βb([],[],B));
bondi_coeff: 2*N + 3*c*diff(c,θ) + 4*c^2*cot(θ);
Nu: block([expr, X, Cu, Cd, Nu, Nd, Ubu, βb, N],
expr: ic_convert(sachs_coeff=X([], [A])),
expr: subst([Cu=Cup, Cd=Cdown, Nu=Nup, Nd=Ndown], expr),
Ubu[1]: -(diff(c, θ) + 2*c*cot(θ)),
βb: -(1/4)*c^2,
expand(ev(expr)),
N: rhs(first(solve(X[1]=bondi_coeff, N))),
subst(diff(P[1],u)=Pu[1], diff(-3*N, u)));
|