v0.14.0
Riemann.hpp
Go to the documentation of this file.
1 /* Declaration for Riemann, which has the symmetries of the
2  Riemann tensor. Antisymmetric on the first two indices and the
3  last two indices, and symmetric on a cyclic permutation of the last
4  three indices. */
5 
6 #pragma once
7 
8 namespace FTensor
9 {
10  template <class T, int Dim> class Riemann
11  {};
12 
13  template <class T> class Riemann<T, 3>
14  {
15  /* The zero variable is because some of the components of a tensor
16  with these symmetries are identically zero. */
17 
18  T zero;
19  T data0101, data0102, data0112, data0202, data0212, data1212;
20 
21  public:
22  Riemann() : zero(0.0) {}
23 
24  /* There are two operator(int,int,int,int)'s, one for non-consts
25  that lets you change the value, and one for consts that doesn't.
26  The non-const one will give you a wrong answer if you aren't
27  careful. The problem is that we only store the minimal set of
28  components, but some have different signs. We can't return the
29  negative of a component, and assign something to it, because that
30  would assign something to a temporary. To get the correct answer
31  if you don't want to change the value, use eval instead. */
32 
33  T operator()(const int N1, const int N2, const int N3, const int N4) const
34  {
35  return N1 == 0
36  ? (N2 == 0
37  ? zero
38  : (N2 == 1
39  ? (N3 == 0
40  ? (N4 == 0 ? zero
41  : (N4 == 1 ? data0101 : data0102))
42  : (N3 == 1
43  ? (N4 == 0 ? -data0101
44  : (N4 == 1 ? zero : data0112))
45  : (N4 == 0 ? -data0102
46  : (N4 == 1 ? -data0112 : zero))))
47  : (N3 == 0
48  ? (N4 == 0 ? zero
49  : (N4 == 1 ? data0102 : data0202))
50  : (N3 == 1
51  ? (N4 == 0 ? -data0102
52  : (N4 == 1 ? zero : data0212))
53  : (N4 == 0
54  ? -data0202
55  : (N4 == 1 ? -data0212 : zero))))))
56  : (N1 == 1
57  ? (N2 == 0
58  ? (N3 == 0
59  ? (N4 == 0 ? zero
60  : (N4 == 1 ? -data0101 : -data0102))
61  : (N3 == 1
62  ? (N4 == 0 ? data0101
63  : (N4 == 1 ? zero : -data0112))
64  : (N4 == 0 ? data0102
65  : (N4 == 1 ? data0112 : zero))))
66  : (N2 == 1
67  ? zero
68  : (N3 == 0
69  ? (N4 == 0 ? zero
70  : (N4 == 1 ? data0112 : data0212))
71  : (N3 == 1
72  ? (N4 == 0
73  ? -data0112
74  : (N4 == 1 ? zero : data1212))
75  : (N4 == 0 ? -data0212
76  : (N4 == 1 ? -data1212
77  : zero))))))
78  : (N2 == 0
79  ? (N3 == 0
80  ? (N4 == 0 ? zero
81  : (N4 == 1 ? -data0102 : -data0202))
82  : (N3 == 1
83  ? (N4 == 0 ? data0102
84  : (N4 == 1 ? zero : -data0212))
85  : (N4 == 0 ? data0202
86  : (N4 == 1 ? data0212 : zero))))
87  : (N2 == 1
88  ? (N3 == 0
89  ? (N4 == 0
90  ? zero
91  : (N4 == 1 ? -data0112 : -data0212))
92  : (N3 == 1
93  ? (N4 == 0
94  ? data0112
95  : (N4 == 1 ? zero : -data1212))
96  : (N4 == 0
97  ? data0212
98  : (N4 == 1 ? data1212 : zero))))
99  : zero)));
100  }
101 
102  T eval(const int N1, const int N2, const int N3, const int N4) const
103  {
104  return N1 == 0
105  ? (N2 == 0
106  ? zero
107  : (N2 == 1
108  ? (N3 == 0
109  ? (N4 == 0 ? zero
110  : (N4 == 1 ? data0101 : data0102))
111  : (N3 == 1
112  ? (N4 == 0 ? -data0101
113  : (N4 == 1 ? zero : data0112))
114  : (N4 == 0 ? -data0102
115  : (N4 == 1 ? -data0112 : zero))))
116  : (N3 == 0
117  ? (N4 == 0 ? zero
118  : (N4 == 1 ? data0102 : data0202))
119  : (N3 == 1
120  ? (N4 == 0 ? -data0102
121  : (N4 == 1 ? zero : data0212))
122  : (N4 == 0
123  ? -data0202
124  : (N4 == 1 ? -data0212 : zero))))))
125  : (N1 == 1
126  ? (N2 == 0
127  ? (N3 == 0
128  ? (N4 == 0 ? zero
129  : (N4 == 1 ? -data0101 : -data0102))
130  : (N3 == 1
131  ? (N4 == 0 ? data0101
132  : (N4 == 1 ? zero : -data0112))
133  : (N4 == 0 ? data0102
134  : (N4 == 1 ? data0112 : zero))))
135  : (N2 == 1
136  ? zero
137  : (N3 == 0
138  ? (N4 == 0 ? zero
139  : (N4 == 1 ? data0112 : data0212))
140  : (N3 == 1
141  ? (N4 == 0
142  ? -data0112
143  : (N4 == 1 ? zero : data1212))
144  : (N4 == 0 ? -data0212
145  : (N4 == 1 ? -data1212
146  : zero))))))
147  : (N2 == 0
148  ? (N3 == 0
149  ? (N4 == 0 ? zero
150  : (N4 == 1 ? -data0102 : -data0202))
151  : (N3 == 1
152  ? (N4 == 0 ? data0102
153  : (N4 == 1 ? zero : -data0212))
154  : (N4 == 0 ? data0202
155  : (N4 == 1 ? data0212 : zero))))
156  : (N2 == 1
157  ? (N3 == 0
158  ? (N4 == 0
159  ? zero
160  : (N4 == 1 ? -data0112 : -data0212))
161  : (N3 == 1
162  ? (N4 == 0
163  ? data0112
164  : (N4 == 1 ? zero : -data1212))
165  : (N4 == 0
166  ? data0212
167  : (N4 == 1 ? data1212 : zero))))
168  : zero)));
169  }
170 
171  T &operator()(const int N1, const int N2, const int N3, const int N4)
172  {
173  return (N1 == 0 && N2 == 1 && N3 == 0 && N4 == 1)
174  ? data0101
175  : ((N1 == 0 && N2 == 1 && N3 == 0 && N4 == 2)
176  ? data0102
177  : ((N1 == 0 && N2 == 1 && N3 == 1 && N4 == 2)
178  ? data0112
179  : ((N1 == 0 && N2 == 2 && N3 == 1 && N4 == 2)
180  ? data0212
181  : ((N1 == 0 && N2 == 2 && N3 == 0 && N4 == 2)
182  ? data0202
183  : ((N1 == 1 && N2 == 2 && N3 == 1 && N4 == 2)
184  ? data1212
185  : zero)))));
186  }
187 
188  /* These operator()'s are the first part in constructing template
189  expressions. */
190 
191  template <char i, char j, char k, char l>
192  Riemann_Expr<Riemann<T, 3>, T, 3, i, j, k, l>
193  operator()(const Index<i, 3>, const Index<j, 3>, const Index<k, 3>,
194  const Index<l, 3>)
195  {
196  return Riemann_Expr<Riemann<T, 3>, T, 3, i, j, k, l>(*this);
197  }
198  };
199 }
200 #include "Riemann/Riemann_Expr.hpp"
FTensor
JSON compatible output.
Definition: Christof_constructor.hpp:6
FTensor::Riemann< T, 3 >::Riemann
Riemann()
Definition: Riemann.hpp:22
FTensor::Riemann< T, 3 >::data1212
T data1212
Definition: Riemann.hpp:19
FTensor::Riemann_Expr
Definition: Riemann_Expr.hpp:15
FTensor::Riemann
Definition: Riemann.hpp:10
FTensor::Riemann< T, 3 >::zero
T zero
Definition: Riemann.hpp:18
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
FTensor::Index
Definition: Index.hpp:23
Riemann_Expr.hpp
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
FTensor::Riemann< T, 3 >::eval
T eval(const int N1, const int N2, const int N3, const int N4) const
Definition: Riemann.hpp:102
FTensor::Riemann< T, 3 >::operator()
Riemann_Expr< Riemann< T, 3 >, T, 3, i, j, k, l > operator()(const Index< i, 3 >, const Index< j, 3 >, const Index< k, 3 >, const Index< l, 3 >)
Definition: Riemann.hpp:193
FTensor::Riemann< T, 3 >::operator()
T & operator()(const int N1, const int N2, const int N3, const int N4)
Definition: Riemann.hpp:171
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
FTensor::Riemann< T, 3 >::operator()
T operator()(const int N1, const int N2, const int N3, const int N4) const
Definition: Riemann.hpp:33
Generated by Doxygen 1.8.17 and hosted at University of Glasgow