aboutsummaryrefslogtreecommitdiff
path: root/src/gr/setup_gr_gfas.maple
blob: ff3968192384c442147c8db29c05d494b3e804bf (plain)
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
# Maple code to set up all gridfn arrays
# $Header$

#
# setup_gr_gfas - setup all the GR gfas
# `Diff/gridfn2` - Diff() simplification based on gridfn properties known here
#

################################################################################

setup_gr_gfas :=
proc()
global
  @include "../maple/coords.minc",
  @include "../maple/gfa.minc",
  @include "../gr/gr_gfas.minc";

# metric and extrinsic curvature
make_gfa('g_dd', {inert}, [1..N, 1..N], symmetric);
make_gfa('K_dd', {inert}, [1..N, 1..N], symmetric);

# index-raised and contracted metric and extrinsic curvature
make_gfa('g_uu', {inert,fnd}, [1..N, 1..N], symmetric);
make_gfa('K_uu', {inert,fnd}, [1..N, 1..N], symmetric);
make_gfa('K', {inert,fnd}, [], none);

# xyz partial derivatives of metric
# ... as far as Maple is concerned, these are indeed separate gridfns;
#     the derivatives are actually taken by the interpolator
make_gfa('partial_d_g_dd', {inert}, [1..N, 1..N, 1..N], symmetric3_23);

# xyz partial derivatives of metric determinant
make_gfa('partial_d_ln_sqrt_g', {inert,fnd}, [1..N], none);

# xyz partial derivatives of inverse metric
# (computed in terms of xyz derivatives of g_dd)
make_gfa('partial_d_g_uu', {inert,fnd}, [1..N, 1..N, 1..N], symmetric3_23);

# radius of horizon
make_gfa('h', {inert,fnd}, [], none);

# outward-pointing *non*-unit normal (covector) to horizon
# and it's xyz-coordinate partial derivatives
make_gfa('s_d', {inert,fnd}, [1..N], none);
make_gfa('partial_d_s_d', {inert,fnd}, [1..N, 1..N], none);

# subexpressions for computing LHS of apparent horizon equation
# ... these are defined by (15) in my 1996 apparent horizon finding paper
make_gfa('HA', {inert,fnd}, [], none);
make_gfa('HB', {inert,fnd}, [], none);
make_gfa('HC', {inert,fnd}, [], none);
make_gfa('HD', {inert,fnd}, [], none);

# LHS of apparent horizon equation
make_gfa('H', {inert,fnd}, [], none);

# 1st derivatives of H[ABCD] and of H
make_gfa('partial_d_HA', {inert,fnd}, [1..N], none);
make_gfa('partial_d_HB', {inert,fnd}, [1..N], none);
make_gfa('partial_d_HC', {inert,fnd}, [1..N], none);
make_gfa('partial_d_HD', {inert,fnd}, [1..N], none);
make_gfa('partial_d_H', {inert,fnd}, [1..N], none);

# Jacobian coefficients for H[ABCD]
# these are the partial derivatives of H[ABCD]
# ... wrt Diff(h,x_rs[x])
make_gfa('partial_HA_wrt_partial_d_h', {inert,fnd}, [1..N_ang], none);
make_gfa('partial_HB_wrt_partial_d_h', {inert,fnd}, [1..N_ang], none);
make_gfa('partial_HC_wrt_partial_d_h', {inert,fnd}, [1..N_ang], none);
make_gfa('partial_HD_wrt_partial_d_h', {inert,fnd}, [1..N_ang], none);
# ... wrt Diff(h,x_rs[x],x_rs[y])
make_gfa('partial_HA_wrt_partial_dd_h', {inert,fnd},
	 [1..N_ang, 1..N_ang], symmetric);
make_gfa('partial_HB_wrt_partial_dd_h', {inert,fnd},
	 [1..N_ang, 1..N_ang], symmetric);

# Jacobian coefficients for H itself
make_gfa('partial_H_wrt_partial_d_h', {inert,fnd}, [1..N_ang], none);
make_gfa('partial_H_wrt_partial_dd_h', {inert,fnd},
	 [1..N_ang, 1..N_ang], symmetric);

NULL;
end proc;


################################################################################

#
# This function implements further simplification rules for  Diff()
# based on gridfn properties which are known here but not in ../maple/.
#
# It currently knows about the following simplifications:
# - Diff(g_dd[i,j], x_xyz[k]) --> partial_d_g_dd[k,i,j]
#
# Anything else is returned unchanged.  (To avoid infinite recursion,
# such a return is *unevaluated*.)
#
# Arguments:
# operand = (in) The thing to be differentiated.
# var_seq = (in) (varargs) An expression sequence of the variables to
#			   differentiate with respect to.
#
`Diff/gridfn2` :=
proc(operand)			# varargs
option remember;		# performance optimization
global
  @include "../maple/coords.minc",
  @include "../maple/gfa.minc",
  @include "../gr/gr_gfas.minc";
local var_list, posn;

var_list := [args[2..nargs]];

if ( type(operand, indexed) and (op(0,operand) = 'g_dd')
     and (nops(var_list) = 1) and member(var_list[1],x_xyz_list,'posn') )
   then return partial_d_g_dd[posn, op(operand)];
end if;

# unevaluated return to avoid infinite recursion
return 'Diff'(operand, op(var_list));
end proc;