1 #ifdef COMPILATION// -*-indent-tabs-mode:t;c-basic-offset:4;tab-width:4;-*-
2 $CXX $0 -o $0x&&$0x&&rm $0x;exit
3 #endif
4 // © Alfredo Correa 2020
5
6 #ifndef MULTI_COMPLEX_HPP
7 #define MULTI_COMPLEX_HPP
8
9 #include "array_ref.hpp"
10
11 #include<complex>
12
13 namespace boost{
14 namespace multi{
15
16 constexpr class adl_conj_fn__{
17 template<class... As> auto _(priority<1>, As&&... as) const JUSTRETURN( std::conj(std::forward<As>(as)...))
18 template<class... As> auto _(priority<2>, As&&... as) const DECLRETURN( conj(std::forward<As>(as)...))
19 template<class T, class... As> auto _(priority<3>, T&& t, As&&... as) const DECLRETURN(std::forward<T>(t).conj(std::forward<As>(as)...))
20 public:
21 template<class... As> auto operator()(As&&... as) const DECLRETURN(_(priority<3>{}, std::forward<As>(as)...))
22 } adl_conj;
23
24 constexpr class adl_real_fn__{
25 template<class... As> auto _(priority<1>, As&&... as) const DECLRETURN( std::real(std::forward<As>(as)...))
26 template<class... As> auto _(priority<2>, As&&... as) const DECLRETURN( real(std::forward<As>(as)...))
27 template<class T, class... As> auto _(priority<3>, T&& t, As&&... as) const DECLRETURN(std::forward<T>(t).real(std::forward<As>(as)...))
28 public:
29 template<class... As> auto operator()(As&&... as) const DECLRETURN(_(priority<3>{}, std::forward<As>(as)...))
30 } adl_real;
31
32 constexpr class adl_imag_fn__{
33 template<class... As> auto _(priority<1>, As&&... as) const DECLRETURN( std::imag(std::forward<As>(as)...))
34 template<class... As> auto _(priority<2>, As&&... as) const DECLRETURN( imag(std::forward<As>(as)...))
35 template<class T, class... As> auto _(priority<3>, T&& t, As&&... as) const DECLRETURN(std::forward<T>(t).imag(std::forward<As>(as)...))
36 public:
37 template<class... As> auto operator()(As&&... as) const DECLRETURN(_(priority<3>{}, std::forward<As>(as)...))
38 } adl_imag;
39
40 struct real_t;
41 struct imag_t;
42
43 template<class ValueType = double>
44 struct complex{
45 using value_type = ValueType;
46 value_type re;
47 value_type im;
48 complex() = default;
49 // cppcheck-suppress noExplicitConstructor ; a real is a special complex without loss
complexboost::multi::complex50 constexpr complex(value_type real) : re{real}, im{value_type{0}}{}
complexboost::multi::complex51 constexpr complex(value_type real, value_type imag) : re{real}, im{imag}{}
52 // cppcheck-suppress noExplicitConstructor ; can be copied from standard type
complexboost::multi::complex53 constexpr complex(std::complex<ValueType> const& other) : re{other.real()}, im{other.imag()}{}
54 /* friend value_type const& real(complex const& c){return c.real;}
55 friend value_type & real(complex & c){return c.real;}
56 friend value_type const& imag(complex const& c){return c.imag;}
57 friend value_type & imag(complex & c){return c.imag;}*/
58 template<
59 class T,
60 std::enable_if_t<
61 sizeof(T)==2*sizeof(value_type) and
62 std::is_assignable<typename T::value_type&, decltype(std::declval<T>().real())>{} and
63 std::is_assignable<typename T::value_type&, decltype(std::declval<T>().imag())>{}, int
64 > =0
65 >
66 constexpr operator T const&() const&{return reinterpret_cast<T const&>(*this);}
67 template<
68 class T,
69 std::enable_if_t<
70 sizeof(T)==2*sizeof(value_type) and
71 std::is_assignable<typename T::value_type&, decltype(std::declval<T>().real())>{} and
72 std::is_assignable<typename T::value_type&, decltype(std::declval<T>().imag())>{}, int
73 > =0
74 >
75 constexpr operator T&()&{return reinterpret_cast<T const&>(*this);}
stdboost::multi::complex76 constexpr std::complex<value_type> const& std() const&{return reinterpret_cast<std::complex<value_type> const&>(*this);}
stdboost::multi::complex77 constexpr std::complex<value_type> & std() &{return reinterpret_cast<std::complex<value_type> &>(*this);}
abs(complex const & self)78 friend constexpr auto abs(complex const& self){return abs(self.std());}
operator -(complex const & self,complex const & other)79 friend constexpr complex operator-(complex const& self, complex const& other){return self.std() - other.std();}
realboost::multi::complex80 constexpr value_type& real()&{return re;}
realboost::multi::complex81 constexpr value_type const& real() const&{return re;}
82 // friend constexpr value_type& real(complex& self){return self.re;}
83 // friend constexpr value_type const& real(complex const& self){return self.re;}
84
imagboost::multi::complex85 constexpr value_type& imag()&{return im;}
imagboost::multi::complex86 constexpr value_type const& imag() const&{return im;}
87 // friend constexpr value_type& imag(complex& self){return self.im;}
88 // friend constexpr value_type const& imag(complex const& self){return self.im;}
89
operator +=boost::multi::complex90 template<class Real> constexpr auto operator+=(Real const& other)&->decltype(re += other, *this){return re += other, *this;}
operator -=boost::multi::complex91 template<class Real> constexpr auto operator-=(Real const& other)&->decltype(re -= other, *this){return re -= other, *this;}
operator *=boost::multi::complex92 template<class Real> constexpr auto operator*=(Real const& other)&->decltype(re *= other, im *= other, *this){return re *= other, im *= other, *this;}
operator /=boost::multi::complex93 template<class Real> constexpr auto operator/=(Real const& other)&->decltype(re /= other, im /= other, *this){return re /= other, im /= other, *this;}
94
operator +=boost::multi::complex95 template<class Complex> constexpr auto operator+=(Complex const& other)&->decltype(re += other.re, im += other.im, *this){return re += other.re, im += other.im, *this;}
operator -=boost::multi::complex96 template<class Complex> constexpr auto operator-=(Complex const& other)&->decltype(re -= other.re, im -= other.im, *this){return re -= other.re, im -= other.im, *this;}
97 };
98
99 struct real_t{
100 template<class Array, typename E = typename std::decay_t<Array>::element, typename ValueType = typename E::value_type>
operator ()boost::multi::real_t101 auto operator()(Array&& a) const
102 ->decltype(std::forward<Array>(a).template reinterpret_array_cast<complex<ValueType>>().template member_cast<ValueType>(&complex<ValueType>::real)){
103 return std::forward<Array>(a).template reinterpret_array_cast<complex<ValueType>>().template member_cast<ValueType>(&complex<ValueType>::real);}
104 template<class T, typename ValueType = typename std::decay_t<T>::value_type,
105 std::enable_if_t<
106 sizeof(T)==2*sizeof(ValueType) and
107 std::is_assignable<ValueType&, decltype(real(std::declval<T>()))>{} and
108 std::is_assignable<ValueType&, decltype(imag(std::declval<T>()))>{}, int
109 > =0
110 >
111 ValueType& operator()(T& t) const{return reinterpret_cast<multi::complex<ValueType>&>(t).real;}
112 template<class T, typename ValueType = typename std::decay_t<T>::value_type,
113 std::enable_if_t<
114 sizeof(T)==2*sizeof(ValueType) and
115 std::is_assignable<ValueType&, decltype(real(std::declval<T>()))>{} and
116 std::is_assignable<ValueType&, decltype(imag(std::declval<T>()))>{}, int
117 > =0
118 >
119 ValueType const& operator()(T const& t) const{return reinterpret_cast<multi::complex<ValueType> const&>(t).real;}
120 };
121
122 struct imag_t{
123 template<class Array, typename E = typename std::decay_t<Array>::element, typename ValueType = typename E::value_type>
operator ()boost::multi::imag_t124 auto operator()(Array&& a) const
125 ->decltype(std::forward<Array>(a).template reinterpret_array_cast<complex<ValueType>>().template member_cast<ValueType>(&complex<ValueType>::imag)){
126 return std::forward<Array>(a).template reinterpret_array_cast<complex<ValueType>>().template member_cast<ValueType>(&complex<ValueType>::imag);}
127 template<class T, typename ValueType = typename std::decay_t<T>::value_type,
128 std::enable_if_t<
129 sizeof(T)==2*sizeof(ValueType) and
130 std::is_assignable<ValueType&, decltype(real(std::declval<T>()))>{} and
131 std::is_assignable<ValueType&, decltype(imag(std::declval<T>()))>{}, int
132 > =0
133 >
134 ValueType& operator()(T& t) const{return reinterpret_cast<multi::complex<ValueType>&>(t).imag;}
135 template<class T, typename ValueType = typename std::decay_t<T>::value_type,
136 std::enable_if_t<
137 sizeof(T)==2*sizeof(ValueType) and
138 std::is_assignable<ValueType&, decltype(real(std::declval<T>()))>{} and
139 std::is_assignable<ValueType&, decltype(imag(std::declval<T>()))>{}, int
140 > =0
141 >
142 ValueType const& operator()(T const& t) const{return reinterpret_cast<multi::complex<ValueType> const&>(t).imag;}
143 };
144
145 /*
146 template<class T2, class Array, class P2 = typename std::pointer_traits<typename std::decay<Array>::type::element_ptr>::template rebind<T2>,
147 typename E = typename std::decay_t<Array>::element,
148 typename R = decltype(real(E{})),
149 std::enable_if_t<sizeof(E)==2*sizeof(typename E::value_type), int> =0
150 >
151 decltype(auto) member_array_cast(Array&& a, real_t const*){
152 struct Complex{double real; double imag;};
153 return multi::member_array_cast<double>(multi::reinterpret_array_cast<Complex>(std::forward<Array>(a)), &Complex::real);
154 }
155 template<class T2, class Array, class P2 = typename std::pointer_traits<typename std::decay<Array>::type::element_ptr>::template rebind<T2>,
156 typename E = typename std::decay_t<Array>::element,
157 typename R = decltype(real(E{})),
158 std::enable_if_t<sizeof(E)==2*sizeof(typename E::value_type), int> =0
159 >
160 decltype(auto) member_array_cast(Array&& a, imag_t const*){
161 struct Complex{double real; double imag;};
162 return multi::member_array_cast<double>(multi::reinterpret_array_cast<Complex>(std::forward<Array>(a)), &Complex::imag);
163 }
164 */
165
166 static constexpr real_t real MAYBE_UNUSED;
167 static constexpr imag_t imag MAYBE_UNUSED;
168
169 }}
170
171 namespace std{
172 template<class T>
173 struct is_trivially_default_constructible<std::complex<T>> :
174 is_trivially_default_constructible<T>{};
175 }
176
177 #if not __INCLUDE_LEVEL__ // _TEST_MULTI_COMPLEX
178
179 #include<cassert>
180 #include<complex>
181 #include "array.hpp"
182
183 namespace multi = boost::multi;
184
185 template<class T> void what(T&&)=delete;
186
main()187 int main(){
188 static_assert( std::is_trivially_default_constructible<std::complex<double>>{}, "!");
189 static_assert( std::is_trivially_copy_constructible<std::complex<double>>{}, "!");
190 static_assert( std::is_trivially_assignable<std::complex<double>&, std::complex<double> const>{}, "!");
191
192 using complex = multi::complex<double>;
193
194 multi::array<complex, 2> A = {
195 { {1.,2.}, {3.,4.} },
196 { {22.,33.}, {5.,9.} }
197 };
198
199 {
200 auto&& Areal = A.member_cast<double>(&multi::complex<double>::re);
201 auto&& Aimag = A.member_cast<double>(&multi::complex<double>::im);
202
203 assert( Areal[1][0] == 22. );
204 assert( Aimag[1][0] == 33. );
205 }
206 {
207 auto&& Areal = A.member_cast<double>(&multi::complex<double>::re);
208 auto&& Aimag = A.member_cast<double>(&multi::complex<double>::im);
209
210 assert( Areal[1][0] == 22. );
211 assert( Aimag[1][0] == 33. );
212 }
213 {
214 // auto&& Areal = multi::real(A); // multi::real(A);
215 // auto&& Aimag = multi::imag(A); // multi::real(A);
216
217 // assert( &Areal[1][0] == &multi::real(A[1][0]) );
218 // assert( &Aimag[1][0] == &multi::imag(A[1][0]) );
219 }
220 }
221
222 #endif
223 #endif
224
225
226