Shamrock 2025.10.0
Astrophysical Code
Loading...
Searching...
No Matches
offset_multipole.hpp
Go to the documentation of this file.
1// -------------------------------------------------------//
2//
3// SHAMROCK code for hydrodynamics
4// Copyright (c) 2021-2026 Timothée David--Cléris <tim.shamrock@proton.me>
5// SPDX-License-Identifier: CeCILL Free Software License Agreement v2.1
6// Shamrock is licensed under the CeCILL 2.1 License, see LICENSE for more information
7//
8// -------------------------------------------------------//
9
10#pragma once
11
21
22namespace shamphys {
23
24 namespace details {
25
26 template<class T>
27 inline shammath::SymTensor3d_1<T> offset_multipole_1(
28 const T &Qt0,
30 const shammath::SymTensor3d_1<T> &dt1) {
31 return Qt1 + Qt0 * dt1;
32 }
33
34 template<class T>
35 inline shammath::SymTensor3d_2<T> offset_multipole_2(
39 const shammath::SymTensor3d_1<T> &dt1) {
40 // mathematica out
41 // auto Qn2 = SymTensor3d_2<T>{
42 // Q.t2.v_00 + 2*Q.t1.v_0*d.t1.v_0 + Q.t0*d.t2.v_00,
43 // Q.t2.v_01 + Q.t1.v_1*d.t1.v_0 + Q.t1.v_0*d.t1.v_1 + Q.t0*d.t2.v_01,
44 // Q.t2.v_02 + Q.t1.v_2*d.t1.v_0 + Q.t1.v_0*d.t1.v_2 + Q.t0*d.t2.v_02,
45 // Q.t2.v_11 + 2*Q.t1.v_1*d.t1.v_1 + Q.t0*d.t2.v_11,
46 // Q.t2.v_12 + Q.t1.v_2*d.t1.v_1 + Q.t1.v_1*d.t1.v_2 + Q.t0*d.t2.v_12,
47 // Q.t2.v_22 + 2*Q.t1.v_2*d.t1.v_2 + Q.t0*d.t2.v_22
48 // };
49
50 // symbolic Qn2 : d_\mu Qn1_\nu + Q1_\mu d_\nu + Q2
52 dt1.v_0 * Qn1.v_0 + Qt1.v_0 * dt1.v_0 + Qt2.v_00,
53 dt1.v_0 * Qn1.v_1 + Qt1.v_0 * dt1.v_1 + Qt2.v_01,
54 dt1.v_0 * Qn1.v_2 + Qt1.v_0 * dt1.v_2 + Qt2.v_02,
55 dt1.v_1 * Qn1.v_1 + Qt1.v_1 * dt1.v_1 + Qt2.v_11,
56 dt1.v_1 * Qn1.v_2 + Qt1.v_1 * dt1.v_2 + Qt2.v_12,
57 dt1.v_2 * Qn1.v_2 + Qt1.v_2 * dt1.v_2 + Qt2.v_22};
58 }
59
60 template<class T>
61 inline shammath::SymTensor3d_3<T> offset_multipole_3(
66 const shammath::SymTensor3d_1<T> &dt1) {
67
68 // mathematica out
69 // auto Qn3 = SymTensor3d_3<T>{
70 // Qt3.v_000 + 3*Qt2.v_00*dt1.v_0 + 3*Q.t1.v_0*d.t2.v_00 + Q.t0*d.t3.v_000,
71 // Qt3.v_001 + 2*Qt2.v_01*dt1.v_0 + Qt2.v_00*dt1.v_1 + Q.t1.v_1*d.t2.v_00 +
72 // 2*Q.t1.v_0*d.t2.v_01 + Q.t0*d.t3.v_001, Qt3.v_002 + 2*Qt2.v_02*dt1.v_0 +
73 // Qt2.v_00*dt1.v_2 + Q.t1.v_2*d.t2.v_00 + 2*Q.t1.v_0*d.t2.v_02 + Q.t0*d.t3.v_002,
74 // Qt3.v_011 + Qt2.v_11*dt1.v_0 + 2*Qt2.v_01*dt1.v_1 + 2*Q.t1.v_1*d.t2.v_01 +
75 // Q.t1.v_0*d.t2.v_11 + Q.t0*d.t3.v_011, Qt3.v_012 + Qt2.v_12*dt1.v_0 +
76 // Qt2.v_02*dt1.v_1 + Qt2.v_01*dt1.v_2 + Q.t1.v_2*d.t2.v_01 + Q.t1.v_1*d.t2.v_02
77 // + Q.t1.v_0*d.t2.v_12 + Q.t0*d.t3.v_012, Qt3.v_022 + Qt2.v_22*dt1.v_0 +
78 // 2*Qt2.v_02*dt1.v_2 + 2*Q.t1.v_2*d.t2.v_02 + Q.t1.v_0*d.t2.v_22 +
79 // Q.t0*d.t3.v_022, Qt3.v_111 + 3*Qt2.v_11*dt1.v_1 + 3*Q.t1.v_1*d.t2.v_11 +
80 // Q.t0*d.t3.v_111, Qt3.v_112 + 2*Qt2.v_12*dt1.v_1 + Qt2.v_11*dt1.v_2 +
81 // Q.t1.v_2*d.t2.v_11 + 2*Q.t1.v_1*d.t2.v_12 + Q.t0*d.t3.v_112, Qt3.v_122 +
82 // Qt2.v_22*dt1.v_1 + 2*Qt2.v_12*dt1.v_2 + 2*Q.t1.v_2*d.t2.v_12 +
83 // Q.t1.v_1*d.t2.v_22 + Q.t0*d.t3.v_122, Qt3.v_222 + 3*Qt2.v_22*dt1.v_2 +
84 // 3*Q.t1.v_2*d.t2.v_22 + Q.t0*d.t3.v_222
85 // };
86
87 // symbolic Qn3 : d_\mu Qn2_\nu\delta + d_\nu * (Qn2_ \mu \delta - d_\mu Qn1_\delta) +
88 // Q2_\mu\nu d_\delta + Q3
90 dt1.v_0 * Qn2.v_00 + dt1.v_0 * (Qn2.v_00 - dt1.v_0 * Qn1.v_0) + Qt2.v_00 * dt1.v_0
91 + Qt3.v_000,
92 dt1.v_0 * Qn2.v_01 + dt1.v_0 * (Qn2.v_01 - dt1.v_0 * Qn1.v_1) + Qt2.v_00 * dt1.v_1
93 + Qt3.v_001,
94 dt1.v_0 * Qn2.v_02 + dt1.v_0 * (Qn2.v_02 - dt1.v_0 * Qn1.v_2) + Qt2.v_00 * dt1.v_2
95 + Qt3.v_002,
96 dt1.v_0 * Qn2.v_11 + dt1.v_1 * (Qn2.v_01 - dt1.v_0 * Qn1.v_1) + Qt2.v_01 * dt1.v_1
97 + Qt3.v_011,
98 dt1.v_0 * Qn2.v_12 + dt1.v_1 * (Qn2.v_02 - dt1.v_0 * Qn1.v_2) + Qt2.v_01 * dt1.v_2
99 + Qt3.v_012,
100 dt1.v_0 * Qn2.v_22 + dt1.v_2 * (Qn2.v_02 - dt1.v_0 * Qn1.v_2) + Qt2.v_02 * dt1.v_2
101 + Qt3.v_022,
102 dt1.v_1 * Qn2.v_11 + dt1.v_1 * (Qn2.v_11 - dt1.v_1 * Qn1.v_1) + Qt2.v_11 * dt1.v_1
103 + Qt3.v_111,
104 dt1.v_1 * Qn2.v_12 + dt1.v_1 * (Qn2.v_12 - dt1.v_1 * Qn1.v_2) + Qt2.v_11 * dt1.v_2
105 + Qt3.v_112,
106 dt1.v_1 * Qn2.v_22 + dt1.v_2 * (Qn2.v_12 - dt1.v_1 * Qn1.v_2) + Qt2.v_12 * dt1.v_2
107 + Qt3.v_122,
108 dt1.v_2 * Qn2.v_22 + dt1.v_2 * (Qn2.v_22 - dt1.v_2 * Qn1.v_2) + Qt2.v_22 * dt1.v_2
109 + Qt3.v_222};
110 }
111
112 template<class T>
113 inline shammath::SymTensor3d_4<T> offset_multipole_4(
120 const shammath::SymTensor3d_2<T> &dt2) {
121
122 // auto Qn4 = SymTensor3d_4<T>{
123 // Q.t4.v_0000 + 4*Q.t3.v_000*d.t1.v_0 + 6*Q.t2.v_00*d.t2.v_00 +
124 // 4*Q.t1.v_0*d.t3.v_000 + Q.t0*d.t4.v_0000, Q.t4.v_0001 + 3*Q.t3.v_001*d.t1.v_0 +
125 // Q.t3.v_000*d.t1.v_1 + 3*Q.t2.v_01*d.t2.v_00 + 3*Q.t2.v_00*d.t2.v_01 +
126 // Q.t1.v_1*d.t3.v_000 + 3*Q.t1.v_0*d.t3.v_001 + Q.t0*d.t4.v_0001, Q.t4.v_0002 +
127 // 3*Q.t3.v_002*d.t1.v_0 + Q.t3.v_000*d.t1.v_2 + 3*Q.t2.v_02*d.t2.v_00 +
128 // 3*Q.t2.v_00*d.t2.v_02 + Q.t1.v_2*d.t3.v_000
129 // + 3*Q.t1.v_0*d.t3.v_002 + Q.t0*d.t4.v_0002, Q.t4.v_0011 + 2*Q.t3.v_011*d.t1.v_0 +
130 // 2*Q.t3.v_001*d.t1.v_1 + Q.t2.v_11*d.t2.v_00 + 4*Q.t2.v_01*d.t2.v_01 +
131 // Q.t2.v_00*d.t2.v_11
132 // + 2*Q.t1.v_1*d.t3.v_001 + 2*Q.t1.v_0*d.t3.v_011 + Q.t0*d.t4.v_0011, Q.t4.v_0012 +
133 // 2*Q.t3.v_012*d.t1.v_0 + Q.t3.v_002*d.t1.v_1 + Q.t3.v_001*d.t1.v_2 +
134 // Q.t2.v_12*d.t2.v_00 + 2*Q.t2.v_02*d.t2.v_01 + 2*Q.t2.v_01*d.t2.v_02 +
135 // Q.t2.v_00*d.t2.v_12 + Q.t1.v_2*d.t3.v_001
136 // + Q.t1.v_1*d.t3.v_002 + 2*Q.t1.v_0*d.t3.v_012 + Q.t0*d.t4.v_0012, Q.t4.v_0022 +
137 // 2*Q.t3.v_022*d.t1.v_0 + 2*Q.t3.v_002*d.t1.v_2 + Q.t2.v_22*d.t2.v_00 +
138 // 4*Q.t2.v_02*d.t2.v_02 + Q.t2.v_00*d.t2.v_22 + 2*Q.t1.v_2*d.t3.v_002 +
139 // 2*Q.t1.v_0*d.t3.v_022 + Q.t0*d.t4.v_0022, Q.t4.v_0111 + Q.t3.v_111*d.t1.v_0 +
140 // 3*Q.t3.v_011*d.t1.v_1 + 3*Q.t2.v_11*d.t2.v_01 + 3*Q.t2.v_01*d.t2.v_11 +
141 // 3*Q.t1.v_1*d.t3.v_011 + Q.t1.v_0*d.t3.v_111 + Q.t0*d.t4.v_0111, Q.t4.v_0112 +
142 // Q.t3.v_112*d.t1.v_0 + 2*Q.t3.v_012*d.t1.v_1 + Q.t3.v_011*d.t1.v_2 +
143 // 2*Q.t2.v_12*d.t2.v_01
144 // + Q.t2.v_11*d.t2.v_02 + Q.t2.v_02*d.t2.v_11 + 2*Q.t2.v_01*d.t2.v_12 +
145 // Q.t1.v_2*d.t3.v_011
146 // + 2*Q.t1.v_1*d.t3.v_012 + Q.t1.v_0*d.t3.v_112 + Q.t0*d.t4.v_0112, Q.t4.v_0122 +
147 // Q.t3.v_122*d.t1.v_0 + Q.t3.v_022*d.t1.v_1 + 2*Q.t3.v_012*d.t1.v_2 +
148 // Q.t2.v_22*d.t2.v_01 + 2*Q.t2.v_12*d.t2.v_02 + 2*Q.t2.v_02*d.t2.v_12 +
149 // Q.t2.v_01*d.t2.v_22 + 2*Q.t1.v_2*d.t3.v_012 + Q.t1.v_1*d.t3.v_022 +
150 // Q.t1.v_0*d.t3.v_122 + Q.t0*d.t4.v_0122, Q.t4.v_0222 + Q.t3.v_222*d.t1.v_0 +
151 // 3*Q.t3.v_022*d.t1.v_2 + 3*Q.t2.v_22*d.t2.v_02 + 3*Q.t2.v_02*d.t2.v_22 +
152 // 3*Q.t1.v_2*d.t3.v_022 + Q.t1.v_0*d.t3.v_222 + Q.t0*d.t4.v_0222, Q.t4.v_1111 +
153 // 4*Q.t3.v_111*d.t1.v_1 + 6*Q.t2.v_11*d.t2.v_11 + 4*Q.t1.v_1*d.t3.v_111 +
154 // Q.t0*d.t4.v_1111, Q.t4.v_1112 + 3*Q.t3.v_112*d.t1.v_1 + Q.t3.v_111*d.t1.v_2 +
155 // 3*Q.t2.v_12*d.t2.v_11 + 3*Q.t2.v_11*d.t2.v_12 + Q.t1.v_2*d.t3.v_111 +
156 // 3*Q.t1.v_1*d.t3.v_112 + Q.t0*d.t4.v_1112, Q.t4.v_1122 + 2*Q.t3.v_122*d.t1.v_1 +
157 // 2*Q.t3.v_112*d.t1.v_2 + Q.t2.v_22*d.t2.v_11 + 4*Q.t2.v_12*d.t2.v_12 +
158 // Q.t2.v_11*d.t2.v_22
159 // + 2*Q.t1.v_2*d.t3.v_112 + 2*Q.t1.v_1*d.t3.v_122 + Q.t0*d.t4.v_1122, Q.t4.v_1222 +
160 // Q.t3.v_222*d.t1.v_1 + 3*Q.t3.v_122*d.t1.v_2 + 3*Q.t2.v_22*d.t2.v_12 +
161 // 3*Q.t2.v_12*d.t2.v_22 + 3*Q.t1.v_2*d.t3.v_122 + Q.t1.v_1*d.t3.v_222 +
162 // Q.t0*d.t4.v_1222, Q.t4.v_2222 + 4*Q.t3.v_222*d.t1.v_2 + 6*Q.t2.v_22*d.t2.v_22 +
163 // 4*Q.t1.v_2*d.t3.v_222 + Q.t0*d.t4.v_2222
164 // };
165
166 // symbolic Qn4 : T4_\mu\nu\delta\epsilon = d_\mu * T3_\nu\delta\epsilon + d_\delta
167 // *(T3_\mu\epsilon\nu - d_\mu * T2_\epsilon\nu ) + d2_\epsilon\nu *Q2_\delta\mu +
168 // d_\epsilon *Q3_\nu\delta\mu + d_\nu *Q3_\delta\epsilon\mu + Q4_\mu\nu\delta\epsilon
169
171 dt1.v_0 * Qn3.v_000 + dt1.v_0 * (Qn3.v_000 - dt1.v_0 * Qn2.v_00)
172 + dt2.v_00 * Qt2.v_00 + dt1.v_0 * Qt3.v_000 + dt1.v_0 * Qt3.v_000 + Qt4.v_0000,
173 dt1.v_0 * Qn3.v_001 + dt1.v_0 * (Qn3.v_001 - dt1.v_0 * Qn2.v_01)
174 + dt2.v_01 * Qt2.v_00 + dt1.v_1 * Qt3.v_000 + dt1.v_0 * Qt3.v_001 + Qt4.v_0001,
175 dt1.v_0 * Qn3.v_002 + dt1.v_0 * (Qn3.v_002 - dt1.v_0 * Qn2.v_02)
176 + dt2.v_02 * Qt2.v_00 + dt1.v_2 * Qt3.v_000 + dt1.v_0 * Qt3.v_002 + Qt4.v_0002,
177 dt1.v_0 * Qn3.v_011 + dt1.v_1 * (Qn3.v_001 - dt1.v_0 * Qn2.v_01)
178 + dt2.v_01 * Qt2.v_01 + dt1.v_1 * Qt3.v_001 + dt1.v_0 * Qt3.v_011 + Qt4.v_0011,
179 dt1.v_0 * Qn3.v_012 + dt1.v_1 * (Qn3.v_002 - dt1.v_0 * Qn2.v_02)
180 + dt2.v_02 * Qt2.v_01 + dt1.v_2 * Qt3.v_001 + dt1.v_0 * Qt3.v_012 + Qt4.v_0012,
181 dt1.v_0 * Qn3.v_022 + dt1.v_2 * (Qn3.v_002 - dt1.v_0 * Qn2.v_02)
182 + dt2.v_02 * Qt2.v_02 + dt1.v_2 * Qt3.v_002 + dt1.v_0 * Qt3.v_022 + Qt4.v_0022,
183 dt1.v_0 * Qn3.v_111 + dt1.v_1 * (Qn3.v_011 - dt1.v_0 * Qn2.v_11)
184 + dt2.v_11 * Qt2.v_01 + dt1.v_1 * Qt3.v_011 + dt1.v_1 * Qt3.v_011 + Qt4.v_0111,
185 dt1.v_0 * Qn3.v_112 + dt1.v_1 * (Qn3.v_012 - dt1.v_0 * Qn2.v_12)
186 + dt2.v_12 * Qt2.v_01 + dt1.v_2 * Qt3.v_011 + dt1.v_1 * Qt3.v_012 + Qt4.v_0112,
187 dt1.v_0 * Qn3.v_122 + dt1.v_2 * (Qn3.v_012 - dt1.v_0 * Qn2.v_12)
188 + dt2.v_12 * Qt2.v_02 + dt1.v_2 * Qt3.v_012 + dt1.v_1 * Qt3.v_022 + Qt4.v_0122,
189 dt1.v_0 * Qn3.v_222 + dt1.v_2 * (Qn3.v_022 - dt1.v_0 * Qn2.v_22)
190 + dt2.v_22 * Qt2.v_02 + dt1.v_2 * Qt3.v_022 + dt1.v_2 * Qt3.v_022 + Qt4.v_0222,
191 dt1.v_1 * Qn3.v_111 + dt1.v_1 * (Qn3.v_111 - dt1.v_1 * Qn2.v_11)
192 + dt2.v_11 * Qt2.v_11 + dt1.v_1 * Qt3.v_111 + dt1.v_1 * Qt3.v_111 + Qt4.v_1111,
193 dt1.v_1 * Qn3.v_112 + dt1.v_1 * (Qn3.v_112 - dt1.v_1 * Qn2.v_12)
194 + dt2.v_12 * Qt2.v_11 + dt1.v_2 * Qt3.v_111 + dt1.v_1 * Qt3.v_112 + Qt4.v_1112,
195 dt1.v_1 * Qn3.v_122 + dt1.v_2 * (Qn3.v_112 - dt1.v_1 * Qn2.v_12)
196 + dt2.v_12 * Qt2.v_12 + dt1.v_2 * Qt3.v_112 + dt1.v_1 * Qt3.v_122 + Qt4.v_1122,
197 dt1.v_1 * Qn3.v_222 + dt1.v_2 * (Qn3.v_122 - dt1.v_1 * Qn2.v_22)
198 + dt2.v_22 * Qt2.v_12 + dt1.v_2 * Qt3.v_122 + dt1.v_2 * Qt3.v_122 + Qt4.v_1222,
199 dt1.v_2 * Qn3.v_222 + dt1.v_2 * (Qn3.v_222 - dt1.v_2 * Qn2.v_22)
200 + dt2.v_22 * Qt2.v_22 + dt1.v_2 * Qt3.v_222 + dt1.v_2 * Qt3.v_222 + Qt4.v_2222};
201 }
202
203 template<class T>
204 inline shammath::SymTensor3d_5<T> offset_multipole_5(
213 const shammath::SymTensor3d_3<T> &dt3) {
214
215 // symbolic Qn5 : {T5}_{\mu\nu\delta\epsilon\sigma} = d_\mu T4_\nu\delta\epsilon\sigma
216 // + d_\nu ( T4_\mu\delta\epsilon\sigma - d_\mu T3_\delta\epsilon\sigma ) + Q2_\mu\nu
217 // d3_\delta\epsilon\sigma +Q3_\mu\nu\sigma d2_\delta\epsilon +Q3_\mu\nu\delta
218 // d2_\epsilon\sigma +Q4_\mu\nu\delta\sigma d_\epsilon +Q3_\mu\nu\epsilon
219 // d2_\delta\sigma +Q4_\mu\nu\epsilon\sigma d_\delta +Q4_\mu\nu\delta\epsilon d_\sigma
220 // +Q5_\mu\nu\delta\epsilon\sigma
221
223 dt1.v_0 * Qn4.v_0000 + dt1.v_0 * (Qn4.v_0000 - dt1.v_0 * Qn3.v_000)
224 + Qt2.v_00 * dt3.v_000 + Qt3.v_000 * dt2.v_00 + Qt3.v_000 * dt2.v_00
225 + Qt4.v_0000 * dt1.v_0 + Qt3.v_000 * dt2.v_00 + Qt4.v_0000 * dt1.v_0
226 + Qt4.v_0000 * dt1.v_0 + Qt5.v_00000,
227 dt1.v_0 * Qn4.v_0001 + dt1.v_0 * (Qn4.v_0001 - dt1.v_0 * Qn3.v_001)
228 + Qt2.v_00 * dt3.v_001 + Qt3.v_001 * dt2.v_00 + Qt3.v_000 * dt2.v_01
229 + Qt4.v_0001 * dt1.v_0 + Qt3.v_000 * dt2.v_01 + Qt4.v_0001 * dt1.v_0
230 + Qt4.v_0000 * dt1.v_1 + Qt5.v_00001,
231 dt1.v_0 * Qn4.v_0002 + dt1.v_0 * (Qn4.v_0002 - dt1.v_0 * Qn3.v_002)
232 + Qt2.v_00 * dt3.v_002 + Qt3.v_002 * dt2.v_00 + Qt3.v_000 * dt2.v_02
233 + Qt4.v_0002 * dt1.v_0 + Qt3.v_000 * dt2.v_02 + Qt4.v_0002 * dt1.v_0
234 + Qt4.v_0000 * dt1.v_2 + Qt5.v_00002,
235 dt1.v_0 * Qn4.v_0011 + dt1.v_0 * (Qn4.v_0011 - dt1.v_0 * Qn3.v_011)
236 + Qt2.v_00 * dt3.v_011 + Qt3.v_001 * dt2.v_01 + Qt3.v_000 * dt2.v_11
237 + Qt4.v_0001 * dt1.v_1 + Qt3.v_001 * dt2.v_01 + Qt4.v_0011 * dt1.v_0
238 + Qt4.v_0001 * dt1.v_1 + Qt5.v_00011,
239 dt1.v_0 * Qn4.v_0012 + dt1.v_0 * (Qn4.v_0012 - dt1.v_0 * Qn3.v_012)
240 + Qt2.v_00 * dt3.v_012 + Qt3.v_002 * dt2.v_01 + Qt3.v_000 * dt2.v_12
241 + Qt4.v_0002 * dt1.v_1 + Qt3.v_001 * dt2.v_02 + Qt4.v_0012 * dt1.v_0
242 + Qt4.v_0001 * dt1.v_2 + Qt5.v_00012,
243 dt1.v_0 * Qn4.v_0022 + dt1.v_0 * (Qn4.v_0022 - dt1.v_0 * Qn3.v_022)
244 + Qt2.v_00 * dt3.v_022 + Qt3.v_002 * dt2.v_02 + Qt3.v_000 * dt2.v_22
245 + Qt4.v_0002 * dt1.v_2 + Qt3.v_002 * dt2.v_02 + Qt4.v_0022 * dt1.v_0
246 + Qt4.v_0002 * dt1.v_2 + Qt5.v_00022,
247 dt1.v_0 * Qn4.v_0111 + dt1.v_0 * (Qn4.v_0111 - dt1.v_0 * Qn3.v_111)
248 + Qt2.v_00 * dt3.v_111 + Qt3.v_001 * dt2.v_11 + Qt3.v_001 * dt2.v_11
249 + Qt4.v_0011 * dt1.v_1 + Qt3.v_001 * dt2.v_11 + Qt4.v_0011 * dt1.v_1
250 + Qt4.v_0011 * dt1.v_1 + Qt5.v_00111,
251 dt1.v_0 * Qn4.v_0112 + dt1.v_0 * (Qn4.v_0112 - dt1.v_0 * Qn3.v_112)
252 + Qt2.v_00 * dt3.v_112 + Qt3.v_002 * dt2.v_11 + Qt3.v_001 * dt2.v_12
253 + Qt4.v_0012 * dt1.v_1 + Qt3.v_001 * dt2.v_12 + Qt4.v_0012 * dt1.v_1
254 + Qt4.v_0011 * dt1.v_2 + Qt5.v_00112,
255 dt1.v_0 * Qn4.v_0122 + dt1.v_0 * (Qn4.v_0122 - dt1.v_0 * Qn3.v_122)
256 + Qt2.v_00 * dt3.v_122 + Qt3.v_002 * dt2.v_12 + Qt3.v_001 * dt2.v_22
257 + Qt4.v_0012 * dt1.v_2 + Qt3.v_002 * dt2.v_12 + Qt4.v_0022 * dt1.v_1
258 + Qt4.v_0012 * dt1.v_2 + Qt5.v_00122,
259 dt1.v_0 * Qn4.v_0222 + dt1.v_0 * (Qn4.v_0222 - dt1.v_0 * Qn3.v_222)
260 + Qt2.v_00 * dt3.v_222 + Qt3.v_002 * dt2.v_22 + Qt3.v_002 * dt2.v_22
261 + Qt4.v_0022 * dt1.v_2 + Qt3.v_002 * dt2.v_22 + Qt4.v_0022 * dt1.v_2
262 + Qt4.v_0022 * dt1.v_2 + Qt5.v_00222,
263 dt1.v_0 * Qn4.v_1111 + dt1.v_1 * (Qn4.v_0111 - dt1.v_0 * Qn3.v_111)
264 + Qt2.v_01 * dt3.v_111 + Qt3.v_011 * dt2.v_11 + Qt3.v_011 * dt2.v_11
265 + Qt4.v_0111 * dt1.v_1 + Qt3.v_011 * dt2.v_11 + Qt4.v_0111 * dt1.v_1
266 + Qt4.v_0111 * dt1.v_1 + Qt5.v_01111,
267 dt1.v_0 * Qn4.v_1112 + dt1.v_1 * (Qn4.v_0112 - dt1.v_0 * Qn3.v_112)
268 + Qt2.v_01 * dt3.v_112 + Qt3.v_012 * dt2.v_11 + Qt3.v_011 * dt2.v_12
269 + Qt4.v_0112 * dt1.v_1 + Qt3.v_011 * dt2.v_12 + Qt4.v_0112 * dt1.v_1
270 + Qt4.v_0111 * dt1.v_2 + Qt5.v_01112,
271 dt1.v_0 * Qn4.v_1122 + dt1.v_1 * (Qn4.v_0122 - dt1.v_0 * Qn3.v_122)
272 + Qt2.v_01 * dt3.v_122 + Qt3.v_012 * dt2.v_12 + Qt3.v_011 * dt2.v_22
273 + Qt4.v_0112 * dt1.v_2 + Qt3.v_012 * dt2.v_12 + Qt4.v_0122 * dt1.v_1
274 + Qt4.v_0112 * dt1.v_2 + Qt5.v_01122,
275 dt1.v_0 * Qn4.v_1222 + dt1.v_1 * (Qn4.v_0222 - dt1.v_0 * Qn3.v_222)
276 + Qt2.v_01 * dt3.v_222 + Qt3.v_012 * dt2.v_22 + Qt3.v_012 * dt2.v_22
277 + Qt4.v_0122 * dt1.v_2 + Qt3.v_012 * dt2.v_22 + Qt4.v_0122 * dt1.v_2
278 + Qt4.v_0122 * dt1.v_2 + Qt5.v_01222,
279 dt1.v_0 * Qn4.v_2222 + dt1.v_2 * (Qn4.v_0222 - dt1.v_0 * Qn3.v_222)
280 + Qt2.v_02 * dt3.v_222 + Qt3.v_022 * dt2.v_22 + Qt3.v_022 * dt2.v_22
281 + Qt4.v_0222 * dt1.v_2 + Qt3.v_022 * dt2.v_22 + Qt4.v_0222 * dt1.v_2
282 + Qt4.v_0222 * dt1.v_2 + Qt5.v_02222,
283 dt1.v_1 * Qn4.v_1111 + dt1.v_1 * (Qn4.v_1111 - dt1.v_1 * Qn3.v_111)
284 + Qt2.v_11 * dt3.v_111 + Qt3.v_111 * dt2.v_11 + Qt3.v_111 * dt2.v_11
285 + Qt4.v_1111 * dt1.v_1 + Qt3.v_111 * dt2.v_11 + Qt4.v_1111 * dt1.v_1
286 + Qt4.v_1111 * dt1.v_1 + Qt5.v_11111,
287 dt1.v_1 * Qn4.v_1112 + dt1.v_1 * (Qn4.v_1112 - dt1.v_1 * Qn3.v_112)
288 + Qt2.v_11 * dt3.v_112 + Qt3.v_112 * dt2.v_11 + Qt3.v_111 * dt2.v_12
289 + Qt4.v_1112 * dt1.v_1 + Qt3.v_111 * dt2.v_12 + Qt4.v_1112 * dt1.v_1
290 + Qt4.v_1111 * dt1.v_2 + Qt5.v_11112,
291 dt1.v_1 * Qn4.v_1122 + dt1.v_1 * (Qn4.v_1122 - dt1.v_1 * Qn3.v_122)
292 + Qt2.v_11 * dt3.v_122 + Qt3.v_112 * dt2.v_12 + Qt3.v_111 * dt2.v_22
293 + Qt4.v_1112 * dt1.v_2 + Qt3.v_112 * dt2.v_12 + Qt4.v_1122 * dt1.v_1
294 + Qt4.v_1112 * dt1.v_2 + Qt5.v_11122,
295 dt1.v_1 * Qn4.v_1222 + dt1.v_1 * (Qn4.v_1222 - dt1.v_1 * Qn3.v_222)
296 + Qt2.v_11 * dt3.v_222 + Qt3.v_112 * dt2.v_22 + Qt3.v_112 * dt2.v_22
297 + Qt4.v_1122 * dt1.v_2 + Qt3.v_112 * dt2.v_22 + Qt4.v_1122 * dt1.v_2
298 + Qt4.v_1122 * dt1.v_2 + Qt5.v_11222,
299 dt1.v_1 * Qn4.v_2222 + dt1.v_2 * (Qn4.v_1222 - dt1.v_1 * Qn3.v_222)
300 + Qt2.v_12 * dt3.v_222 + Qt3.v_122 * dt2.v_22 + Qt3.v_122 * dt2.v_22
301 + Qt4.v_1222 * dt1.v_2 + Qt3.v_122 * dt2.v_22 + Qt4.v_1222 * dt1.v_2
302 + Qt4.v_1222 * dt1.v_2 + Qt5.v_12222,
303 dt1.v_2 * Qn4.v_2222 + dt1.v_2 * (Qn4.v_2222 - dt1.v_2 * Qn3.v_222)
304 + Qt2.v_22 * dt3.v_222 + Qt3.v_222 * dt2.v_22 + Qt3.v_222 * dt2.v_22
305 + Qt4.v_2222 * dt1.v_2 + Qt3.v_222 * dt2.v_22 + Qt4.v_2222 * dt1.v_2
306 + Qt4.v_2222 * dt1.v_2 + Qt5.v_22222};
307
308 // auto Qn5 = SymTensor3d_5<T>{
309 // Q.t5.v_00000 + 5*Q.t4.v_0000*d.t1.v_0 + 10*Q.t3.v_000*d.t2.v_00 +
310 // 10*Q.t2.v_00*d.t3.v_000
311 // + 5*Q.t1.v_0*d.t4.v_0000 + Q.t0*d.t5.v_00000, Q.t5.v_00001 +
312 // 4*Q.t4.v_0001*d.t1.v_0 + Q.t4.v_0000*d.t1.v_1 + 6*Q.t3.v_001*d.t2.v_00 +
313 // 4*Q.t3.v_000*d.t2.v_01 + 4*Q.t2.v_01*d.t3.v_000 + 6*Q.t2.v_00*d.t3.v_001 +
314 // Q.t1.v_1*d.t4.v_0000 + 4*Q.t1.v_0*d.t4.v_0001 + Q.t0*d.t5.v_00001, Q.t5.v_00002 +
315 // 4*Q.t4.v_0002*d.t1.v_0 + Q.t4.v_0000*d.t1.v_2 + 6*Q.t3.v_002*d.t2.v_00 +
316 // 4*Q.t3.v_000*d.t2.v_02 + 4*Q.t2.v_02*d.t3.v_000 + 6*Q.t2.v_00*d.t3.v_002 +
317 // Q.t1.v_2*d.t4.v_0000 + 4*Q.t1.v_0*d.t4.v_0002 + Q.t0*d.t5.v_00002, Q.t5.v_00011 +
318 // 3*Q.t4.v_0011*d.t1.v_0 + 2*Q.t4.v_0001*d.t1.v_1 + 3*Q.t3.v_011*d.t2.v_00 +
319 // 6*Q.t3.v_001*d.t2.v_01 + Q.t3.v_000*d.t2.v_11 + Q.t2.v_11*d.t3.v_000 +
320 // 6*Q.t2.v_01*d.t3.v_001 + 3*Q.t2.v_00*d.t3.v_011 + 2*Q.t1.v_1*d.t4.v_0001 +
321 // 3*Q.t1.v_0*d.t4.v_0011 + Q.t0*d.t5.v_00011, Q.t5.v_00012 + 3*Q.t4.v_0012*d.t1.v_0
322 // + Q.t4.v_0002*d.t1.v_1 + Q.t4.v_0001*d.t1.v_2 + 3*Q.t3.v_012*d.t2.v_00 +
323 // 3*Q.t3.v_002*d.t2.v_01 + 3*Q.t3.v_001*d.t2.v_02 + Q.t3.v_000*d.t2.v_12 +
324 // Q.t2.v_12*d.t3.v_000 + 3*Q.t2.v_02*d.t3.v_001 + 3*Q.t2.v_01*d.t3.v_002 +
325 // 3*Q.t2.v_00*d.t3.v_012 + Q.t1.v_2*d.t4.v_0001 + Q.t1.v_1*d.t4.v_0002 +
326 // 3*Q.t1.v_0*d.t4.v_0012 + Q.t0*d.t5.v_00012, Q.t5.v_00022 + 3*Q.t4.v_0022*d.t1.v_0
327 // + 2*Q.t4.v_0002*d.t1.v_2 + 3*Q.t3.v_022*d.t2.v_00 + 6*Q.t3.v_002*d.t2.v_02 +
328 // Q.t3.v_000*d.t2.v_22 + Q.t2.v_22*d.t3.v_000 + 6*Q.t2.v_02*d.t3.v_002 +
329 // 3*Q.t2.v_00*d.t3.v_022 + 2*Q.t1.v_2*d.t4.v_0002 + 3*Q.t1.v_0*d.t4.v_0022 +
330 // Q.t0*d.t5.v_00022, Q.t5.v_00111 + 2*Q.t4.v_0111*d.t1.v_0 + 3*Q.t4.v_0011*d.t1.v_1
331 // + Q.t3.v_111*d.t2.v_00 + 6*Q.t3.v_011*d.t2.v_01 + 3*Q.t3.v_001*d.t2.v_11 +
332 // 3*Q.t2.v_11*d.t3.v_001 + 6*Q.t2.v_01*d.t3.v_011 + Q.t2.v_00*d.t3.v_111 +
333 // 3*Q.t1.v_1*d.t4.v_0011 + 2*Q.t1.v_0*d.t4.v_0111 + Q.t0*d.t5.v_00111, Q.t5.v_00112
334 // + 2*Q.t4.v_0112*d.t1.v_0 + 2*Q.t4.v_0012*d.t1.v_1 + Q.t4.v_0011*d.t1.v_2 +
335 // Q.t3.v_112*d.t2.v_00 + 4*Q.t3.v_012*d.t2.v_01 + 2*Q.t3.v_011*d.t2.v_02 +
336 // Q.t3.v_002*d.t2.v_11 + 2*Q.t3.v_001*d.t2.v_12 + 2*Q.t2.v_12*d.t3.v_001 +
337 // Q.t2.v_11*d.t3.v_002 + 2*Q.t2.v_02*d.t3.v_011 + 4*Q.t2.v_01*d.t3.v_012 +
338 // Q.t2.v_00*d.t3.v_112 + Q.t1.v_2*d.t4.v_0011 + 2*Q.t1.v_1*d.t4.v_0012 +
339 // 2*Q.t1.v_0*d.t4.v_0112 + Q.t0*d.t5.v_00112, Q.t5.v_00122 + 2*Q.t4.v_0122*d.t1.v_0
340 // + Q.t4.v_0022*d.t1.v_1 + 2*Q.t4.v_0012*d.t1.v_2 + Q.t3.v_122*d.t2.v_00 +
341 // 2*Q.t3.v_022*d.t2.v_01 + 4*Q.t3.v_012*d.t2.v_02 + 2*Q.t3.v_002*d.t2.v_12 +
342 // Q.t3.v_001*d.t2.v_22 + Q.t2.v_22*d.t3.v_001 + 2*Q.t2.v_12*d.t3.v_002 +
343 // 4*Q.t2.v_02*d.t3.v_012 + 2*Q.t2.v_01*d.t3.v_022 + Q.t2.v_00*d.t3.v_122 +
344 // 2*Q.t1.v_2*d.t4.v_0012 + Q.t1.v_1*d.t4.v_0022 + 2*Q.t1.v_0*d.t4.v_0122 +
345 // Q.t0*d.t5.v_00122, Q.t5.v_00222 + 2*Q.t4.v_0222*d.t1.v_0 + 3*Q.t4.v_0022*d.t1.v_2
346 // + Q.t3.v_222*d.t2.v_00 + 6*Q.t3.v_022*d.t2.v_02 + 3*Q.t3.v_002*d.t2.v_22 +
347 // 3*Q.t2.v_22*d.t3.v_002 + 6*Q.t2.v_02*d.t3.v_022 + Q.t2.v_00*d.t3.v_222 +
348 // 3*Q.t1.v_2*d.t4.v_0022 + 2*Q.t1.v_0*d.t4.v_0222 + Q.t0*d.t5.v_00222, Q.t5.v_01111
349 // + Q.t4.v_1111*d.t1.v_0 + 4*Q.t4.v_0111*d.t1.v_1 + 4*Q.t3.v_111*d.t2.v_01 +
350 // 6*Q.t3.v_011*d.t2.v_11 + 6*Q.t2.v_11*d.t3.v_011 + 4*Q.t2.v_01*d.t3.v_111 +
351 // 4*Q.t1.v_1*d.t4.v_0111 + Q.t1.v_0*d.t4.v_1111 + Q.t0*d.t5.v_01111, Q.t5.v_01112 +
352 // Q.t4.v_1112*d.t1.v_0 + 3*Q.t4.v_0112*d.t1.v_1 + Q.t4.v_0111*d.t1.v_2 +
353 // 3*Q.t3.v_112*d.t2.v_01 + Q.t3.v_111*d.t2.v_02 + 3*Q.t3.v_012*d.t2.v_11 +
354 // 3*Q.t3.v_011*d.t2.v_12 + 3*Q.t2.v_12*d.t3.v_011 + 3*Q.t2.v_11*d.t3.v_012 +
355 // Q.t2.v_02*d.t3.v_111 + 3*Q.t2.v_01*d.t3.v_112 + Q.t1.v_2*d.t4.v_0111 +
356 // 3*Q.t1.v_1*d.t4.v_0112 + Q.t1.v_0*d.t4.v_1112 + Q.t0*d.t5.v_01112, Q.t5.v_01122 +
357 // Q.t4.v_1122*d.t1.v_0 + 2*Q.t4.v_0122*d.t1.v_1 + 2*Q.t4.v_0112*d.t1.v_2 +
358 // 2*Q.t3.v_122*d.t2.v_01 + 2*Q.t3.v_112*d.t2.v_02 + Q.t3.v_022*d.t2.v_11 +
359 // 4*Q.t3.v_012*d.t2.v_12 + Q.t3.v_011*d.t2.v_22 + Q.t2.v_22*d.t3.v_011 +
360 // 4*Q.t2.v_12*d.t3.v_012 + Q.t2.v_11*d.t3.v_022 + 2*Q.t2.v_02*d.t3.v_112 +
361 // 2*Q.t2.v_01*d.t3.v_122 + 2*Q.t1.v_2*d.t4.v_0112 + 2*Q.t1.v_1*d.t4.v_0122 +
362 // Q.t1.v_0*d.t4.v_1122 + Q.t0*d.t5.v_01122, Q.t5.v_01222 + Q.t4.v_1222*d.t1.v_0 +
363 // Q.t4.v_0222*d.t1.v_1 + 3*Q.t4.v_0122*d.t1.v_2 + Q.t3.v_222*d.t2.v_01 +
364 // 3*Q.t3.v_122*d.t2.v_02 + 3*Q.t3.v_022*d.t2.v_12 + 3*Q.t3.v_012*d.t2.v_22 +
365 // 3*Q.t2.v_22*d.t3.v_012 + 3*Q.t2.v_12*d.t3.v_022 + 3*Q.t2.v_02*d.t3.v_122 +
366 // Q.t2.v_01*d.t3.v_222 + 3*Q.t1.v_2*d.t4.v_0122 + Q.t1.v_1*d.t4.v_0222 +
367 // Q.t1.v_0*d.t4.v_1222 + Q.t0*d.t5.v_01222, Q.t5.v_02222 + Q.t4.v_2222*d.t1.v_0 +
368 // 4*Q.t4.v_0222*d.t1.v_2 + 4*Q.t3.v_222*d.t2.v_02 + 6*Q.t3.v_022*d.t2.v_22 +
369 // 6*Q.t2.v_22*d.t3.v_022 + 4*Q.t2.v_02*d.t3.v_222 + 4*Q.t1.v_2*d.t4.v_0222 +
370 // Q.t1.v_0*d.t4.v_2222 + Q.t0*d.t5.v_02222, Q.t5.v_11111 + 5*Q.t4.v_1111*d.t1.v_1 +
371 // 10*Q.t3.v_111*d.t2.v_11 + 10*Q.t2.v_11*d.t3.v_111 + 5*Q.t1.v_1*d.t4.v_1111 +
372 // Q.t0*d.t5.v_11111, Q.t5.v_11112 + 4*Q.t4.v_1112*d.t1.v_1 + Q.t4.v_1111*d.t1.v_2 +
373 // 6*Q.t3.v_112*d.t2.v_11 + 4*Q.t3.v_111*d.t2.v_12 + 4*Q.t2.v_12*d.t3.v_111 +
374 // 6*Q.t2.v_11*d.t3.v_112 + Q.t1.v_2*d.t4.v_1111 + 4*Q.t1.v_1*d.t4.v_1112 +
375 // Q.t0*d.t5.v_11112, Q.t5.v_11122 + 3*Q.t4.v_1122*d.t1.v_1 + 2*Q.t4.v_1112*d.t1.v_2
376 // + 3*Q.t3.v_122*d.t2.v_11 + 6*Q.t3.v_112*d.t2.v_12 + Q.t3.v_111*d.t2.v_22 +
377 // Q.t2.v_22*d.t3.v_111 + 6*Q.t2.v_12*d.t3.v_112 + 3*Q.t2.v_11*d.t3.v_122 +
378 // 2*Q.t1.v_2*d.t4.v_1112 + 3*Q.t1.v_1*d.t4.v_1122 + Q.t0*d.t5.v_11122, Q.t5.v_11222
379 // + 2*Q.t4.v_1222*d.t1.v_1 + 3*Q.t4.v_1122*d.t1.v_2 + Q.t3.v_222*d.t2.v_11 +
380 // 6*Q.t3.v_122*d.t2.v_12 + 3*Q.t3.v_112*d.t2.v_22 + 3*Q.t2.v_22*d.t3.v_112 +
381 // 6*Q.t2.v_12*d.t3.v_122 + Q.t2.v_11*d.t3.v_222 + 3*Q.t1.v_2*d.t4.v_1122 +
382 // 2*Q.t1.v_1*d.t4.v_1222 + Q.t0*d.t5.v_11222, Q.t5.v_12222 + Q.t4.v_2222*d.t1.v_1 +
383 // 4*Q.t4.v_1222*d.t1.v_2 + 4*Q.t3.v_222*d.t2.v_12 + 6*Q.t3.v_122*d.t2.v_22 +
384 // 6*Q.t2.v_22*d.t3.v_122 + 4*Q.t2.v_12*d.t3.v_222 + 4*Q.t1.v_2*d.t4.v_1222 +
385 // Q.t1.v_1*d.t4.v_2222 + Q.t0*d.t5.v_12222, Q.t5.v_22222 + 5*Q.t4.v_2222*d.t1.v_2 +
386 // 10*Q.t3.v_222*d.t2.v_22 + 10*Q.t2.v_22*d.t3.v_222 + 5*Q.t1.v_2*d.t4.v_2222 +
387 // Q.t0*d.t5.v_22222
388 // };
389 }
390
391 } // namespace details
392
394 template<class T, u32 low_order, u32 high_order>
397 const sycl::vec<T, 3> &offset) {
398
399 using namespace shammath;
400 using namespace shamphys::details;
401
402 if constexpr (low_order == 0 && high_order == 5) {
404
405 auto Qn1 = offset_multipole_1(Q.t0, Q.t1, d.t1);
406
407 auto Qn2 = offset_multipole_2(Qn1, Q.t1, Q.t2, d.t1);
408
409 auto Qn3 = offset_multipole_3(Qn1, Qn2, Q.t2, Q.t3, d.t1);
410
411 auto Qn4 = offset_multipole_4(Qn2, Qn3, Q.t2, Q.t3, Q.t4, d.t1, d.t2);
412
413 auto Qn5 = offset_multipole_5(Qn3, Qn4, Q.t2, Q.t3, Q.t4, Q.t5, d.t1, d.t2, d.t3);
414
415 return {Q.t0, Qn1, Qn2, Qn3, Qn4, Qn5};
416 } else if constexpr (low_order == 0 && high_order == 4) {
417
419
420 auto Qn1 = offset_multipole_1(Q.t0, Q.t1, d.t1);
421
422 auto Qn2 = offset_multipole_2(Qn1, Q.t1, Q.t2, d.t1);
423
424 auto Qn3 = offset_multipole_3(Qn1, Qn2, Q.t2, Q.t3, d.t1);
425
426 auto Qn4 = offset_multipole_4(Qn2, Qn3, Q.t2, Q.t3, Q.t4, d.t1, d.t2);
427
428 return {Q.t0, Qn1, Qn2, Qn3, Qn4};
429 } else if constexpr (low_order == 0 && high_order == 3) {
430
432
433 auto Qn1 = offset_multipole_1(Q.t0, Q.t1, d.t1);
434
435 auto Qn2 = offset_multipole_2(Qn1, Q.t1, Q.t2, d.t1);
436
437 auto Qn3 = offset_multipole_3(Qn1, Qn2, Q.t2, Q.t3, d.t1);
438
439 return {Q.t0, Qn1, Qn2, Qn3};
440 } else if constexpr (low_order == 0 && high_order == 2) {
441
443
444 auto Qn1 = offset_multipole_1(Q.t0, Q.t1, d.t1);
445
446 auto Qn2 = offset_multipole_2(Qn1, Q.t1, Q.t2, d.t1);
447
448 return {Q.t0, Qn1, Qn2};
449 } else if constexpr (low_order == 0 && high_order == 1) {
450
452
453 auto Qn1 = offset_multipole_1(Q.t0, Q.t1, d.t1);
454
455 return {Q.t0, Qn1};
456 } else if constexpr (low_order == 0 && high_order == 0) {
457 return {Q.t0};
458 } else {
459 static_assert(shambase::always_false_v<T>, "This combination of orders is not valid");
460 }
461 }
462
463 template<class T, u32 low_order, u32 high_order>
466 const sycl::vec<T, 3> &from,
467 const sycl::vec<T, 3> &to) {
468 return offset_multipole_delta(Q_old, from - to);
469 }
470
471} // namespace shamphys
Namespace for internal details of the logs module.
void throw_with_loc(std::string message, SourceLocation loc=SourceLocation{})
Throw an exception and append the source location to it.
namespace for math utility
Definition AABB.hpp:26
shammath::SymTensorCollection< T, low_order, high_order > offset_multipole_delta(const shammath::SymTensorCollection< T, low_order, high_order > &Q, const sycl::vec< T, 3 > &offset)
utility to offset a multipole, see PHD
Traits for C++ types.