AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
Arithmetic.hpp
1#pragma once
2
3#include <algorithm>
4
5#include <amdis/common/Math.hpp>
6#include <amdis/common/Concepts.hpp>
7#include <amdis/operations/Basic.hpp>
8#include <amdis/operations/Composer.hpp>
9
10namespace AMDiS
11{
12 namespace Operation
13 {
19 struct Plus
20 {
21 template <class... Ts>
22 constexpr auto operator()(Ts const&... ts) const
23 {
24 return Math::sum(ts...);
25 }
26 };
27
28#ifndef DOXYGEN
29 // [0] + g => g
30 template <class G>
32 : ComposerBuilder<Id, G> {};
33
34 // g + [0] => g
35 template <class G>
37 : ComposerBuilder<Id, G> {};
38
39 // [0] + [0] => [0]
40 template <>
42 : ComposerBuilder<Id, Zero> {};
43#endif
44
45 template <class... Int>
46 constexpr int order(Plus const&, Int... orders)
47 {
48 return Math::max(int(orders)...);
49 }
50
51 template <std::size_t I>
52 constexpr auto partial(Plus const&, index_t<I>)
53 {
54 static_assert((I < 2), "Derivatives of `Plus` only defined for the binary case.");
55 return One{};
56 }
57
58 // -------------------------------------------------------------------------
59
61 struct Minus
62 {
63 template <class T, class S>
64 constexpr auto operator()(T const& lhs, S const& rhs) const
65 {
66 return lhs - rhs;
67 }
68
69 friend constexpr int order(Minus const&, int lhs, int rhs)
70 {
71 return Math::max(lhs, rhs);
72 }
73
74 friend constexpr auto partial(Minus const&, index_t<0>)
75 {
76 return One{};
77 }
78
79 friend constexpr auto partial(Minus const&, index_t<1>)
80 {
81 return StaticConstant<int,-1>{};
82 }
83 };
84
85#ifndef DOXYGEN
86 // g - [0] => g
87 template <class G>
89 : ComposerBuilder<Id, G> {};
90
91 // [0] - [0] => [0]
92 template <>
94 : ComposerBuilder<Id, Zero> {};
95#endif
96
97 // -------------------------------------------------------------------------
98
101 {
102 template <class... Ts>
103 constexpr auto operator()(Ts const&... ts) const
104 {
105 return (ts * ...);
106 }
107 };
108
109
110#ifndef DOXYGEN
111 // [0] * g => [0]
112 template <class G>
114 : ComposerBuilder<Id, Zero> {};
115
116 // g * [0] => [0]
117 template <class G>
119 : ComposerBuilder<Id, Zero> {};
120
121 // [0] * [0] => [0]
122 template <>
124 : ComposerBuilder<Id, Zero> {};
125#endif
126
127
128 template <class... Int>
129 constexpr int order(Multiplies const&, Int... orders)
130 {
131 return Math::sum(int(orders)...);
132 }
133
134 // only for binary *
135 // d_0 (x * y) = y, d_1 (x * y) = x
136 template <std::size_t I>
137 constexpr auto partial(Multiplies const&, index_t<I>)
138 {
139 static_assert((I < 2), "Derivatives of `Multiplies` only defined for the binary case.");
140 return Arg<1-I>{};
141 }
142
143 // -------------------------------------------------------------------------
144
146 template <class Factor>
147 struct Scale
148 {
149 Factor const factor_;
150
151 template <class T>
152 constexpr auto operator()(T const& value) const
153 {
154 return factor_ * value;
155 }
156
157 friend constexpr int order(Scale const&, int d)
158 {
159 return d;
160 }
161
162 // d_x f*x = f
163 friend constexpr auto partial(Scale const& s, index_t<0>)
164 {
165 return Constant<Factor>{s.factor_};
166 }
167 };
168
169#ifndef DOXYGEN
170 // f*(g*A) = (f*g)*A
171 template <class F, class G>
173 {
175
176 template <class F_, class G_>
177 static constexpr type build(F_&& f, G_&& g)
178 {
179 return type{f.factor_ * g.factor_};
180 }
181 };
182#endif
183
184 // -------------------------------------------------------------------------
185
187 struct Divides
188 {
189 template <class T, class S>
190 constexpr auto operator()(T const& lhs, S const& rhs) const
191 {
192 return lhs / rhs;
193 }
194
195 // d_0 f(x,y) = 1 / y
196 friend constexpr auto partial(Divides const&, index_t<0>)
197 {
198 return compose(Divides{}, One{}, Arg<1>{});
199 }
200
201 // d_1 f(x,y) = (y - x)/y^2
202 friend constexpr auto partial(Divides const&, index_t<1>);
203 };
204
205 // -------------------------------------------------------------------------
206
208 struct Negate
209 {
210 template <class T>
211 constexpr auto operator()(T const& x) const
212 {
213 return -x;
214 }
215
216 friend constexpr int order(Negate const&, int d)
217 {
218 return d;
219 }
220
221 friend constexpr auto partial(Negate const&, index_t<0>)
222 {
223 return StaticConstant<int,-1>{};
224 }
225 };
226
227#ifndef DOXYGEN
228 // g + -h => g - h
229 template <class G>
231 : ComposerBuilder<Minus, G, Id> {};
232
233 // [0] - g => -g
234 template <class G>
236 : ComposerBuilder<Id, Negate> {};
237
238 // -(g - h) => (h - g)
239 template <>
241 : ComposerBuilder<Minus, Arg<1>, Arg<0>> {};
242#endif
243
244 // -------------------------------------------------------------------------
245
246 // forward declaration
247 template <int p, bool positive>
248 struct PowImpl;
249
250 template <int p>
251 struct PowType
252 {
253 using type = PowImpl<p, (p>0)>;
254 };
255
256 template <> struct PowType<1> { using type = Id; };
257 template <> struct PowType<0> { using type = One; };
258
259 template <int p>
260 using Pow = typename PowType<p>::type;
261
262 using Sqr = Pow<2>;
263
265 template <int p>
266 struct PowImpl<p, true>
267 {
268 template <class T>
269 constexpr auto operator()(T const& x) const
270 {
271 return Math::pow<p>(x);
272 }
273
274 friend constexpr int order(PowImpl const&, int d)
275 {
276 return p*d;
277 }
278
279 friend constexpr auto partial(PowImpl const&, index_t<0>)
280 {
281 return compose(Multiplies{}, StaticConstant<int,p>{}, Pow<p-1>{});
282 }
283 };
284
285 template <int p>
286 struct PowImpl<p, false>
287 : public Composer<Divides, One, Pow<-p>>
288 {
289 constexpr PowImpl()
290 : Composer<Divides, One, Pow<-p>>{}
291 {}
292
293 template <class T>
294 constexpr auto operator()(T const& x) const
295 {
296 return T(1) / Math::pow<-p>(x);
297 }
298 };
299
300
301#ifndef DOXYGEN
302 // (x^p1)^p2 => x^(p1*p2)
303 template <int p1, bool pos1, int p2, bool pos2>
304 struct ComposerBuilder<PowImpl<p1,pos1>, PowImpl<p2,pos2>>
305 : ComposerBuilder<Id, Pow<p1*p2>> {};
306
307 // x^p * y^p => (x*y)^p
308 template <int p, bool pos>
310 : ComposerBuilder<Pow<p>, Multiplies> {};
311#endif
312
313 // d_1 f(x,y) = (y - x)/y^2
314 inline constexpr auto partial(Divides const&, index_t<1>)
315 {
316 return compose(Divides{}, compose(Minus{}, Arg<1>{}, Arg<0>{}),
317 compose(Pow<2>{}, Arg<1>{}));
318 }
319
321 struct Power
322 {
323 constexpr Power(double p)
324 : p_(p)
325 {}
326
327 template <class T>
328 auto operator()(T const& x) const
329 {
330 return std::pow(x, p_);
331 }
332
333 friend constexpr auto partial(Power const& P, index_t<0>)
334 {
335 return compose(Multiplies{}, Constant<double>{P.p_}, Power{P.p_-1});
336 }
337
338 double p_;
339 };
340
343 } // end namespace Operation
344} // end namespace AMDiS
Definition: Basic.hpp:165
Definition: Composer.hpp:52
Composition of Functors.
Definition: Composer.hpp:30
Functor representing a constant value.
Definition: Basic.hpp:88
Functor that represents A/B.
Definition: Arithmetic.hpp:188
(Unary-)Functor representing the identity
Definition: Basic.hpp:65
Functor that represents A-B.
Definition: Arithmetic.hpp:62
Functor that represents A*B.
Definition: Arithmetic.hpp:101
Functor that represents A-B.
Definition: Arithmetic.hpp:209
Functor that represents A+B.
Definition: Arithmetic.hpp:20
Definition: Arithmetic.hpp:248
Definition: Arithmetic.hpp:252
Functor that represents x^p,.
Definition: Arithmetic.hpp:322
Functor that represents f*A.
Definition: Arithmetic.hpp:148
Functor representing a static constant value.
Definition: Basic.hpp:38