00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #ifndef __MMX_COMPLEX_HPP
00014 #define __MMX_COMPLEX_HPP
00015 #include <basix/double.hpp>
00016 #include <basix/function.hpp>
00017 #include <numerix/floating.hpp>
00018 #include <numerix/rounded.hpp>
00019 namespace mmx {
00020 #define TMPL_DEF template<typename C>
00021 #define TMPL template<typename C>
00022 #define Complex complex<C>
00023 TMPL class complex;
00024 TMPL inline C Re (const Complex& z);
00025 TMPL inline C Im (const Complex& z);
00026 TMPL inline C& Re (Complex& z);
00027 TMPL inline C& Im (Complex& z);
00028
00029
00030
00031
00032
00033 TMPL_DEF
00034 class complex {
00035 MMX_ALLOCATORS;
00036 private:
00037 C re;
00038 C im;
00039 public:
00040 inline complex () {}
00041 inline complex (const format<C>& fm):
00042 re (get_sample (fm)), im (get_sample (fm)) {}
00043 template<typename T> inline complex (const T& x):
00044 re (as<C> (x)), im (as<C> (promote (0, x))) {}
00045 template<typename T> inline complex (const complex<T>& z):
00046 re (Re (z)), im (Im (z)) {}
00047 template<typename XT,typename YT> inline
00048 complex (const XT& x, const YT& y):
00049 re (x), im (y) {}
00050 friend C Re LESSGTR (const Complex& z);
00051 friend C Im LESSGTR (const Complex& z);
00052 friend C& Re LESSGTR (Complex& z);
00053 friend C& Im LESSGTR (Complex& z);
00054
00055 inline Complex& operator += (const Complex& z);
00056 inline Complex& operator -= (const Complex& z);
00057 inline Complex& operator *= (const C& x);
00058 inline Complex& operator /= (const C& x);
00059 inline Complex& operator *= (const Complex& z);
00060 inline Complex& operator /= (const Complex& z);
00061 inline Complex& operator <<= (const xint& shift);
00062 inline Complex& operator >>= (const xint& shift);
00063 };
00064 DEFINE_UNARY_FORMAT_1 (complex)
00065
00066 STYPE_TO_TYPE(TMPL,scalar_type,Complex,C);
00067 UNARY_RETURN_TYPE(TMPL,Re_op,Complex,C);
00068 UNARY_RETURN_TYPE(TMPL,abs_op,Complex,C);
00069 UNARY_RETURN_TYPE(TMPL,center_op,Complex,complex<Center_type(C) >);
00070 UNARY_RETURN_TYPE(TMPL,radius_op,Complex,Radius_type(C));
00071 STYPE_TO_TYPE(TMPL,default_radius_type,Complex,Default_radius_type(C));
00072 TMPL struct rounding_helper<Complex >: public rounding_helper<C> {};
00073 TMPL struct projective_helper<Complex > { static const bool val= true; };
00074
00075 TMPL
00076 struct fast_helper<Complex > {
00077 typedef complex<Fast_type(C) > fast_type;
00078 static inline fast_type dd (const Complex& z) {
00079 return complex<Fast_type (C) > (fast (Re (z)), fast (Im (z))); }
00080 static inline Complex uu (const fast_type& z) {
00081 return complex<C> (slow<C> (Re (z)), slow<C> (Im (z))); }
00082 };
00083
00084 template<typename T,typename F>
00085 struct as_helper<complex<T>,complex<F> > {
00086 static inline complex<T> cv (const complex<F>& z) {
00087 return complex<T> (as<T> (Re (z)), as<T> (Im (z))); }
00088 };
00089
00090 template<typename T, typename F> inline void
00091 set_as (complex<T>& r, const complex<F>& z) {
00092 set_as (Re(r), Re(z));
00093 set_as (Im(r), Im(z));
00094 }
00095
00096 template<typename T, typename F> inline void
00097 set_as (complex<T>& r, const F& x) {
00098 set_as (Re(r), x);
00099 set_as (Im(r), promote (0, x));
00100 }
00101
00102
00103
00104
00105
00106 TMPL inline C Re (const Complex& z) { return z.re; }
00107 TMPL inline C Im (const Complex& z) { return z.im; }
00108 TMPL inline C& Re (Complex& z) { return z.re; }
00109 TMPL inline C& Im (Complex& z) { return z.im; }
00110 TMPL inline format<C> CF (const Complex& z) { return get_format (Re(z)); }
00111
00112 TMPL inline Complex copy (const Complex& z) {
00113 return Complex (copy (Re (z)), copy (Im (z))); }
00114 TMPL inline Complex duplicate (const Complex& z) {
00115 return Complex (duplicate (Re (z)), duplicate (Im (z))); }
00116 TMPL inline bool is_exact_zero (const Complex& z) {
00117 return is_exact_zero (Re (z)) && is_exact_zero (Im (z)); }
00118
00119 template<typename Op, typename C> inline nat
00120 unary_hash (const Complex& z) {
00121 nat h= Op::op (Re (z));
00122 return (h<<1) ^ (h<<5) ^ (h>>29) ^ Op::op (Im (z));
00123 }
00124
00125 template<typename Op, typename C> bool
00126 binary_test (const Complex& z1, const Complex& z2) {
00127 return Op::op (Re (z1), Re (z2)) && Op::op (Im (z1), Im (z2));
00128 }
00129
00130 template<typename S1, typename D> complex<D>
00131 map (const function_1<D,Argument(S1) >& fun, const complex<S1>& z) {
00132 return complex<D> (fun (Re (z)), fun (Im (z)));
00133 }
00134
00135 template<typename S1, typename D> complex<D>
00136 map (D (*fun) (const S1&), const complex<S1>& z) {
00137 return complex<D> (fun (Re (z)), fun (Im (z)));
00138 }
00139
00140 TRUE_IDENTITY_OP_SUGAR(TMPL,Complex)
00141 EXACT_IDENTITY_OP_SUGAR(TMPL,Complex)
00142 HARD_IDENTITY_OP_SUGAR(TMPL,Complex)
00143 EQUAL_INT_SUGAR(TMPL,Complex);
00144 EQUAL_SCALAR_SUGAR(TMPL,Complex,C);
00145 EQUAL_SCALAR_SUGAR_BIS(TMPL,Complex,C);
00146
00147 TMPL syntactic
00148 flatten (const Complex& z) {
00149 if (is_nan (z)) return Nan (syntactic);
00150 if (is_fuzz (z)) return Fuzz (syntactic);
00151 if (is_infinite (z)) return Infinity (syntactic);
00152 if (is_exact_zero (z)) return 0;
00153 if (!is_zero (Re (z)) && is_zero (Im (z)))
00154 return flatten (Re (z));
00155 if (is_zero (Re (z)) && !is_zero (Im (z)))
00156 return flatten (Im (z)) * Imaginary (syntactic);
00157 return flatten (Re (z)) + flatten (Im (z)) * Imaginary (syntactic);
00158 }
00159
00160 TMPL syntactic
00161 flatten (const Complex& z, xnat dd) {
00162 if (is_nan (z)) return Nan (syntactic);
00163 if (is_fuzz (z)) return Fuzz (syntactic);
00164 if (is_infinite (z)) return Infinity (syntactic);
00165 if (is_exact_zero (z)) return 0;
00166 if (!is_zero (Re (z)) && is_zero (Im (z)))
00167 return flatten (Re (z), dd);
00168 if (is_zero (Re (z)) && !is_zero (Im (z)))
00169 return flatten (Im (z), dd) * Imaginary (syntactic);
00170 return flatten (Re (z), dd) + flatten (Im (z), dd) * Imaginary (syntactic);
00171 }
00172
00173 TMPL
00174 struct binary_helper<Complex >: public void_binary_helper<Complex > {
00175 static inline string short_type_name () {
00176 return "Cc" * Short_type_name (C); }
00177 static inline generic full_type_name () {
00178 return gen ("Complex", Full_type_name (C)); }
00179 static inline generic disassemble (const Complex& z) {
00180 return gen_vec (as<generic> (Re (z)), as<generic> (Im (z))); }
00181 static inline Complex assemble (const generic& v) {
00182 return Complex (as<C> (vector_access (v, 0)),
00183 as<C> (vector_access (v, 1))); }
00184 static inline void write (const port& out, const Complex& z) {
00185 binary_write<C> (out, Re (z));
00186 binary_write<C> (out, Im (z)); }
00187 static inline Complex read (const port& in) {
00188 C r= binary_read<C> (in);
00189 C i= binary_read<C> (in);
00190 return Complex (r, i); }
00191 };
00192
00193
00194
00195
00196
00197 TMPL inline Complex
00198 gaussian (const C& re) {
00199 return Complex (re, promote (0, re));
00200 }
00201
00202 TMPL inline Complex
00203 gaussian (const C& re, const C& im) {
00204 return Complex (re, im);
00205 }
00206
00207 template<typename C, typename R> inline void
00208 gaussian_as (C& z, const R& x, const R& y) {
00209 z= gaussian (x, y);
00210 }
00211
00212 template<typename R, typename C> inline C
00213 gaussian (const R& re, const format<C>& fm) {
00214 C z= promote (0, fm); gaussian_as (z, re, promote (0, re)); return z;
00215 }
00216
00217 template<typename R, typename C> inline C
00218 gaussian (const R& re, const R& im, const format<C>& fm) {
00219 C z= promote (0, fm); gaussian_as (z, re, im); return z;
00220 }
00221
00222 TMPL inline Complex
00223 polar (const C& r, const C& t) {
00224 return Complex (r * cos (t), r * sin (t));
00225 }
00226
00227 TMPL inline C
00228 abs (const Complex& z) {
00229 return hypot (Re (z), Im (z));
00230 }
00231
00232 TMPL inline C
00233 arg (const Complex& z) {
00234 return atan2 (Im (z), Re (z));
00235 }
00236
00237 TMPL inline C
00238 norm (const Complex& z) {
00239 return square (Re (z)) + square (Im (z));
00240 }
00241
00242 TMPL inline Complex
00243 conj (const Complex& z) {
00244 return Complex (Re (z), -Im (z));
00245 }
00246
00247 TMPL inline Complex
00248 times_i (const Complex& z) {
00249 return Complex (-Im (z), Re (z));
00250 }
00251
00252 TMPL inline Complex
00253 over_i (const Complex& z) {
00254 return Complex (Im (z), -Re (z));
00255 }
00256
00257
00258
00259
00260
00261 TMPL inline void set_default (Complex& z) {
00262 set_default (Re (z)); set_zero (Im (z)); }
00263 TMPL inline void set_nan (Complex& z) {
00264 set_nan (Re (z)); set_nan (Im (z)); }
00265 TMPL inline void set_infinity (Complex& z) {
00266 set_infinity (Re (z)); set_infinity (Im (z)); }
00267 TMPL inline void set_fuzz (Complex& z) {
00268 set_fuzz (Re (z)); set_fuzz (Im (z)); }
00269 TMPL inline void set_smallest (Complex& z) {
00270 set_smallest (Re (z)); set_zero (Im (z)); }
00271 TMPL inline void set_largest (Complex& z) {
00272 set_largest (Re (z)); set_zero (Im (z)); }
00273 TMPL inline void set_accuracy (Complex& z) {
00274 set_accuracy (Re (z)); set_zero (Im (z)); }
00275 TMPL inline void set_log2 (Complex& z) {
00276 set_log2 (Re (z)); set_zero (Im (z)); }
00277 TMPL inline void set_pi (Complex& z) {
00278 set_pi (Re (z)); set_zero (Im (z)); }
00279 TMPL inline void set_euler (Complex& z) {
00280 set_euler (Re (z)); set_zero (Im (z)); }
00281 TMPL inline void set_imaginary (Complex& z) {
00282 set_zero (Re (z)); set_one (Im (z)); }
00283 TMPL inline Complex times_infinity (const Complex& z) {
00284 return z * infinity_cst<Complex > (); }
00285
00286
00287
00288
00289
00290 TMPL inline Complex
00291 operator + (const Complex& z1, const Complex& z2) {
00292 return Complex (Re (z1) + Re (z2), Im (z1) + Im (z2));
00293 }
00294
00295 TMPL inline Complex
00296 operator + (const Complex& z1, const C& z2) {
00297 return Complex (Re (z1) + z2, Im (z1));
00298 }
00299
00300 TMPL inline Complex
00301 operator + (const C& z1, const Complex& z2) {
00302 return Complex (z1 + Re (z2), Im (z2));
00303 }
00304
00305 TMPL inline Complex
00306 operator - (const Complex& z) {
00307 return Complex (-Re (z), -Im (z));
00308 }
00309
00310 TMPL inline Complex
00311 operator - (const Complex& z1, const Complex& z2) {
00312 return Complex (Re (z1) - Re (z2), Im (z1) - Im (z2));
00313 }
00314
00315 TMPL inline Complex
00316 operator - (const Complex& z1, const C& z2) {
00317 return Complex (Re (z1) - z2, Im (z1));
00318 }
00319
00320 TMPL inline Complex
00321 operator - (const C& z1, const Complex& z2) {
00322 return Complex (z1 - Re (z2), -Im (z2));
00323 }
00324
00325 TMPL inline Complex
00326 operator * (const C& x, const Complex& z) {
00327 return Complex (x * Re (z), x * Im (z));
00328 }
00329
00330 TMPL inline Complex
00331 operator * (const Complex& z, const C& x) {
00332 return Complex (Re (z) * x, Im (z) * x);
00333 }
00334
00335 TMPL inline Complex
00336 operator * (const Complex& z1, const Complex& z2) {
00337 return Complex (Re (z1) * Re (z2) - Im (z1) * Im (z2),
00338 Re (z1) * Im (z2) + Im (z1) * Re (z2));
00339 }
00340
00341 TMPL inline Complex
00342 square (const Complex& z) {
00343 C h= Re (z) * Im (z);
00344 return Complex (square (Re (z)) - square (Im (z)), h + h);
00345 }
00346
00347 TMPL inline Complex
00348 invert (const Complex& z) {
00349 C a= norm (z);
00350 return Complex (Re (z) / a, -Im (z) / a);
00351 }
00352
00353 TMPL inline Complex
00354 operator / (const Complex& z, const C& x) {
00355 return Complex (Re (z) / x, Im (z) / x);
00356 }
00357
00358 TMPL inline Complex
00359 operator / (const C& x, const Complex& z) {
00360 C a= x / norm (z);
00361 return Complex (a * Re (z), -a * Im (z));
00362 }
00363
00364 TMPL inline Complex
00365 operator / (const Complex& z1, const Complex& z2) {
00366 C a= norm (z2);
00367 return Complex ((Re (z1) * Re (z2) + Im (z1) * Im (z2)) / a,
00368 (Im (z1) * Re (z2) - Re (z1) * Im (z2)) / a);
00369 }
00370
00371 ARITH_SCALAR_INT_SUGAR(TMPL,Complex);
00372 GCD_SUGAR(TMPL,Complex);
00373
00374
00375
00376
00377
00378 TMPL inline Complex&
00379 Complex::operator += (const Complex& z) {
00380 re += z.re;
00381 im += z.im;
00382 return *this;
00383 }
00384
00385 TMPL inline Complex&
00386 Complex::operator -= (const Complex& z) {
00387 re -= z.re;
00388 im -= z.im;
00389 return *this;
00390 }
00391
00392 TMPL inline Complex&
00393 Complex::operator *= (const C& x) {
00394 re *= x;
00395 im *= x;
00396 return *this;
00397 }
00398
00399 TMPL inline Complex&
00400 Complex::operator /= (const C& x) {
00401 re /= x;
00402 im /= x;
00403 return *this;
00404 }
00405
00406 TMPL inline Complex&
00407 Complex::operator *= (const Complex& z) {
00408 C RE= re, zRE= z.re;
00409 re= re * z.re - im * z.im;
00410 im= RE * z.im + im * zRE;
00411 return *this;
00412 }
00413
00414 TMPL inline Complex&
00415 Complex::operator /= (const Complex& z) {
00416 C RE= re;
00417 C a = norm (z);
00418 re= (re * z.re + im * z.im) / a;
00419 im= (im * z.re - RE * z.im) / a;
00420 return *this;
00421 }
00422
00423 TMPL inline void
00424 neg (Complex& r) {
00425 neg (Re (r));
00426 neg (Im (r));
00427 }
00428
00429 TMPL inline void
00430 neg (Complex& r, const Complex& z1) {
00431 neg (Re (r), Re (z1));
00432 neg (Im (r), Im (z1));
00433 }
00434
00435 TMPL inline void
00436 add (Complex& r, const Complex& z1, const Complex& z2) {
00437 add (Re (r), Re (z1), Re (z2));
00438 add (Im (r), Im (z1), Im (z2));
00439 }
00440
00441 TMPL inline void
00442 sub (Complex& r, const Complex& z1, const Complex& z2) {
00443 sub (Re (r), Re (z1), Re (z2));
00444 sub (Im (r), Im (z1), Im (z2));
00445 }
00446
00447 TMPL inline void
00448 mul (Complex& r, const C& c, const Complex& z) {
00449 mul (Re (r), c, Re (z));
00450 mul (Im (r), c, Im (z));
00451 }
00452
00453 TMPL inline void
00454 mul (Complex& r, const Complex& z, const C& c) {
00455 mul (Re (r), Re (z), c);
00456 mul (Im (r), Im (z), c);
00457 }
00458
00459 TMPL inline void
00460 div (Complex& r, const Complex& z, const C& c) {
00461 div (Re (r), Re (z), c);
00462 div (Im (r), Im (z), c);
00463 }
00464
00465
00466
00467
00468
00469 TMPL Complex
00470 sqrt (const Complex& z) {
00471 C a= hypot (Re (z), Im (z));
00472 if (Re (z) > promote (0, a)) {
00473 C r= sqrt ((a + Re (z)) / promote (2, a));
00474 return Complex (r, (Im (z) / r) / promote (2, a));
00475 }
00476 else if (is_zero (a)) return promote (0, z);
00477 else {
00478 C i = sqrt ((a - Re (z)) / promote (2, a));
00479 if (Im (z) < promote (0, a)) i= -i;
00480 return Complex ((Im (z) / i) / promote (2, a), i);
00481 }
00482 }
00483
00484 TMPL inline Complex
00485 exp (const Complex& z) {
00486 return polar (exp (Re (z)), Im (z));
00487 }
00488
00489 TMPL inline Complex
00490 log (const Complex& z) {
00491 return Complex (log (abs (z)), arg (z));
00492 }
00493
00494 TMPL inline Complex
00495 pow (const Complex& z, const Complex& y) {
00496 C r= log (abs (z));
00497 C t= arg (z);
00498 return polar (exp (r * Re (y) - Im (y) * t),
00499 Im (y) * r + Re (y) * t);
00500 }
00501
00502 TMPL inline Complex
00503 pow (const Complex& z, const C& y) {
00504 return exp (y * log (z));
00505 }
00506
00507 TMPL inline Complex
00508 pow (const C& z, const Complex& y) {
00509 return exp (y * log (z));
00510 }
00511
00512 TMPL inline Complex
00513 pow (const Complex& z, const int& y) {
00514 return exp (y * log (z));
00515 }
00516
00517 TMPL inline Complex
00518 pow (const int& z, const Complex& y) {
00519 return exp (y * log (z));
00520 }
00521
00522 TMPL inline Complex
00523 cos (const Complex& z) {
00524 return Complex (cos (Re (z)) * cosh (Im (z)),
00525 -sin (Re (z)) * sinh (Im (z)));
00526 }
00527
00528 TMPL inline Complex
00529 sin (const Complex& z) {
00530 return Complex (sin (Re (z)) * cosh (Im (z)),
00531 cos (Re (z)) * sinh (Im (z)));
00532 }
00533
00534 TMPL inline Complex
00535 tan (const Complex& z) {
00536 return sin (z) / cos (z);
00537 }
00538
00539 TMPL inline Complex
00540 acos (const C& z) {
00541 return over_i (acosh (z));
00542 }
00543
00544 TMPL inline Complex
00545 asin (const C& z) {
00546 return over_i (asinh (times_i (z)));
00547 }
00548
00549 TMPL inline Complex
00550 atan (const C& z) {
00551 return over_i (atanh (times_i (z)));
00552 }
00553
00554 TMPL inline Complex
00555 acos (const Complex& z) {
00556 return over_i (acosh (z));
00557 }
00558
00559 TMPL inline Complex
00560 asin (const Complex& z) {
00561 return over_i (asinh (times_i (z)));
00562 }
00563
00564 TMPL inline Complex
00565 atan (const Complex& z) {
00566 return over_i (atanh (times_i (z)));
00567 }
00568
00569 TMPL inline Complex
00570 cosh (const Complex& z) {
00571 return Complex (cosh (Re (z)) * cos (Im (z)),
00572 sinh (Re (z)) * sin (Im (z)));
00573 }
00574
00575 TMPL inline Complex
00576 sinh (const Complex& z) {
00577 return Complex (sinh (Re (z)) * cos (Im (z)),
00578 cosh (Re (z)) * sin (Im (z)));
00579 }
00580
00581 TMPL inline Complex
00582 tanh (const Complex& z) {
00583 return sinh (z) / cosh (z);
00584 }
00585
00586 INV_TRIGO_SUGAR(TMPL,Complex)
00587 INV_HYPER_SUGAR(TMPL,Complex)
00588 ARG_HYPER_SUGAR(TMPL,Complex)
00589
00590
00591
00592
00593
00594 TMPL inline bool is_finite (const Complex& z) {
00595 return is_finite (Re (z)) && is_finite (Im (z)); }
00596 TMPL inline bool is_nan (const Complex& z) {
00597 return is_nan (Re (z)) || is_nan (Im (z)); }
00598 TMPL inline bool is_infinite (const Complex& z) {
00599 return !is_nan (z) && (is_infinite (Re (z)) || is_infinite (Im (z))); }
00600 TMPL inline bool is_fuzz (const Complex& z) {
00601 return !is_nan (z) && (is_fuzz (Re (z)) || is_fuzz (Im (z))); }
00602 TMPL inline bool is_reliable (const Complex& z) {
00603 return is_reliable (Re (z)); }
00604
00605 TMPL inline Complex change_precision (const Complex& z, xnat prec) {
00606 return Complex (change_precision (Re (z), prec),
00607 change_precision (Im (z), prec)); }
00608 TMPL inline xnat precision (const Complex& z) {
00609 return min (precision (Re (z)), precision (Im (z))); }
00610
00611 TMPL inline Abs_type(C) rounding_error (const complex<C>& z) {
00612 typedef Round_up(Abs_type(C) ) Up;
00613 return Up::add (rounding_error (Re (z)), rounding_error (Im (z))); }
00614 TMPL inline Abs_type(C) additive_error (const complex<C>& z) {
00615 typedef Round_up(Abs_type(C) ) Up;
00616 return Up::add (additive_error (Re (z)), additive_error (Im (z))); }
00617 TMPL inline Abs_type(C) multiplicative_error (const complex<C>& z) {
00618 typedef Round_up(Abs_type(C) ) Up;
00619 return incexp2 (Up::add (multiplicative_error (Re (z)),
00620 multiplicative_error (Im (z))), 3); }
00621 TMPL inline Abs_type(C) elementary_error (const complex<C>& z) {
00622 typedef Round_up(Abs_type(C) ) Up;
00623 return incexp2 (Up::add (elementary_error (Re (z)),
00624 elementary_error (Im (z))), 3); }
00625
00626
00627
00628
00629
00630
00631
00632
00633
00634
00635
00636
00637
00638 TMPL inline xint exponent (const Complex& z) {
00639 return max (exponent (Re (z)), exponent (Im (z))); }
00640 TMPL inline double magnitude (const Complex& z) {
00641 return magnitude (square (Re(z)) + square (Im(z))) / 2.0; }
00642 TMPL inline Complex operator << (const Complex& z, const xint& shift) {
00643 return Complex (incexp2 (Re (z), shift), incexp2 (Im (z), shift)); }
00644 TMPL inline Complex operator >> (const Complex& z, const xint& shift) {
00645 return Complex (decexp2 (Re (z), shift), decexp2 (Im (z), shift)); }
00646
00647 TMPL inline Complex&
00648 Complex::operator <<= (const xint& shift) {
00649 incexp2_assign (re, shift);
00650 incexp2_assign (im, shift);
00651 return *this;
00652 }
00653
00654 TMPL inline Complex&
00655 Complex::operator >>= (const xint& shift) {
00656 decexp2_assign (re, shift);
00657 decexp2_assign (im, shift);
00658 return *this;
00659 }
00660
00661 TMPL inline Complex sharpen (const Complex& z) {
00662 return Complex (sharpen (Re (z)), sharpen (Im (z))); }
00663 template<typename C, typename R> inline Complex
00664 blur (const Complex& z, const R& r) {
00665 return Complex (blur (Re (z), r), blur (Im (z), r)); }
00666
00667 #undef TMPL_DEF
00668 #undef TMPL
00669 #undef Complex
00670 }
00671 #endif // __MMX_COMPLEX_HPP