AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
QuadMath.hpp
1#pragma once
2
3#if HAVE_QUADMATH
4#include <functional>
5#include <type_traits>
6
7#include <fmt/format.h>
8
9#include <dune/common/hash.hh>
10#include <dune/common/quadmath.hh>
11
12#include <amdis/common/DerivativeTraits.hpp>
13#include <amdis/common/StaticSize.hpp>
14#include <amdis/common/ValueCategory.hpp>
15#include <amdis/gridfunctions/ConstantGridFunction.hpp>
16
17namespace std
18{
19 template <class T>
20 struct common_type<Dune::Float128, T>
21 {
22 using type = Dune::Float128;
23 };
24
25 template <class T>
26 struct common_type<T, Dune::Float128>
27 {
28 using type = Dune::Float128;
29 };
30
31 template <>
32 struct common_type<Dune::Float128, Dune::Float128>
33 {
34 using type = Dune::Float128;
35 };
36
37 template <>
38 struct hash<Dune::Float128>
39 {
40 typedef Dune::Float128 argument_type;
41 typedef std::size_t result_type;
42
43 std::size_t operator()(const Dune::Float128& arg) const
44 {
45 hash<long double> hasher_ld;
46 return hasher_ld((long double)(arg));
47 }
48 };
49
50 template <>
51 struct hash<const Dune::Float128>
52 : public hash<Dune::Float128>
53 {};
54
55} // end namespace std
56
57
58namespace Dune
59{
60 namespace Impl
61 {
62 // specialization for float arguments due to ambiguity
63 template <class T,
64 std::enable_if_t<not std::is_integral_v<T> && std::is_arithmetic_v<T>, int> = 0>
65 inline Float128 pow(const Float128& x, const T& p)
66 {
67 return powq(float128_t(x), float128_t(p));
68 }
69
70 } // end namespace Impl
71} // end namespace Dune
72
73
74namespace AMDiS
75{
76 namespace Concepts
77 {
78 namespace Definition
79 {
80 template <>
81 struct ConstantToGridFunction<Dune::Float128>
82 : std::true_type {};
83
84 } // end namespace Definition
85 } // end namespace Concepts
86
87 namespace Impl
88 {
89 template <>
90 struct SizeImpl<Dune::Float128>
91 {
92 static constexpr auto eval(Dune::Float128)
93 -> std::integral_constant<std::size_t, 1> { return {}; }
94 };
95
96 template <>
97 struct NumRowsImpl<Dune::Float128>
98 {
99 static constexpr auto eval(Dune::Float128)
100 -> std::integral_constant<std::size_t, 1> { return {}; }
101 };
102
103 template <>
104 struct NumColsImpl<Dune::Float128>
105 {
106 static constexpr auto eval(Dune::Float128)
107 -> std::integral_constant<std::size_t, 1> { return {}; }
108 };
109
110 } // end namespace Impl
111
112 template <class K, int N>
113 struct DerivativeTraits<Dune::Float128(Dune::FieldVector<K,N>), tag::gradient>
114 {
115 using Range = Dune::FieldVector<Dune::Float128,N>;
116 };
117
118 template <>
119 struct ValueCategory<Dune::Float128>
120 {
121 using type = tag::scalar;
122 };
123
124} // end namespace AMDiS
125
126template <>
127struct fmt::formatter<Dune::Float128>
128 : public formatter<long double>
129{
130 template <class FormatContext>
131 auto format(Dune::Float128 const& f, FormatContext& ctx) const
132 {
133 return formatter<long double>::format(static_cast<long double>(f), ctx);
134 }
135};
136
137
138#endif // HAVE_QUADMATH
constexpr bool ConstantToGridFunction
Concepts that is true for all ''simple'' types that can be converted automatically to a GridFunction,...
Definition: ConstantGridFunction.hpp:187