Merge of itasc branch. Project files, scons and cmake should be working. Makefile...
[blender.git] / extern / Eigen2 / Eigen / src / Core / MathFunctions.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra. Eigen itself is part of the KDE project.
3 //
4 // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
5 //
6 // Eigen is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 3 of the License, or (at your option) any later version.
10 //
11 // Alternatively, you can redistribute it and/or
12 // modify it under the terms of the GNU General Public License as
13 // published by the Free Software Foundation; either version 2 of
14 // the License, or (at your option) any later version.
15 //
16 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
17 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU Lesser General Public
22 // License and a copy of the GNU General Public License along with
23 // Eigen. If not, see <http://www.gnu.org/licenses/>.
24
25 #ifndef EIGEN_MATHFUNCTIONS_H
26 #define EIGEN_MATHFUNCTIONS_H
27
28 template<typename T> inline typename NumTraits<T>::Real precision();
29 template<typename T> inline typename NumTraits<T>::Real machine_epsilon();
30 template<typename T> inline T ei_random(T a, T b);
31 template<typename T> inline T ei_random();
32 template<typename T> inline T ei_random_amplitude()
33 {
34   if(NumTraits<T>::HasFloatingPoint) return static_cast<T>(1);
35   else return static_cast<T>(10);
36 }
37
38 template<typename T> inline T ei_hypot(T x, T y)
39 {
40   T _x = ei_abs(x);
41   T _y = ei_abs(y);
42   T p = std::max(_x, _y);
43   T q = std::min(_x, _y);
44   T qp = q/p;
45   return p * ei_sqrt(T(1) + qp*qp);
46 }
47
48 /**************
49 ***   int   ***
50 **************/
51
52 template<> inline int precision<int>() { return 0; }
53 template<> inline int machine_epsilon<int>() { return 0; }
54 inline int ei_real(int x)  { return x; }
55 inline int ei_imag(int)    { return 0; }
56 inline int ei_conj(int x)  { return x; }
57 inline int ei_abs(int x)   { return abs(x); }
58 inline int ei_abs2(int x)  { return x*x; }
59 inline int ei_sqrt(int)  { ei_assert(false); return 0; }
60 inline int ei_exp(int)  { ei_assert(false); return 0; }
61 inline int ei_log(int)  { ei_assert(false); return 0; }
62 inline int ei_sin(int)  { ei_assert(false); return 0; }
63 inline int ei_cos(int)  { ei_assert(false); return 0; }
64 inline int ei_atan2(int, int)  { ei_assert(false); return 0; }
65 inline int ei_pow(int x, int y) { return int(std::pow(double(x), y)); }
66
67 template<> inline int ei_random(int a, int b)
68 {
69   // We can't just do rand()%n as only the high-order bits are really random
70   return a + static_cast<int>((b-a+1) * (rand() / (RAND_MAX + 1.0)));
71 }
72 template<> inline int ei_random()
73 {
74   return ei_random<int>(-ei_random_amplitude<int>(), ei_random_amplitude<int>());
75 }
76 inline bool ei_isMuchSmallerThan(int a, int, int = precision<int>())
77 {
78   return a == 0;
79 }
80 inline bool ei_isApprox(int a, int b, int = precision<int>())
81 {
82   return a == b;
83 }
84 inline bool ei_isApproxOrLessThan(int a, int b, int = precision<int>())
85 {
86   return a <= b;
87 }
88
89 /**************
90 *** float   ***
91 **************/
92
93 template<> inline float precision<float>() { return 1e-5f; }
94 template<> inline float machine_epsilon<float>() { return 1.192e-07f; }
95 inline float ei_real(float x)  { return x; }
96 inline float ei_imag(float)    { return 0.f; }
97 inline float ei_conj(float x)  { return x; }
98 inline float ei_abs(float x)   { return std::abs(x); }
99 inline float ei_abs2(float x)  { return x*x; }
100 inline float ei_sqrt(float x)  { return std::sqrt(x); }
101 inline float ei_exp(float x)   { return std::exp(x); }
102 inline float ei_log(float x)   { return std::log(x); }
103 inline float ei_sin(float x)   { return std::sin(x); }
104 inline float ei_cos(float x)   { return std::cos(x); }
105 inline float ei_atan2(float y, float x) { return std::atan2(y,x); }
106 inline float ei_pow(float x, float y)  { return std::pow(x, y); }
107
108 template<> inline float ei_random(float a, float b)
109 {
110 #ifdef EIGEN_NICE_RANDOM
111   int i;
112   do { i = ei_random<int>(256*int(a),256*int(b));
113   } while(i==0);
114   return float(i)/256.f;
115 #else
116   return a + (b-a) * float(std::rand()) / float(RAND_MAX);
117 #endif
118 }
119 template<> inline float ei_random()
120 {
121   return ei_random<float>(-ei_random_amplitude<float>(), ei_random_amplitude<float>());
122 }
123 inline bool ei_isMuchSmallerThan(float a, float b, float prec = precision<float>())
124 {
125   return ei_abs(a) <= ei_abs(b) * prec;
126 }
127 inline bool ei_isApprox(float a, float b, float prec = precision<float>())
128 {
129   return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec;
130 }
131 inline bool ei_isApproxOrLessThan(float a, float b, float prec = precision<float>())
132 {
133   return a <= b || ei_isApprox(a, b, prec);
134 }
135
136 /**************
137 *** double  ***
138 **************/
139
140 template<> inline double precision<double>() { return 1e-11; }
141 template<> inline double machine_epsilon<double>() { return 2.220e-16; }
142
143 inline double ei_real(double x)  { return x; }
144 inline double ei_imag(double)    { return 0.; }
145 inline double ei_conj(double x)  { return x; }
146 inline double ei_abs(double x)   { return std::abs(x); }
147 inline double ei_abs2(double x)  { return x*x; }
148 inline double ei_sqrt(double x)  { return std::sqrt(x); }
149 inline double ei_exp(double x)   { return std::exp(x); }
150 inline double ei_log(double x)   { return std::log(x); }
151 inline double ei_sin(double x)   { return std::sin(x); }
152 inline double ei_cos(double x)   { return std::cos(x); }
153 inline double ei_atan2(double y, double x) { return std::atan2(y,x); }
154 inline double ei_pow(double x, double y) { return std::pow(x, y); }
155
156 template<> inline double ei_random(double a, double b)
157 {
158 #ifdef EIGEN_NICE_RANDOM
159   int i;
160   do { i= ei_random<int>(256*int(a),256*int(b));
161   } while(i==0);
162   return i/256.;
163 #else
164   return a + (b-a) * std::rand() / RAND_MAX;
165 #endif
166 }
167 template<> inline double ei_random()
168 {
169   return ei_random<double>(-ei_random_amplitude<double>(), ei_random_amplitude<double>());
170 }
171 inline bool ei_isMuchSmallerThan(double a, double b, double prec = precision<double>())
172 {
173   return ei_abs(a) <= ei_abs(b) * prec;
174 }
175 inline bool ei_isApprox(double a, double b, double prec = precision<double>())
176 {
177   return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec;
178 }
179 inline bool ei_isApproxOrLessThan(double a, double b, double prec = precision<double>())
180 {
181   return a <= b || ei_isApprox(a, b, prec);
182 }
183
184 /*********************
185 *** complex<float> ***
186 *********************/
187
188 template<> inline float precision<std::complex<float> >() { return precision<float>(); }
189 template<> inline float machine_epsilon<std::complex<float> >() { return machine_epsilon<float>(); }
190 inline float ei_real(const std::complex<float>& x) { return std::real(x); }
191 inline float ei_imag(const std::complex<float>& x) { return std::imag(x); }
192 inline std::complex<float> ei_conj(const std::complex<float>& x) { return std::conj(x); }
193 inline float ei_abs(const std::complex<float>& x) { return std::abs(x); }
194 inline float ei_abs2(const std::complex<float>& x) { return std::norm(x); }
195 inline std::complex<float> ei_exp(std::complex<float> x)  { return std::exp(x); }
196 inline std::complex<float> ei_sin(std::complex<float> x)  { return std::sin(x); }
197 inline std::complex<float> ei_cos(std::complex<float> x)  { return std::cos(x); }
198 inline std::complex<float> ei_atan2(std::complex<float>, std::complex<float> )  { ei_assert(false); return 0; }
199
200 template<> inline std::complex<float> ei_random()
201 {
202   return std::complex<float>(ei_random<float>(), ei_random<float>());
203 }
204 inline bool ei_isMuchSmallerThan(const std::complex<float>& a, const std::complex<float>& b, float prec = precision<float>())
205 {
206   return ei_abs2(a) <= ei_abs2(b) * prec * prec;
207 }
208 inline bool ei_isMuchSmallerThan(const std::complex<float>& a, float b, float prec = precision<float>())
209 {
210   return ei_abs2(a) <= ei_abs2(b) * prec * prec;
211 }
212 inline bool ei_isApprox(const std::complex<float>& a, const std::complex<float>& b, float prec = precision<float>())
213 {
214   return ei_isApprox(ei_real(a), ei_real(b), prec)
215       && ei_isApprox(ei_imag(a), ei_imag(b), prec);
216 }
217 // ei_isApproxOrLessThan wouldn't make sense for complex numbers
218
219 /**********************
220 *** complex<double> ***
221 **********************/
222
223 template<> inline double precision<std::complex<double> >() { return precision<double>(); }
224 template<> inline double machine_epsilon<std::complex<double> >() { return machine_epsilon<double>(); }
225 inline double ei_real(const std::complex<double>& x) { return std::real(x); }
226 inline double ei_imag(const std::complex<double>& x) { return std::imag(x); }
227 inline std::complex<double> ei_conj(const std::complex<double>& x) { return std::conj(x); }
228 inline double ei_abs(const std::complex<double>& x) { return std::abs(x); }
229 inline double ei_abs2(const std::complex<double>& x) { return std::norm(x); }
230 inline std::complex<double> ei_exp(std::complex<double> x)  { return std::exp(x); }
231 inline std::complex<double> ei_sin(std::complex<double> x)  { return std::sin(x); }
232 inline std::complex<double> ei_cos(std::complex<double> x)  { return std::cos(x); }
233 inline std::complex<double> ei_atan2(std::complex<double>, std::complex<double>)  { ei_assert(false); return 0; }
234
235 template<> inline std::complex<double> ei_random()
236 {
237   return std::complex<double>(ei_random<double>(), ei_random<double>());
238 }
239 inline bool ei_isMuchSmallerThan(const std::complex<double>& a, const std::complex<double>& b, double prec = precision<double>())
240 {
241   return ei_abs2(a) <= ei_abs2(b) * prec * prec;
242 }
243 inline bool ei_isMuchSmallerThan(const std::complex<double>& a, double b, double prec = precision<double>())
244 {
245   return ei_abs2(a) <= ei_abs2(b) * prec * prec;
246 }
247 inline bool ei_isApprox(const std::complex<double>& a, const std::complex<double>& b, double prec = precision<double>())
248 {
249   return ei_isApprox(ei_real(a), ei_real(b), prec)
250       && ei_isApprox(ei_imag(a), ei_imag(b), prec);
251 }
252 // ei_isApproxOrLessThan wouldn't make sense for complex numbers
253
254
255 /******************
256 *** long double ***
257 ******************/
258
259 template<> inline long double precision<long double>() { return precision<double>(); }
260 template<> inline long double machine_epsilon<long double>() { return 1.084e-19l; }
261 inline long double ei_real(long double x)  { return x; }
262 inline long double ei_imag(long double)    { return 0.; }
263 inline long double ei_conj(long double x)  { return x; }
264 inline long double ei_abs(long double x)   { return std::abs(x); }
265 inline long double ei_abs2(long double x)  { return x*x; }
266 inline long double ei_sqrt(long double x)  { return std::sqrt(x); }
267 inline long double ei_exp(long double x)   { return std::exp(x); }
268 inline long double ei_log(long double x)   { return std::log(x); }
269 inline long double ei_sin(long double x)   { return std::sin(x); }
270 inline long double ei_cos(long double x)   { return std::cos(x); }
271 inline long double ei_atan2(long double y, long double x) { return std::atan2(y,x); }
272 inline long double ei_pow(long double x, long double y)  { return std::pow(x, y); }
273
274 template<> inline long double ei_random(long double a, long double b)
275 {
276   return ei_random<double>(static_cast<double>(a),static_cast<double>(b));
277 }
278 template<> inline long double ei_random()
279 {
280   return ei_random<double>(-ei_random_amplitude<double>(), ei_random_amplitude<double>());
281 }
282 inline bool ei_isMuchSmallerThan(long double a, long double b, long double prec = precision<long double>())
283 {
284   return ei_abs(a) <= ei_abs(b) * prec;
285 }
286 inline bool ei_isApprox(long double a, long double b, long double prec = precision<long double>())
287 {
288   return ei_abs(a - b) <= std::min(ei_abs(a), ei_abs(b)) * prec;
289 }
290 inline bool ei_isApproxOrLessThan(long double a, long double b, long double prec = precision<long double>())
291 {
292   return a <= b || ei_isApprox(a, b, prec);
293 }
294
295 #endif // EIGEN_MATHFUNCTIONS_H