DUNE-ACFEM (unstable)

fractionconstant.hh
1#ifndef __DUNE_ACFEM_COMMON_FRACTIONCONSTANT_HH__
2#define __DUNE_ACFEM_COMMON_FRACTIONCONSTANT_HH__
3
4#include <string>
5#include <numeric>
6#include <ratio>
7
8#include <dune/acfem/common/types.hh>
9namespace Dune {
10
11 namespace ACFem {
12
13 namespace TypedValue {
14
25 template<class Int, Int N, Int D>
26 struct FractionConstant
27 {
28 using IntegerType = Int;
29 static constexpr IntegerType numerator = N;
30 static constexpr IntegerType denominator = D;
31
32 static_assert(D > 0 && (std::gcd(N, D) == 1 || std::gcd(N, D) == -1),
33 "GCD should be 1 here ... oops");
34
35 constexpr FractionConstant()
36 {}
37
38 template<class Int2, Int2 N2, Int2 D2>
39 constexpr FractionConstant(Constant<Int2, N2>, Constant<Int2, D2>)
40 {
41 static_assert((N == (Int)N2 / std::gcd((Int)N2, (Int)D2))
42 &&
43 (D == (Int)D2 / std::gcd((Int)N2, (Int)D2)),
44 "Assignment from integral constants only possible if their quotient matches the value of the fraction");
45 }
46
47 template<class Int2, Int2 N2>
48 constexpr FractionConstant(Constant<Int2, N2>)
49 : FractionConstant(Constant<Int, (Int)N2>{}, Constant<Int, (Int)1>{})
50 {}
51
52 constexpr FractionConstant(const FractionConstant&)
53 {}
54
55 constexpr FractionConstant(FractionConstant&&)
56 {}
57
58 constexpr FractionConstant& operator=(const FractionConstant&)
59 {
60 return *this;
61 }
62
63 template<class T, std::enable_if_t<std::numeric_limits<std::decay_t<T> >::is_specialized, int> = 0>
64 constexpr operator T() const
65 {
66 return static_cast<T>(N) / static_cast<T>(D);
67 }
68
69 };
70
72 template<class T, bool RequireDecay = false, class SFINAE = void>
73 struct IsFractionConstant
74 : FalseType
75 {};
76
77 template<class T>
78 struct IsFractionConstant<T, false, std::enable_if_t<!IsDecay<T>::value> >
79 : IsFractionConstant<std::decay_t<T> >
80 {};
81
82 template<class I, I N, I D, bool RequireDecay>
83 struct IsFractionConstant<FractionConstant<I, N, D>, RequireDecay>
84 : TrueType
85 {};
86
87 template<class T>
88 constexpr T sign(T val) {
89 return (T(0) < val) - (val < T(0));
90 }
91
92 template<class T>
93 constexpr T abs(T val) {
94 return sign(val) * val;
95 }
96
100 template<class I, I N, I D = 1>
101 using Fraction =
102 FractionConstant<I,
103 sign(D) * N / std::gcd(N, D),
104 abs(D) / std::gcd(N, D)>;
105
107 template<std::ptrdiff_t N, std::ptrdiff_t D = 1>
108 using IntFraction = Fraction<std::ptrdiff_t, N, D>;
109
123 template<std::ptrdiff_t N, std::ptrdiff_t D = 1>
124 constexpr auto intFraction(IntConstant<N> = IntConstant<N>{}, IntConstant<D> = IntConstant<D>{})
125 {
126 return IntFraction<N, D>{};
127 }
128
129 template<std::size_t N, std::size_t D = 1>
130 constexpr auto intFraction(IndexConstant<N>, IndexConstant<D> = IndexConstant<D>{})
131 {
132 return IntFraction<(std::ptrdiff_t)N, (std::ptrdiff_t)D>{};
133 }
134
135 template<class T>
136 constexpr std::ptrdiff_t numerator(T&& t = T{})
137 {
138 static_assert(IsFractionConstant<T>::value, "T must fulfil IsFractionConstant<T>");
139 return std::decay_t<T>::numerator;
140 }
141
142 template<class T>
143 constexpr std::ptrdiff_t denominator(T&& t = T{})
144 {
145 static_assert(IsFractionConstant<T>::value, "T must fulfil IsFractionConstant<T>");
146 return std::decay_t<T>::denominator;
147 }
148
150 template<std::ptrdiff_t S>
151 using Sign = IntFraction<(S > 0 ? 1 : (S < 0 ? -1 : 0))>;
152
154 template<std::ptrdiff_t S>
155 constexpr auto sign(IntConstant<S> = IntConstant<S>{})
156 {
157 return Sign<S>{};
158 }
159
161 template<class T>
162 using SignOf = Sign<T{}>;
163
165 template<class T>
166 constexpr auto signOf(T&&)
167 {
168 return SignOf<T>{};
169 }
170
172 template<class T, class SFINAE = void>
173 struct IsSign;
174
176 template<class T>
177 struct IsSign<T, std::enable_if_t<!IsFractionConstant<T>::value> >
178 : FalseType
179 {};
180
184 template<class T>
185 struct IsSign<T, std::enable_if_t<IsFractionConstant<T>::value> >
186 : BoolConstant<(numerator<T>()*numerator<T>() <= 1 && denominator<T>() == 1)>
187 {};
188
190 template<class T1, class T2>
191 using FieldPromotionType = std::decay_t<decltype(std::declval<T1>() * std::declval<T2>())>;
192
193 template<class I, I N, I D>
194 constexpr auto operator-(FractionConstant<I, N, D>)
195 {
196 return Fraction<I, -N, D>{};
197 }
198
199 template<class I, I N, I D>
200 constexpr auto reciprocal(FractionConstant<I, N, D>)
201 {
202 static_assert(N != 0, "Reciprocal of zero attempted.");
203 return Fraction<I, D, N>{};
204 }
205
206 template<class I1, I1 N1, I1 D1,
207 class I2, I2 N2, I2 D2>
208 constexpr auto operator*(FractionConstant<I1, N1, D1>, FractionConstant<I2, N2, D2>)
209 {
210 // normalize
211 using T1 = Fraction<I1, N1, D1>;
212 using T2 = Fraction<I2, N2, D2>;
213
214 // cross normalize to reduce the chance of overflow
215 using C1 = Fraction<FieldPromotionType<I1, I2>, T1::numerator, T2::denominator>;
216 using C2 = Fraction<FieldPromotionType<I1, I2>, T2::numerator, T1::denominator>;
217
218 // calculate
219 return Fraction<FieldPromotionType<I1, I2>,
220 C1::numerator * C2::numerator,
221 C1::denominator * C2::denominator>{};
222 }
223
224 template<class I1, I1 N1, I1 D1,
225 class I2, I2 N2, I2 D2>
226 constexpr auto operator+(FractionConstant<I1, N1, D1>, FractionConstant<I2, N2, D2>)
227 {
228 // normalize
229 using Int = FieldPromotionType<I1, I2>;
230 using T1 = Fraction<Int, N1, D1>;
231 using T2 = Fraction<Int, N2, D2>;
232
233 // calculate
234
235 // reduce the chance of overflow
236 constexpr Int denomGcd = std::gcd(T1::denominator, T2::denominator);
237 constexpr Int denom1 = T1::denominator/denomGcd;
238 constexpr Int denom2 = T2::denominator/denomGcd;
239
240 return Fraction<Int,
241 T1::numerator*denom2 + T2::numerator*denom1,
242 denom1*denom2*denomGcd>{};
243 }
244
245 template<class I1, I1 N1, I1 D1,
246 class I2, I2 N2, I2 D2>
247 constexpr auto operator-(FractionConstant<I1, N1, D1>, FractionConstant<I2, N2, D2>)
248 {
249 // normalize
250 using Int = FieldPromotionType<I1, I2>;
251 using T1 = Fraction<Int, N1, D1>;
252 using T2 = Fraction<Int, N2, D2>;
253
254 // calculate
255
256 // reduce the chance of overflow
257 constexpr Int denomGcd = std::gcd(T1::denominator, T2::denominator);
258 constexpr Int denom1 = T1::denominator/denomGcd;
259 constexpr Int denom2 = T2::denominator/denomGcd;
260
261 return Fraction<Int,
262 T1::numerator*denom2 - T2::numerator*denom1,
263 denom1*denom2*denomGcd>{};
264 }
265
266 template<class I1, I1 N1, I1 D1,
267 class I2, I2 N2, I2 D2>
268 constexpr auto operator/(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
269 {
270 static_assert(N2 != 0, "Division by Zero.");
271
272 return f1 * reciprocal(f2);
273 }
274
275 template<class I1, I1 N1, class I2, I2 N2>
276 constexpr auto operator%(FractionConstant<I1, N1, 1>, FractionConstant<I2, N2, 1>)
277 {
278 return Fraction<FieldPromotionType<I1, I2>, (N1 % N2), 1>{};
279 }
280
281 template<class I1, I1 N1, class I2, I2 N2>
282 constexpr auto div(FractionConstant<I1, N1, 1>, FractionConstant<I2, N2, 1>)
283 {
284 return Fraction<FieldPromotionType<I1, I2>, N1 / N2, 1>{};
285 }
286
288 template<class I1, I1 N1, I1 D1, class I2>
289 constexpr auto pow(FractionConstant<I1, N1, D1>, FractionConstant<I2, 0, 1>)
290 {
291 return FractionConstant<I1, 1, 1>{};
292 }
293
295 template<class I1, I1 N1, I1 D1, class I2>
296 constexpr auto pow(FractionConstant<I1, N1, D1>, FractionConstant<I2, 1, 1>)
297 {
298 return FractionConstant<I1, N1, D1>{};
299 }
300
302 template<class I1, I1 N1, I1 D1, class I2, I2 N2, std::enable_if_t<(N2 > 1), int> = 0>
303 constexpr auto pow(FractionConstant<I1, N1, D1>, FractionConstant<I2, N2, 1>)
304 {
305 return FractionConstant<I1, N1, D1>{}*pow(FractionConstant<I1, N1, D1>{}, FractionConstant<I2, N2-1, 1>{});
306 }
307
308 template<class I, I N, I D>
309 constexpr auto sqr(FractionConstant<I, N, D>)
310 {
311 return pow(FractionConstant<I, N, D>{}, FractionConstant<I, 2, 1>{});
312 }
313
315 template<class I1, I1 N1, I1 D1, class I2, I2 N2, std::enable_if_t<(N2 < 0 && N1 != 0), int> = 0>
316 constexpr auto pow(FractionConstant<I1, N1, D1>, FractionConstant<I2, N2, 1>)
317 {
318 return pow(Fraction<I1, D1, N1>{}, FractionConstant<I2, N2, 1>{});
319 }
320
321 template<class I, I N, I D>
322 constexpr auto floor(FractionConstant<I, N, D>)
323 {
324 if constexpr (D == 1) {
325 return FractionConstant<I, N, 1>{};
326 } else if constexpr (N < 0) {
327 // -1.5 = -3/2 -> -2
328 // -0.5 = -1/2 -> -1
329 return FractionConstant<I, N/D-1, 1>{};
330 } else {
331 // 1.5 = 3/2 -> 1
332 return FractionConstant<I, N/D, 1>{};
333 }
334 }
335
336 template<class I, I N, I D>
337 constexpr auto ceil(FractionConstant<I, N, D>)
338 {
339 if constexpr (D == 1) {
340 return FractionConstant<I, N, 1>{};
341 } else {
342 // -1.5 = -3/2 -> (-3+1)/2 = -1
343 // 1.5 = 3/2 -> (3+1)/2 = 2
344 return FractionConstant<I, (N+1)/D, 1>{};
345 }
346 }
347
349 template<class I, I N, I D>
350 constexpr auto round(FractionConstant<I, N, D>)
351 {
352 // -9/14, (-9+14)/14 = 5/14 = 0
353 return FractionConstant<I, (N < 0 ? (N-D/2)/D : (N-D/2)/D) , 1>{};
354 }
355
356 template<class I, I N, I D>
357 constexpr auto abs(FractionConstant<I, N, D>)
358 {
359 if constexpr (N < 0) {
360 return FractionConstant<I, -N, D>{};
361 } else {
362 return FractionConstant<I, N, D>{};
363 }
364 }
365
366 template<class I1, I1 N1, I1 D1,
367 class I2, I2 N2, I2 D2>
368 constexpr bool operator==(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
369 {
370 // closure type
371 using Int = FieldPromotionType<I1, I2>;
372
373 // normalize
374 using T1 = Fraction<Int, N1, D1>;
375 using T2 = Fraction<Int, N2, D2>;
376
377 // reduce the chance of overflow
378 constexpr Int denomGcd = std::gcd(T1::denominator, T2::denominator);
379 constexpr Int denom1 = T1::denominator/denomGcd;
380 constexpr Int denom2 = T2::denominator/denomGcd;
381
382 return T1::numerator*denom2 == T2::numerator*denom1;
383 }
384
385 template<class I1, I1 N1, I1 D1,
386 class I2, I2 N2, I2 D2>
387 constexpr bool operator!=(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
388 {
389 return !(f1 == f2);
390 }
391
392 template<class I1, I1 N1, I1 D1,
393 class I2, I2 N2, I2 D2>
394 constexpr bool operator<(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
395 {
396 // closure type
397 using Int = FieldPromotionType<I1, I2>;
398
399 // normalize
400 using T1 = Fraction<Int, N1, D1>;
401 using T2 = Fraction<Int, N2, D2>;
402
403 // reduce the chance of overflow
404 constexpr Int denomGcd = std::gcd(T1::denominator, T2::denominator);
405 constexpr Int denom1 = T1::denominator/denomGcd;
406 constexpr Int denom2 = T2::denominator/denomGcd;
407
408 return T1::numerator*denom2 < T2::numerator*denom1;
409 }
410
411 template<class I1, I1 N1, I1 D1,
412 class I2, I2 N2, I2 D2>
413 constexpr bool operator<=(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
414 {
415 return f1 < f2 || f1 == f2;
416 }
417
418 template<class I1, I1 N1, I1 D1,
419 class I2, I2 N2, I2 D2>
420 constexpr bool operator>(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
421 {
422 return f2 < f1;
423 }
424
425 template<class I1, I1 N1, I1 D1,
426 class I2, I2 N2, I2 D2>
427 constexpr bool operator>=(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
428 {
429 return f2 <= f1;
430 }
431
432 using std::min;
433
434 template<class I1, I1 N1, I1 D1,
435 class I2, I2 N2, I2 D2,
436 std::enable_if_t<!std::is_same<FractionConstant<I1, N1, D1>,
437 FractionConstant<I2, N2, D2> >::value, int> = 0>
438 constexpr auto min(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
439 {
440 if constexpr (f1 <= f2) {
441 return f1;
442 } else {
443 return f2;
444 }
445 }
446
447 using std::max;
448
449 template<class I1, I1 N1, I1 D1,
450 class I2, I2 N2, I2 D2,
451 std::enable_if_t<!std::is_same<FractionConstant<I1, N1, D1>,
452 FractionConstant<I2, N2, D2> >::value, int> = 0>
453 constexpr auto max(FractionConstant<I1, N1, D1> f1, FractionConstant<I2, N2, D2> f2)
454 {
455 if constexpr (f1 >= f2) {
456 return f1;
457 } else {
458 return f2;
459 }
460 }
461
462 template<class I, I N, I D>
463 std::string toString(FractionConstant<I, N, D>)
464 {
465 // normalize
466 using T = Fraction<I, N, D>;
467 if constexpr (T::denominator == 1) {
468 return std::to_string(T::numerator) + "_c";
469 } else {
470 return "(" + std::to_string(T::numerator) + "/" + std::to_string(T::denominator) + ")_c";
471 }
472 }
473
475 template<class I, I N, I D>
476 std::ostream& operator<<(std::ostream& out, FractionConstant<I, N, D> frac)
477 {
478 return out << toString(frac);
479 }
480
481 } // NS TypedValue
482
483 template<class T>
484 using IsFractionConstant = TypedValue::IsFractionConstant<T>;
485
486 } // namespace ACFem
487
488} // Namespace Dune
489
490#endif // __DUNE_ACFEM_COMMON_FRACTIONCONSTANT_HH__
Dune::ACFem::operator<<
std::ostream & operator<<(std::ostream &out, TypeString< T > &&t)
Output operator for TypePrint tag-structure.
Definition: ostream.hh:39
Dune::ACFem::Tensor::pow
auto pow(T1 &&t1, T2 &&t2)
Power operations with scalar exponents are promoted to component-wise power operations by multiplying...
Definition: expressions.hh:525
Dune::ACFem::Constant
integral_constant< T, V > Constant
Short-cut for any integral constant.
Definition: types.hh:40
Dune::ACFem::BoolConstant
Constant< bool, V > BoolConstant
Short-cut for integral constant of type bool.
Definition: types.hh:48
Dune::ACFem::FalseType
BoolConstant< false > FalseType
Alias for std::false_type.
Definition: types.hh:110
Dune::ACFem::TrueType
BoolConstant< true > TrueType
Alias for std::true_type.
Definition: types.hh:107
std
STL namespace.
Dune::ACFem::TypedValue::FractionConstant
A class implementing compile-time constant fractions of integers.
Definition: fractionconstant.hh:27
Dune::ACFem::TypedValue::IsFractionConstant
Definition: fractionconstant.hh:75
Dune::ACFem::TypedValue::IsSign
Definition: fractionconstant.hh:173
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Jul 15, 22:36, 2024)