v0.14.0
Loading...
Searching...
No Matches
Dg_pointer.hpp
Go to the documentation of this file.
1/* A version for pointers. */
2
3#pragma once
4
5namespace FTensor
6{
7 template <class T, int Tensor_Dim01, int Tensor_Dim2>
8 class Dg<T *, Tensor_Dim01, Tensor_Dim2>
9 {
10
11 protected:
12
13 mutable T *restrict
14 data[(Tensor_Dim01 * (Tensor_Dim01 + 1)) / 2][Tensor_Dim2];
15
16 public:
17
18 template <class... U> Dg(U *... d) : data{d...}
19 {
20 static_assert(sizeof...(d) == sizeof(data) / sizeof(T),
21 "Incorrect number of Arguments. Constructor should "
22 "initialize the entire Tensor");
23 }
24
25 Dg() {}
26
27 /* There are two operator(int,int,int)'s, one for non-consts that lets you
28 change the value, and one for consts that doesn't. */
29
30 T &operator()(const int N1, const int N2, const int N3)
31 {
32#ifdef FTENSOR_DEBUG
33 if(N1 >= Tensor_Dim01 || N1 < 0 || N2 >= Tensor_Dim01 || N2 < 0
34 || N3 >= Tensor_Dim2 || N3 < 0)
35 {
36 std::stringstream s;
37 s << "Bad index in Dg<T*," << Tensor_Dim01 << "," << Tensor_Dim2
38 << ">.operator(" << N1 << "," << N2 << "," << N3 << ")"
39 << std::endl;
40 throw std::out_of_range(s.str());
41 }
42#endif
43 return N1 > N2 ? *data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2][N3]
44 : *data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2][N3];
45 }
46
47 T operator()(const int N1, const int N2, const int N3) const
48 {
49#ifdef FTENSOR_DEBUG
50 if(N1 >= Tensor_Dim01 || N1 < 0 || N2 >= Tensor_Dim01 || N2 < 0
51 || N3 >= Tensor_Dim2 || N3 < 0)
52 {
53 std::stringstream s;
54 s << "Bad index in Dg<T*," << Tensor_Dim01 << "," << Tensor_Dim2
55 << ">.operator(" << N1 << "," << N2 << "," << N3 << ") const"
56 << std::endl;
57 throw std::out_of_range(s.str());
58 }
59#endif
60 return N1 > N2 ? *data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2][N3]
61 : *data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2][N3];
62 }
63
64 T *ptr(const int N1, const int N2, const int N3) const
65 {
66#ifdef FTENSOR_DEBUG
67 if(N1 >= Tensor_Dim01 || N1 < 0 || N2 >= Tensor_Dim01 || N2 < 0
68 || N3 >= Tensor_Dim2 || N3 < 0)
69 {
70 std::stringstream s;
71 s << "Bad index in Dg<T," << Tensor_Dim01 << "," << Tensor_Dim2
72 << ">.ptr(" << N1 << "," << N2 << "," << N3 << ")" << std::endl;
73 throw std::out_of_range(s.str());
74 }
75#endif
76 return N1 > N2 ? data[N1 + (N2 * (2 * Tensor_Dim01 - N2 - 1)) / 2][N3]
77 : data[N2 + (N1 * (2 * Tensor_Dim01 - N1 - 1)) / 2][N3];
78 }
79
80 /* These operator()'s are the first part in constructing template
81 expressions. */
82
83 template <char i, char j, char k, int Dim01, int Dim2>
85 operator()(const Index<i, Dim01> index1, const Index<j, Dim01> index2,
86 const Index<k, Dim2> index3)
87 {
88 return Dg_Expr<Dg<T *, Tensor_Dim01, Tensor_Dim2>, T, Dim01, Dim2, i, j,
89 k>(*this);
90 }
91
92 template <char i, char j, char k, int Dim01, int Dim2>
94 operator()(const Index<i, Dim01> index1, const Index<j, Dim01> index2,
95 const Index<k, Dim2> index3) const
96 {
98 i, j, k>(*this);
99 }
100
101 /* These operators are for internal contractions. The commented out
102 versions are more general, but are ambiguous. This means,
103 unfortunately, that you can't do something like A(i,j,j) where i
104 and j have different dimensions. */
105
106 template <char i, char j, int Dim>
109 T, Dim, i>
110 operator()(const Index<i, Dim> index1, const Index<j, Dim> index2,
111 const Index<j, Dim> index3) const
112 {
113 using TensorExpr
115 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this));
116 }
117
118 template <char i, char j, int Dim>
121 T, Dim, i>
122 operator()(const Index<j, Dim> index1, const Index<i, Dim> index2,
123 const Index<j, Dim> index3) const
124 {
125 using TensorExpr
127 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this));
128 }
129
130 // template<char i, char j, int Dim0, int Dim12>
131 // Tensor1_Expr<const Tensor3_contracted_12<Dg
132 // <T*,Tensor_Dim01,Tensor_Dim2>,T,Dim12>,T,Dim0,i>
133 // operator()(const Index<i,Dim0> index1, const Index<j,Dim12> index2,
134 // const Index<j,Dim12> index3) const
135 // {
136 // typedef Tensor3_contracted_12<Dg<T*,Tensor_Dim01,Tensor_Dim2>,
137 // T,Dim12> TensorExpr;
138 // return Tensor1_Expr<TensorExpr,T,Dim0,i>(TensorExpr(*this));
139 // }
140
141 // template<char i, char j, int Dim02, int Dim1>
142 // Tensor1_Expr<const Tensor3_contracted_02<Dg
143 // <T*,Tensor_Dim01,Tensor_Dim2>,T,Dim02>,T,Dim1,i>
144 // operator()(const Index<j,Dim02> index1, const Index<i,Dim1> index2,
145 // const Index<j,Dim02> index3) const
146 // {
147 // typedef Tensor3_contracted_02<Dg<T*,Tensor_Dim01,Tensor_Dim2>,
148 // T,Dim02> TensorExpr;
149 // return Tensor1_Expr<TensorExpr,T,Dim1,i>(TensorExpr(*this));
150 // }
151
152 template <char i, char j, int Dim01, int Dim2>
155 T, Dim2, i>
156 operator()(const Index<j, Dim01> index1, const Index<j, Dim01> index2,
157 const Index<i, Dim2> index3) const
158 {
159 using TensorExpr
161 return Tensor1_Expr<TensorExpr, T, Dim2, i>(TensorExpr(*this));
162 }
163
164 /* This is for expressions where a number is used for one slot, and
165 indices for the others, yielding a Tensor2_Expr or
166 Tensor2_symmetric_Expr. The non-const versions don't actually
167 create a Dg_number_rhs_* object, while the const versions
168 do create a Dg_number_*. */
169
170 /* First slot. */
171
172 template <char i, char j, int N, int Dim1, int Dim2>
174 Dim1, Dim2, i, j>
175 operator()(const Number<N> n1, const Index<i, Dim1> index1,
176 const Index<j, Dim2> index2)
177 {
178 using TensorExpr
181 }
182
183 template <char i, char j, int N, int Dim1, int Dim2>
184 const Tensor2_Expr<
186 Dim1, Dim2, i, j>
187 operator()(const Number<N> n1, const Index<i, Dim1> index1,
188 const Index<j, Dim2> index2) const
189 {
190 using TensorExpr
192 return Tensor2_Expr<TensorExpr, T, Dim1, Dim2, i, j>(TensorExpr(*this));
193 }
194
195 /* Second slot. */
196
197 template <char i, char j, int N, int Dim0, int Dim2>
199 Dim0, Dim2, i, j>
200 operator()(const Index<i, Dim0> index1, const Number<N> n1,
201 const Index<j, Dim2> index2)
202 {
203 using TensorExpr
206 }
207
208 template <char i, char j, int N, int Dim0, int Dim2>
209 const Tensor2_Expr<
211 Dim0, Dim2, i, j>
212 operator()(const Index<i, Dim0> index1, const Number<N> n1,
213 const Index<j, Dim2> index2) const
214 {
215 using TensorExpr
217 return Tensor2_Expr<TensorExpr, T, Dim0, Dim2, i, j>(TensorExpr(*this));
218 }
219
220 /* Third slot. */
221
222 template <char i, char j, int N, int Dim>
225 operator()(const Index<i, Dim> index1, const Index<j, Dim> index2,
226 const Number<N> n1)
227 {
228 using TensorExpr
231 }
232
233 template <char i, char j, int N, int Dim>
236 Dim, i, j>
237 operator()(const Index<i, Dim> index1, const Index<j, Dim> index2,
238 const Number<N> n1) const
239 {
240 using TensorExpr
243 TensorExpr(*this));
244 }
245
246 /* This is for expressions where a number is used for two slots, and
247 an Index for the other, yielding a Tensor1_Expr. The non-const
248 versions don't actually create a Dg_number_rhs_* object,
249 while the const versions do create a Dg_number_*. */
250
251 /* Index in first slot. */
252
253 template <char i, int N1, int N2, int Dim>
255 T, Dim, i>
256 operator()(const Index<i, Dim> index, const Number<N1> n1,
257 const Number<N2> n2)
258 {
259 using TensorExpr
262 }
263
264 template <char i, int N1, int N2, int Dim>
265 const Tensor1_Expr<
267 T, Dim, i>
268 operator()(const Index<i, Dim> index, const Number<N1> n1,
269 const Number<N2> n2) const
270 {
271 using TensorExpr
273 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this));
274 }
275
276 /* Index in second slot. I use the same structures as for the Index
277 in the first slot since the tensor is symmetric on the first two
278 indices. */
279
280 template <char i, int N1, int N2, int Dim>
282 T, Dim, i>
283 operator()(const Number<N1> n1, const Index<i, Dim> index,
284 const Number<N2> n2)
285 {
286 using TensorExpr
289 }
290
291 template <char i, int N1, int N2, int Dim>
292 const Tensor1_Expr<
294 T, Dim, i>
295 operator()(const Number<N1> n1, const Index<i, Dim> index,
296 const Number<N2> n2) const
297 {
298 using TensorExpr
300 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this));
301 }
302
303 /* Index in third slot. */
304
305 template <char i, int N1, int N2, int Dim>
307 T, Dim, i>
308 operator()(const Number<N1> n1, const Number<N2> n2,
309 const Index<i, Dim> index)
310 {
311 using TensorExpr
314 }
315
316 template <char i, int N1, int N2, int Dim>
317 const Tensor1_Expr<
319 T, Dim, i>
320 operator()(const Number<N1> n1, const Number<N2> n2,
321 const Index<i, Dim> index) const
322 {
323 using TensorExpr
325 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this));
326 }
327
328 /* Specializations for using actual numbers instead of Number<>.
329 This is for expressions where an actual number is used for one
330 slot, and indices for the others, yielding a Tensor2_Expr or
331 Tensor2_symmetric_Expr. I only define the const versions. */
332
333 /* First slot. */
334
335 template <char i, char j, int Dim1, int Dim2>
336 const Tensor2_Expr<
338 Dim2, i, j>
339 operator()(const int N, const Index<i, Dim1> index1,
340 const Index<j, Dim2> index2) const
341 {
342 using TensorExpr
345 TensorExpr(*this, N));
346 }
347
348 /* Second slot. */
349
350 template <char i, char j, int Dim0, int Dim2>
351 const Tensor2_Expr<
353 Dim2, i, j>
354 operator()(const Index<i, Dim0> index1, const int N,
355 const Index<j, Dim2> index2) const
356 {
357 using TensorExpr
360 TensorExpr(*this, N));
361 }
362
363 /* Third slot. */
364
365 template <char i, char j, int Dim>
368 i, j>
369 operator()(const Index<i, Dim> index1, const Index<j, Dim> index2,
370 const int N) const
371 {
372 using TensorExpr
375 TensorExpr(*this, N));
376 }
377
378 /* This is for expressions where a numeral is used for two slots, and
379 an Index for the other, yielding a Tensor1_Expr. */
380
381 /* Index in first slot. */
382
383 template <char i, int Dim>
384 const Tensor1_Expr<
386 i>
387 operator()(const Index<i, Dim> index, const int N1, const int N2) const
388 {
389 using TensorExpr
391 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this, N1, N2));
392 }
393
394 /* Index in second slot. I use the same structures as for the Index
395 in the first slot since the tensor is symmetric on the first two
396 indices. */
397
398 template <char i, int Dim>
399 const Tensor1_Expr<
401 i>
402 operator()(const int N1, const Index<i, Dim> index, const int N2) const
403 {
404 using TensorExpr
406 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this, N1, N2));
407 }
408
409 /* Index in third slot. */
410
411 template <char i, int Dim>
412 const Tensor1_Expr<
414 i>
415 operator()(const int N1, const int N2, const Index<i, Dim> index) const
416 {
417 using TensorExpr
419 return Tensor1_Expr<TensorExpr, T, Dim, i>(TensorExpr(*this, N1, N2));
420 }
421
422 /* The ++ operator increments the pointer, not the number that the
423 pointer points to. This allows iterating over a grid. */
424
426 {
427 for(int i = 0; i < (Tensor_Dim01 * (Tensor_Dim01 + 1)) / 2; ++i)
428 for(int j = 0; j < Tensor_Dim2; ++j)
429 ++data[i][j];
430 return *this;
431 }
432
433 private:
434 template <int I>
435 Dg(const Dg<PackPtr<T *, I>, Tensor_Dim01, Tensor_Dim2> &) = delete;
436 };
437
438 template <class T, int Tensor_Dim01, int Tensor_Dim2, int I>
439 class Dg<PackPtr<T *, I>, Tensor_Dim01, Tensor_Dim2>
441
442 public:
443 template <class... U>
444 Dg(U *... d) : Dg<T *, Tensor_Dim01, Tensor_Dim2>(d...) {}
445
446 Dg() : Dg<T *, Tensor_Dim01, Tensor_Dim2>() {}
447
448 /* The ++ operator increments the pointer, not the number that the
449 pointer points to. This allows iterating over a grid. */
450
451 const Dg &operator++() const {
452 for (int i = 0; i < (Tensor_Dim01 * (Tensor_Dim01 + 1)) / 2; ++i)
453 for (int j = 0; j < Tensor_Dim2; ++j)
455 return *this;
456 }
457 };
458 } // namespace FTensor
static Number< 2 > N2
static Number< 1 > N1
Tensor2_symmetric_Expr< Dg_number_rhs_2< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N >, T, Dim, i, j > operator()(const Index< i, Dim > index1, const Index< j, Dim > index2, const Number< N > n1)
Definition: Dg_pointer.hpp:225
T & operator()(const int N1, const int N2, const int N3)
Definition: Dg_pointer.hpp:30
const Tensor2_Expr< const Dg_number_0< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N >, T, Dim0, Dim2, i, j > operator()(const Index< i, Dim0 > index1, const Number< N > n1, const Index< j, Dim2 > index2) const
Definition: Dg_pointer.hpp:212
T * ptr(const int N1, const int N2, const int N3) const
Definition: Dg_pointer.hpp:64
Tensor1_Expr< const Tensor3_contracted_12< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, Dim >, T, Dim, i > operator()(const Index< i, Dim > index1, const Index< j, Dim > index2, const Index< j, Dim > index3) const
Definition: Dg_pointer.hpp:110
const Tensor2_Expr< const Dg_number_0< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N >, T, Dim1, Dim2, i, j > operator()(const Number< N > n1, const Index< i, Dim1 > index1, const Index< j, Dim2 > index2) const
Definition: Dg_pointer.hpp:187
T operator()(const int N1, const int N2, const int N3) const
Definition: Dg_pointer.hpp:47
const Tensor2_Expr< const Dg_numeral_0< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T >, T, Dim1, Dim2, i, j > operator()(const int N, const Index< i, Dim1 > index1, const Index< j, Dim2 > index2) const
Definition: Dg_pointer.hpp:339
Tensor1_Expr< const Tensor3_contracted_02< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, Dim >, T, Dim, i > operator()(const Index< j, Dim > index1, const Index< i, Dim > index2, const Index< j, Dim > index3) const
Definition: Dg_pointer.hpp:122
Tensor1_Expr< Dg_number_rhs_12< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N1, N2 >, T, Dim, i > operator()(const Number< N1 > n1, const Index< i, Dim > index, const Number< N2 > n2)
Definition: Dg_pointer.hpp:283
const Tensor1_Expr< const Dg_number_12< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N1, N2 >, T, Dim, i > operator()(const Number< N1 > n1, const Index< i, Dim > index, const Number< N2 > n2) const
Definition: Dg_pointer.hpp:295
const Tensor1_Expr< const Dg_numeral_01< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T >, T, Dim, i > operator()(const int N1, const int N2, const Index< i, Dim > index) const
Definition: Dg_pointer.hpp:415
const Tensor1_Expr< const Dg_numeral_12< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T >, T, Dim, i > operator()(const Index< i, Dim > index, const int N1, const int N2) const
Definition: Dg_pointer.hpp:387
const Tensor1_Expr< const Dg_number_01< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N1, N2 >, T, Dim, i > operator()(const Number< N1 > n1, const Number< N2 > n2, const Index< i, Dim > index) const
Definition: Dg_pointer.hpp:320
const Tensor1_Expr< const Dg_number_12< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N1, N2 >, T, Dim, i > operator()(const Index< i, Dim > index, const Number< N1 > n1, const Number< N2 > n2) const
Definition: Dg_pointer.hpp:268
Dg(const Dg< PackPtr< T *, I >, Tensor_Dim01, Tensor_Dim2 > &)=delete
const Tensor2_Expr< const Dg_numeral_0< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T >, T, Dim0, Dim2, i, j > operator()(const Index< i, Dim0 > index1, const int N, const Index< j, Dim2 > index2) const
Definition: Dg_pointer.hpp:354
Tensor1_Expr< Dg_number_rhs_12< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N1, N2 >, T, Dim, i > operator()(const Index< i, Dim > index, const Number< N1 > n1, const Number< N2 > n2)
Definition: Dg_pointer.hpp:256
const Tensor2_symmetric_Expr< const Dg_numeral_2< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T >, T, Dim, i, j > operator()(const Index< i, Dim > index1, const Index< j, Dim > index2, const int N) const
Definition: Dg_pointer.hpp:369
const Tensor2_symmetric_Expr< const Dg_number_2< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N >, T, Dim, i, j > operator()(const Index< i, Dim > index1, const Index< j, Dim > index2, const Number< N > n1) const
Definition: Dg_pointer.hpp:237
Tensor1_Expr< Dg_number_rhs_01< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N1, N2 >, T, Dim, i > operator()(const Number< N1 > n1, const Number< N2 > n2, const Index< i, Dim > index)
Definition: Dg_pointer.hpp:308
Tensor2_Expr< Dg_number_rhs_0< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N >, T, Dim0, Dim2, i, j > operator()(const Index< i, Dim0 > index1, const Number< N > n1, const Index< j, Dim2 > index2)
Definition: Dg_pointer.hpp:200
Tensor1_Expr< const Tensor3_contracted_01< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, Dim01 >, T, Dim2, i > operator()(const Index< j, Dim01 > index1, const Index< j, Dim01 > index2, const Index< i, Dim2 > index3) const
Definition: Dg_pointer.hpp:156
const Dg< T *, Tensor_Dim01, Tensor_Dim2 > & operator++() const
Definition: Dg_pointer.hpp:425
Tensor2_Expr< Dg_number_rhs_0< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, N >, T, Dim1, Dim2, i, j > operator()(const Number< N > n1, const Index< i, Dim1 > index1, const Index< j, Dim2 > index2)
Definition: Dg_pointer.hpp:175
Dg_Expr< Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, Dim01, Dim2, i, j, k > operator()(const Index< i, Dim01 > index1, const Index< j, Dim01 > index2, const Index< k, Dim2 > index3)
Definition: Dg_pointer.hpp:85
const Tensor1_Expr< const Dg_numeral_12< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T >, T, Dim, i > operator()(const int N1, const Index< i, Dim > index, const int N2) const
Definition: Dg_pointer.hpp:402
Dg_Expr< const Dg< T *, Tensor_Dim01, Tensor_Dim2 >, T, Dim01, Dim2, i, j, k > operator()(const Index< i, Dim01 > index1, const Index< j, Dim01 > index2, const Index< k, Dim2 > index3) const
Definition: Dg_pointer.hpp:94
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
const double T
Tensors class implemented by Walter Landry.
Definition: FTensor.hpp:51
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
constexpr IntegrationType I
const int N
Definition: speed_test.cpp:3