v0.14.0
Loading...
Searching...
No Matches
dd_boundary_Tensor0.hpp
Go to the documentation of this file.
1/* Takes a second derivative of a Tensor0<T*> yielding a Tensor1.
2 This is primarily useful at boundaries, where you have to take
3 one-sided derivatives. */
4
5#pragma once
6
7namespace FTensor
8{
9 template <class T, int Dim, char i, char j> class dd_boundary_Tensor0
10 {
15
16 public:
17 typename promote<T, double>::V operator()(const int N1, const int N2) const
18 {
19 return N1 == N2
20 ? (boundary(N1, 0)
21 ? (*(&a + 2 * d_ijk(N1)) - 2 * *(&a + d_ijk(N1)) + a)
22 : (boundary(N1, 1)
23 ? (a - 2 * *(&a - d_ijk(N1)) + *(&a - 2 * d_ijk(N1)))
24 : (*(&a + d_ijk(N1)) - 2 * a + *(&a - d_ijk(N1)))))
25 * d_xyz(N1) * d_xyz(N1)
26 : (boundary(N1, 0)
27 ? (boundary(N2, 0)
28 ? ((*(&a + d_ijk(N1) + d_ijk(N2)) - *(&a + d_ijk(N2))
29 - *(&a + d_ijk(N1)) + a)
30 * d_xyz(N1) * d_xyz(N2))
31 : (boundary(N2, 1) ? ((*(&a + d_ijk(N1)) - a
32 - *(&a + d_ijk(N1) - d_ijk(N2))
33 + *(&a - d_ijk(N2)))
34 * d_xyz(N1) * d_xyz(N2))
35 : (*(&a + d_ijk(N1) + d_ijk(N2))
36 - *(&a + d_ijk(N2))
37 - *(&a + d_ijk(N1) - d_ijk(N2))
38 + *(&a - d_ijk(N2)))
39 * d_xyz(N1) * d_xyz(N2) * 0.5))
40 : (boundary(N1, 1)
41 ? (boundary(N2, 0)
42 ? ((*(&a + d_ijk(N2))
43 - *(&a - d_ijk(N1) + d_ijk(N2)) - a
44 + *(&a - d_ijk(N1)))
45 * d_xyz(N1) * d_xyz(N2))
46 : (boundary(N2, 1)
47 ? ((a - *(&a - d_ijk(N1)) - *(&a - d_ijk(N2))
48 + *(&a - d_ijk(N1) - d_ijk(N2)))
49 * d_xyz(N1) * d_xyz(N2))
50 : (*(&a + d_ijk(N2))
51 - *(&a - d_ijk(N1) + d_ijk(N2))
52 - *(&a - d_ijk(N2))
53 + *(&a - d_ijk(N1) - d_ijk(N2)))
54 * d_xyz(N1) * d_xyz(N2) * 0.5))
55 : (boundary(N2, 0)
56 ? ((*(&a + d_ijk(N1) + d_ijk(N2))
57 - *(&a - d_ijk(N1) + d_ijk(N2))
58 - *(&a + d_ijk(N1)) + *(&a - d_ijk(N1)))
59 * d_xyz(N1) * d_xyz(N2) * 0.5)
60 : (boundary(N2, 1)
61 ? ((*(&a + d_ijk(N1)) - *(&a - d_ijk(N1))
62 - *(&a + d_ijk(N1) - d_ijk(N2))
63 + *(&a - d_ijk(N1) - d_ijk(N2)))
64 * d_xyz(N1) * d_xyz(N2) * 0.5)
65 : ((*(&a + d_ijk(N1) + d_ijk(N2))
66 - *(&a - d_ijk(N1) + d_ijk(N2))
67 - *(&a + d_ijk(N1) - d_ijk(N2))
68 + *(&a - d_ijk(N1) - d_ijk(N2)))
69 * d_xyz(N1) * d_xyz(N2) * 0.25)))));
70 }
72 const Tensor1<double, Dim> &D_xyz,
73 const Tensor2<bool, Dim, 2> &Boundary)
74 : a(A), d_ijk(D_ijk), d_xyz(D_xyz), boundary(Boundary)
75 {}
76 };
77
78 template <class T, int Dim, char i, char j>
79 const Tensor2_symmetric_Expr<const dd_boundary_Tensor0<T, Dim, i, j>,
80 typename promote<T, double>::V, Dim, i, j>
81 dd_boundary(const Tensor0<T *> &a, const Index<i, Dim> index3,
82 const Index<j, Dim> index4, const Tensor1<int, Dim> &d_ijk,
83 const Tensor1<double, Dim> &d_xyz,
84 const Tensor2<bool, Dim, 2> &boundary)
85 {
86 using Tensor_Expr = dd_boundary_Tensor0<T, Dim, i, j>;
88 Dim, i, j>(
89 Tensor_Expr(a, d_ijk, d_xyz, boundary));
90 }
91}
static Number< 2 > N2
static Number< 1 > N1
constexpr double a
const Tensor1< int, Dim > & d_ijk
promote< T, double >::V operator()(const int N1, const int N2) const
const Tensor1< double, Dim > & d_xyz
dd_boundary_Tensor0(const Tensor0< T * > &A, const Tensor1< int, Dim > &D_ijk, const Tensor1< double, Dim > &D_xyz, const Tensor2< bool, Dim, 2 > &Boundary)
const Tensor2< bool, Dim, 2 > & boundary
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
Tensors class implemented by Walter Landry.
Definition: FTensor.hpp:51
const Tensor2_symmetric_Expr< const dd_boundary_Tensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd_boundary(const Tensor0< T * > &a, const Index< i, Dim > index3, const Index< j, Dim > index4, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz, const Tensor2< bool, Dim, 2 > &boundary)
constexpr AssemblyType A