v0.13.1
Loading...
Searching...
No Matches
Ddg_value.hpp
Go to the documentation of this file.
1/* A general version, not for pointers. */
2
3#pragma once
4
5namespace FTensor
6{
7 template <class T, int Tensor_Dim01, int Tensor_Dim23> class Ddg
8 {
9 T data[(Tensor_Dim01 * (Tensor_Dim01 + 1)) / 2]
10 [(Tensor_Dim23 * (Tensor_Dim23 + 1)) / 2];
11
12 public:
13 /* There are two operator(int,int,int,int)'s, one for non-consts
14 that lets you change the value, and one for consts that
15 doesn't. */
16
17 T &operator()(const int N1, const int N2, const int N3, const int N4)
18 {
19#ifdef FTENSOR_DEBUG
20 if(N1 >= Tensor_Dim01 || N1 < 0 || N2 >= Tensor_Dim01 || N2 < 0
21 || N3 >= Tensor_Dim23 || N3 < 0 || N4 >= Tensor_Dim23 || N4 < 0)
22 {
23 std::stringstream s;
24 s << "Bad index in Dg<T," << Tensor_Dim01 << "," << Tensor_Dim23
25 << ">.operator(" << N1 << "," << N2 << "," << N3 << "," << N4
26 << ")" << std::endl;
27 throw std::out_of_range(s.str());
28 }
29#endif
30 return N1 > N2
31 ? (N3 > N4 ? data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
32 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2]
33 : data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
34 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2])
35 : (N3 > N4 ? data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
36 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2]
37 : data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
38 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2]);
39 }
40
41 T operator()(const int N1, const int N2, const int N3, const int N4) const
42 {
43#ifdef FTENSOR_DEBUG
44 if(N1 >= Tensor_Dim01 || N1 < 0 || N2 >= Tensor_Dim01 || N2 < 0
45 || N3 >= Tensor_Dim23 || N3 < 0 || N4 >= Tensor_Dim23 || N4 < 0)
46 {
47 std::stringstream s;
48 s << "Bad index in Dg<T," << Tensor_Dim01 << "," << Tensor_Dim23
49 << ">.operator(" << N1 << "," << N2 << "," << N3 << "," << N4
50 << ") const" << std::endl;
51 throw std::out_of_range(s.str());
52 }
53#endif
54 return N1 > N2
55 ? (N3 > N4 ? data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
56 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2]
57 : data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
58 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2])
59 : (N3 > N4 ? data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
60 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2]
61 : data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
62 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2]);
63 }
64
65 template <int N1, int N2, int N3, int N4>
66 inline auto &operator()(const Number<N1> &, const Number<N2> &,
67 const Number<N3> &, const Number<N4> &) {
68
69 static_assert(N1 < Tensor_Dim01, "Bad index N1");
70 static_assert(N2 < Tensor_Dim01, "Bad index N2");
71 static_assert(N3 < Tensor_Dim23, "Bad index N3");
72 static_assert(N4 < Tensor_Dim23, "Bad index N4");
73
74 if constexpr (N1 > N2) {
75 if constexpr (N3 > N4)
76 return data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
77 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2];
78 else
79 return data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
80 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2];
81 } else {
82 if constexpr (N3 > N4)
83 return data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
84 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2];
85 else
86 return data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
87 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2];
88 }
89 }
90
91 template <int N1, int N2, int N3, int N4>
92 inline auto operator()(const Number<N1> &, const Number<N2> &,
93 const Number<N3> &, const Number<N4> &) const {
94
95 static_assert(N1 < Tensor_Dim01, "Bad index N1");
96 static_assert(N2 < Tensor_Dim01, "Bad index N2");
97 static_assert(N3 < Tensor_Dim23, "Bad index N3");
98 static_assert(N4 < Tensor_Dim23, "Bad index N4");
99
100 if constexpr (N1 > N2) {
101 if constexpr (N3 > N4)
102 return data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
103 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2];
104 else
105 return data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2]
106 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2];
107 } else {
108 if constexpr (N3 > N4)
109 return data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
110 [N3 + (N4 * (2 * Tensor_Dim23 - N4 - 1)) / 2];
111 else
112 return data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2]
113 [N4 + (N3 * (2 * Tensor_Dim23 - N3 - 1)) / 2];
114 }
115 }
116
117 /* These operator()'s are the first part in constructing template
118 expressions. They can be used to slice off lower dimensional
119 parts. They are not entirely safe, since you can accidently use a
120 higher dimension than what is really allowed (like Dim=5). */
121
122 template <char i, char j, char k, char l, int Dim01, int Dim23>
123 typename std::enable_if<(Tensor_Dim01 >= Dim01 && Tensor_Dim23 >= Dim23),
125 Dim01, Dim23, i, j, k, l>>::type
126 operator()(const Index<i, Dim01>, const Index<j, Dim01>,
127 const Index<k, Dim23>, const Index<l, Dim23>)
128 {
129 return Ddg_Expr<Ddg<T, Tensor_Dim01, Tensor_Dim23>, T, Dim01, Dim23, i,
130 j, k, l>(*this);
131 }
132
133 template <char i, char j, char k, char l, int Dim01, int Dim23>
134 typename std::enable_if<(Tensor_Dim01 >= Dim01 && Tensor_Dim23 >= Dim23),
136 T, Dim01, Dim23, i, j, k, l>>::type
137 operator()(const Index<i, Dim01>, const Index<j, Dim01>,
138 const Index<k, Dim23>, const Index<l, Dim23>) const
139 {
141 Dim23, i, j, k, l>(*this);
142 }
143
144 /* This is for expressions where a number is used for two slots, and
145 an index for the other two, yielding a Tensor2_symmetric_Expr. */
146
147 template <char i, char j, int N0, int N1, int Dim>
148 typename std::enable_if<
149 (Tensor_Dim01 > N0 && Tensor_Dim01 > N1 && Tensor_Dim23 >= Dim),
152 Dim, i, j>>::type
153 operator()(const Number<N0>, const Number<N1>, const Index<i, Dim>,
154 const Index<j, Dim>) const
155 {
156 using TensorExpr
159 TensorExpr(*this));
160 }
161
162 template <char i, char j, int N0, int N1, int Dim>
163 typename std::enable_if<
164 (Tensor_Dim01 > N0 && Tensor_Dim01 > N1 && Tensor_Dim23 >= Dim),
167 Dim, i, j>>::type
168 operator()(const Number<N0>, const Number<N1>, const Index<i, Dim>,
169 const Index<j, Dim>)
170 {
171 using TensorExpr
174 }
175
176 /* This is for expressions where a number is used for one slot, and
177 an index for the other three, yielding a Dg_Expr. */
178
179 template <char i, char j, char k, int N0, int Dim1, int Dim23>
180 typename std::enable_if<
181 (Tensor_Dim01 > N0 && Tensor_Dim01 >= Dim1 && Tensor_Dim23 >= Dim23),
183 Dim23, Dim1, i, j, k>>::type
184 operator()(const Number<N0>, const Index<k, Dim1>, const Index<i, Dim23>,
185 const Index<j, Dim23>) const
186 {
187 using TensorExpr
189 return Dg_Expr<TensorExpr, T, Dim23, Dim1, i, j, k>(TensorExpr(*this));
190 }
191
192 template <char i, char j, char k, int N0, int Dim1, int Dim23>
193 typename std::enable_if<
194 (Tensor_Dim01 > N0 && Tensor_Dim01 >= Dim1 && Tensor_Dim23 >= Dim23),
196 Dim23, Dim1, i, j, k>>::type
197 operator()(const Number<N0>, const Index<k, Dim1>, const Index<i, Dim23>,
198 const Index<j, Dim23>)
199 {
200 using TensorExpr
203 }
204
205 /* This is for expressions where an int (not a Number) is used for
206 two slots, and an index for the other two, yielding a
207 Tensor2_symmetric_Expr. */
208
209 template <char i, char j, int Dim>
210 typename std::enable_if<
211 (Tensor_Dim23 >= Dim),
214 j>>::type
215 operator()(const int N0, const int N1, const Index<i, Dim>,
216 const Index<j, Dim>) const
217 {
218 using TensorExpr
221 TensorExpr(*this, N0, N1));
222 }
223
224 template <char i, char j, int Dim>
225 typename std::enable_if<
226 (Tensor_Dim01 >= Dim),
229 j>>::type
230 operator()(const Index<i, Dim>, const Index<j, Dim>, const int N2,
231 const int N3) const
232 {
233 using TensorExpr
236 TensorExpr(*this, N2, N3));
237 }
238
239 /* int in two slots but yielding a Tensor2 */
240
241 template <char i, char j, int Dim1, int Dim3>
242 typename std::enable_if<
243 (Tensor_Dim01 >= Dim1 && Tensor_Dim23 >= Dim3),
245 T, Dim1, Dim3, i, j>>::type
246 operator()(const int N0, const Index<i, Dim1>, const int N2,
247 const Index<j, Dim3>) const
248 {
249 using TensorExpr
252 TensorExpr(*this, N0, N2));
253 }
254
255 /* int in three slots, Index in the other yielding a Tensor1_Expr. */
256
257 template <char i, int Dim>
258 typename std::enable_if<
259 (Tensor_Dim01 >= Dim),
261 T, Dim, i>>::type
262 operator()(const Index<i, Dim>, const int N1, const int N2, const int N3)
263 {
264 using TensorExpr
267 TensorExpr(*this, N1, N2, N3));
268 }
269
270 template <char i, int Dim>
271 typename std::enable_if<
272 (Tensor_Dim01 >= Dim),
274 T, Dim, i>>::type
275 operator()(const int N1, const Index<i, Dim>, const int N2, const int N3)
276 {
277 using TensorExpr
280 TensorExpr(*this, N1, N2, N3));
281 }
282
283 /* This is for expressions where an int (not a Number) is used for
284 one slot, and an index for the other three, yielding a
285 Dg_Expr. */
286
287 template <char i, char j, char k, int Dim1, int Dim23>
288 typename std::enable_if<
289 (Tensor_Dim01 >= Dim1 && Tensor_Dim23 >= Dim23),
291 Dim23, Dim1, i, j, k>>::type
292 operator()(const int N0, const Index<k, Dim1>, const Index<i, Dim23>,
293 const Index<j, Dim23>) const
294 {
295 using TensorExpr
298 TensorExpr(*this, N0));
299 }
300 };
301}
static Number< 2 > N2
static Number< 1 > N1
static Number< 0 > N0
auto & operator()(const Number< N1 > &, const Number< N2 > &, const Number< N3 > &, const Number< N4 > &)
Definition: Ddg_value.hpp:66
T data[(Tensor_Dim01 *(Tensor_Dim01+1))/2][(Tensor_Dim23 *(Tensor_Dim23+1))/2]
Definition: Ddg_value.hpp:10
T operator()(const int N1, const int N2, const int N3, const int N4) const
Definition: Ddg_value.hpp:41
auto operator()(const Number< N1 > &, const Number< N2 > &, const Number< N3 > &, const Number< N4 > &) const
Definition: Ddg_value.hpp:92
T & operator()(const int N1, const int N2, const int N3, const int N4)
Definition: Ddg_value.hpp:17
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
const double T
Tensors class implemented by Walter Landry.
Definition: FTensor.hpp:51