11#include " bout/christoffel_symbols.hxx"
22#include " bout/coordinates.hxx"
3- #include " bout/derivs.hxx"
43#include " bout/mesh.hxx"
54
65ChristoffelSymbols::ChristoffelSymbols (Coordinates& coordinates) {
@@ -25,63 +24,83 @@ ChristoffelSymbols::ChristoffelSymbols(Coordinates& coordinates) {
2524 const auto & g_13 = covariantMetricTensor.g13 ();
2625 const auto & g_23 = covariantMetricTensor.g23 ();
2726
28- G1_11_m = 0.5 * g11 * DDX (g_11) + g12 * (DDX (g_12) - 0.5 * DDY (g_11))
29- + g13 * (DDX (g_13) - 0.5 * DDZ (g_11));
30- G1_22_m = g11 * (DDY (g_12) - 0.5 * DDX (g_22)) + 0.5 * g12 * DDY (g_22)
31- + g13 * (DDY (g_23) - 0.5 * DDZ (g_22));
32- G1_33_m = g11 * (DDZ (g_13) - 0.5 * DDX (g_33)) + g12 * (DDZ (g_23) - 0.5 * DDY (g_33))
33- + 0.5 * g13 * DDZ (g_33);
34- G1_12_m = 0.5 * g11 * DDY (g_11) + 0.5 * g12 * DDX (g_22)
35- + 0.5 * g13 * (DDY (g_13) + DDX (g_23) - DDZ (g_12));
36- G1_13_m = 0.5 * g11 * DDZ (g_11) + 0.5 * g12 * (DDZ (g_12) + DDX (g_23) - DDY (g_13))
37- + 0.5 * g13 * DDX (g_33);
38- G1_23_m = 0.5 * g11 * (DDZ (g_12) + DDY (g_13) - DDX (g_23))
39- + 0.5 * g12 * (DDZ (g_22) + DDY (g_23) - DDY (g_23))
40- // + 0.5 *g13*(DDZ(g_32) + DDY(g_33) - DDZ(g_23));
41- // which equals
42- + 0.5 * g13 * DDY (g_33);
27+ G1_11_m = 0.5 * g11 * coordinates.DDX (g_11)
28+ + g12 * (coordinates.DDX (g_12) - 0.5 * coordinates.DDY (g_11))
29+ + g13 * (coordinates.DDX (g_13) - 0.5 * coordinates.DDZ (g_11));
30+ G1_22_m = g11 * (coordinates.DDY (g_12) - 0.5 * coordinates.DDX (g_22))
31+ + 0.5 * g12 * coordinates.DDY (g_22)
32+ + g13 * (coordinates.DDY (g_23) - 0.5 * coordinates.DDZ (g_22));
33+ G1_33_m = g11 * (coordinates.DDZ (g_13) - 0.5 * coordinates.DDX (g_33))
34+ + g12 * (coordinates.DDZ (g_23) - 0.5 * coordinates.DDY (g_33))
35+ + 0.5 * g13 * coordinates.DDZ (g_33);
36+ G1_12_m =
37+ 0.5 * g11 * coordinates.DDY (g_11) + 0.5 * g12 * coordinates.DDX (g_22)
38+ + 0.5 * g13
39+ * (coordinates.DDY (g_13) + coordinates.DDX (g_23) - coordinates.DDZ (g_12));
40+ G1_13_m =
41+ 0.5 * g11 * coordinates.DDZ (g_11)
42+ + 0.5 * g12
43+ * (coordinates.DDZ (g_12) + coordinates.DDX (g_23) - coordinates.DDY (g_13))
44+ + 0.5 * g13 * coordinates.DDX (g_33);
45+ G1_23_m =
46+ 0.5 * g11 * (coordinates.DDZ (g_12) + coordinates.DDY (g_13) - coordinates.DDX (g_23))
47+ + 0.5 * g12
48+ * (coordinates.DDZ (g_22) + coordinates.DDY (g_23) - coordinates.DDY (g_23))
49+ // + 0.5 *g13*(coordinates.DDZ(g_32) + coordinates.DDY(g_33) - coordinates.DDZ(g_23));
50+ // which equals
51+ + 0.5 * g13 * coordinates.DDY (g_33);
4352
44- G2_11_m = 0.5 * g12 * DDX (g_11) + g22 * (DDX (g_12) - 0.5 * DDY (g_11))
45- + g23 * (DDX (g_13) - 0.5 * DDZ (g_11));
46- G2_22_m = g12 * (DDY (g_12) - 0.5 * DDX (g_22)) + 0.5 * g22 * DDY (g_22)
47- + g23 * (DDY (g23) - 0.5 * DDZ (g_22));
48- G2_33_m = g12 * (DDZ (g_13) - 0.5 * DDX (g_33)) + g22 * (DDZ (g_23) - 0.5 * DDY (g_33))
49- + 0.5 * g23 * DDZ (g_33);
50- G2_12_m = 0.5 * g12 * DDY (g_11) + 0.5 * g22 * DDX (g_22)
51- + 0.5 * g23 * (DDY (g_13) + DDX (g_23) - DDZ (g_12));
53+ G2_11_m = 0.5 * g12 * coordinates.DDX (g_11)
54+ + g22 * (coordinates.DDX (g_12) - 0.5 * coordinates.DDY (g_11))
55+ + g23 * (coordinates.DDX (g_13) - 0.5 * coordinates.DDZ (g_11));
56+ G2_22_m = g12 * (coordinates.DDY (g_12) - 0.5 * coordinates.DDX (g_22))
57+ + 0.5 * g22 * coordinates.DDY (g_22)
58+ + g23 * (coordinates.DDY (g23) - 0.5 * coordinates.DDZ (g_22));
59+ G2_33_m = g12 * (coordinates.DDZ (g_13) - 0.5 * coordinates.DDX (g_33))
60+ + g22 * (coordinates.DDZ (g_23) - 0.5 * coordinates.DDY (g_33))
61+ + 0.5 * g23 * coordinates.DDZ (g_33);
62+ G2_12_m =
63+ 0.5 * g12 * coordinates.DDY (g_11) + 0.5 * g22 * coordinates.DDX (g_22)
64+ + 0.5 * g23
65+ * (coordinates.DDY (g_13) + coordinates.DDX (g_23) - coordinates.DDZ (g_12));
5266 G2_13_m =
53- // 0.5 *g21*(DDZ(g_11) + DDX(g13 ()) - DDX(g_13))
67+ // 0.5 *g21*(coordinates. DDZ(g_11) + coordinates. DDX(covariantMetricTensor.Getg13 ()) - coordinates. DDX(g_13))
5468 // which equals
55- 0.5 * g12 * (DDZ (g_11) + DDX (g_13) - DDX (g_13))
56- // + 0.5 *g22*(DDZ(g21 ()) + DDX(g_23) - DDY(g_13))
69+ 0.5 * g12 * (coordinates. DDZ (g_11) + coordinates. DDX (g_13) - coordinates. DDX (g_13))
70+ // + 0.5 *g22*(coordinates. DDZ(covariantMetricTensor.Getg21 ()) + coordinates. DDX(g_23) - coordinates. DDY(g_13))
5771 // which equals
58- + 0.5 * g22 * (DDZ (g_12) + DDX (g_23) - DDY (g_13))
59- // + 0.5 *g23*(DDZ(g31()) + DDX(g_33) - DDZ(g_13));
72+ + 0.5 * g22
73+ * (coordinates.DDZ (g_12) + coordinates.DDX (g_23) - coordinates.DDY (g_13))
74+ // + 0.5 *g23*(coordinates.DDZ(covariantMetricTensor.Getg31()) + coordinates.DDX(g_33) - coordinates.DDZ(g_13));
6075 // which equals
61- + 0.5 * g23 * DDX (g_33);
62- G2_23_m = 0.5 * g12 * (DDZ (g_12) + DDY (g_13) - DDX (g_23)) + 0.5 * g22 * DDZ (g_22)
63- + 0.5 * g23 * DDY (g_33);
76+ + 0.5 * g23 * coordinates.DDX (g_33);
77+ G2_23_m =
78+ 0.5 * g12 * (coordinates.DDZ (g_12) + coordinates.DDY (g_13) - coordinates.DDX (g_23))
79+ + 0.5 * g22 * coordinates.DDZ (g_22) + 0.5 * g23 * coordinates.DDY (g_33);
6480
65- G3_11_m = 0.5 * g13 * DDX (g_11) + g23 * (DDX (g_12) - 0.5 * DDY (g_11))
66- + g33 * (DDX (g_13) - 0.5 * DDZ (g_11));
67- G3_22_m = g13 * (DDY (g_12) - 0.5 * DDX (g_22)) + 0.5 * g23 * DDY (g_22)
68- + g33 * (DDY (g_23) - 0.5 * DDZ (g_22));
69- G3_33_m = g13 * (DDZ (g_13) - 0.5 * DDX (g_33)) + g23 * (DDZ (g_23) - 0.5 * DDY (g_33))
70- + 0.5 * g33 * DDZ (g_33);
81+ G3_11_m = 0.5 * g13 * coordinates.DDX (g_11)
82+ + g23 * (coordinates.DDX (g_12) - 0.5 * coordinates.DDY (g_11))
83+ + g33 * (coordinates.DDX (g_13) - 0.5 * coordinates.DDZ (g_11));
84+ G3_22_m = g13 * (coordinates.DDY (g_12) - 0.5 * coordinates.DDX (g_22))
85+ + 0.5 * g23 * coordinates.DDY (g_22)
86+ + g33 * (coordinates.DDY (g_23) - 0.5 * coordinates.DDZ (g_22));
87+ G3_33_m = g13 * (coordinates.DDZ (g_13) - 0.5 * coordinates.DDX (g_33))
88+ + g23 * (coordinates.DDZ (g_23) - 0.5 * coordinates.DDY (g_33))
89+ + 0.5 * g33 * coordinates.DDZ (g_33);
7190 G3_12_m =
72- // 0.5 *g31*(DDY(g_11) + DDX(g12 ()) - DDX(g_12))
91+ // 0.5 *g31*(coordinates. DDY(g_11) + coordinates. DDX(covariantMetricTensor.Getg12 ()) - coordinates. DDX(g_12))
7392 // which equals to
74- 0.5 * g13 * DDY (g_11)
75- // + 0.5 *g32*(DDY(g21 ()) + DDX(g_22) - DDY(g_12))
93+ 0.5 * g13 * coordinates. DDY (g_11)
94+ // + 0.5 *g32*(coordinates. DDY(covariantMetricTensor.Getg21 ()) + coordinates. DDX(g_22) - coordinates. DDY(g_12))
7695 // which equals to
77- + 0.5 * g23 * DDX (g_22)
78- // + 0.5 *g33*(DDY(g31 ()) + DDX(g32 ()) - DDZ(g_12));
96+ + 0.5 * g23 * coordinates. DDX (g_22)
97+ // + 0.5 *g33*(coordinates. DDY(covariantMetricTensor.Getg31 ()) + coordinates. DDX(covariantMetricTensor.Getg32 ()) - coordinates. DDZ(g_12));
7998 // which equals to
80- + 0.5 * g33 * (DDY (g_13)) + DDX (g_23) - DDZ (g_12);
81- G3_13_m = 0.5 * g13 * DDZ (g_11) + 0.5 * g23 * (DDZ (g_12) + DDX (g_23) - DDY (g_13))
82- + 0.5 * g33 * DDX (g_33);
83- G3_23_m = 0.5 * g13 * (DDZ (g_12) + DDY (g_13)) - DDX (g_23) + 0.5 * g23 * DDZ (g_22)
84- + 0.5 * g33 * DDY (g_33);
99+ + 0.5 * g33 * (coordinates. DDY (g_13)) + coordinates. DDX (g_23) - coordinates. DDZ (g_12);
100+ G3_13_m = 0.5 * g13 * coordinates. DDZ (g_11) + 0.5 * g23 * (coordinates. DDZ (g_12) + coordinates. DDX (g_23) - coordinates. DDY (g_13))
101+ + 0.5 * g33 * coordinates. DDX (g_33);
102+ G3_23_m = 0.5 * g13 * (coordinates. DDZ (g_12) + coordinates. DDY (g_13)) - coordinates. DDX (g_23) + 0.5 * g23 * coordinates. DDZ (g_22)
103+ + 0.5 * g33 * coordinates. DDY (g_33);
85104
86105 output_progress.write (" \t Communicating connection terms\n " );
87106
0 commit comments