1// The template and inlines for the -*- C++ -*- complex number classes. 2 3// Copyright (C) 1997, 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 4// Free Software Foundation, Inc. 5// 6// This file is part of the GNU ISO C++ Library. This library is free 7// software; you can redistribute it and/or modify it under the 8// terms of the GNU General Public License as published by the 9// Free Software Foundation; either version 2, or (at your option) 10// any later version. 11 12// This library is distributed in the hope that it will be useful, 13// but WITHOUT ANY WARRANTY; without even the implied warranty of 14// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15// GNU General Public License for more details. 16 17// You should have received a copy of the GNU General Public License along 18// with this library; see the file COPYING. If not, write to the Free 19// Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, 20// USA. 21 22// As a special exception, you may use this file as part of a free software 23// library without restriction. Specifically, if other files instantiate 24// templates or use macros or inline functions from this file, or you compile 25// this file and link it with other files to produce an executable, this 26// file does not by itself cause the resulting executable to be covered by 27// the GNU General Public License. This exception does not however 28// invalidate any other reasons why the executable file might be covered by 29// the GNU General Public License. 30 31/** @file complex 32 * This is a Standard C++ Library header. 33 */ 34 35// 36// ISO C++ 14882: 26.2 Complex Numbers 37// Note: this is not a conforming implementation. 38// Initially implemented by Ulrich Drepper <drepper@cygnus.com> 39// Improved by Gabriel Dos Reis <dosreis@cmla.ens-cachan.fr> 40// 41 42#ifndef _GLIBCXX_COMPLEX 43#define _GLIBCXX_COMPLEX 1 44 45#pragma GCC system_header 46 47#include <bits/c++config.h> 48#include <bits/cpp_type_traits.h> 49#include <cmath> 50#include <sstream> 51 52_GLIBCXX_BEGIN_NAMESPACE(std) 53 54 // Forward declarations. 55 template<typename _Tp> class complex; 56 template<> class complex<float>; 57 template<> class complex<double>; 58 template<> class complex<long double>; 59 60 /// Return magnitude of @a z. 61 template<typename _Tp> _Tp abs(const complex<_Tp>&); 62 /// Return phase angle of @a z. 63 template<typename _Tp> _Tp arg(const complex<_Tp>&); 64 /// Return @a z magnitude squared. 65 template<typename _Tp> _Tp norm(const complex<_Tp>&); 66 67 /// Return complex conjugate of @a z. 68 template<typename _Tp> complex<_Tp> conj(const complex<_Tp>&); 69 /// Return complex with magnitude @a rho and angle @a theta. 70 template<typename _Tp> complex<_Tp> polar(const _Tp&, const _Tp& = 0); 71 72 // Transcendentals: 73 /// Return complex cosine of @a z. 74 template<typename _Tp> complex<_Tp> cos(const complex<_Tp>&); 75 /// Return complex hyperbolic cosine of @a z. 76 template<typename _Tp> complex<_Tp> cosh(const complex<_Tp>&); 77 /// Return complex base e exponential of @a z. 78 template<typename _Tp> complex<_Tp> exp(const complex<_Tp>&); 79 /// Return complex natural logarithm of @a z. 80 template<typename _Tp> complex<_Tp> log(const complex<_Tp>&); 81 /// Return complex base 10 logarithm of @a z. 82 template<typename _Tp> complex<_Tp> log10(const complex<_Tp>&); 83 /// Return complex cosine of @a z. 84 template<typename _Tp> complex<_Tp> pow(const complex<_Tp>&, int); 85 /// Return @a x to the @a y'th power. 86 template<typename _Tp> complex<_Tp> pow(const complex<_Tp>&, const _Tp&); 87 /// Return @a x to the @a y'th power. 88 template<typename _Tp> complex<_Tp> pow(const complex<_Tp>&, 89 const complex<_Tp>&); 90 /// Return @a x to the @a y'th power. 91 template<typename _Tp> complex<_Tp> pow(const _Tp&, const complex<_Tp>&); 92 /// Return complex sine of @a z. 93 template<typename _Tp> complex<_Tp> sin(const complex<_Tp>&); 94 /// Return complex hyperbolic sine of @a z. 95 template<typename _Tp> complex<_Tp> sinh(const complex<_Tp>&); 96 /// Return complex square root of @a z. 97 template<typename _Tp> complex<_Tp> sqrt(const complex<_Tp>&); 98 /// Return complex tangent of @a z. 99 template<typename _Tp> complex<_Tp> tan(const complex<_Tp>&); 100 /// Return complex hyperbolic tangent of @a z. 101 template<typename _Tp> complex<_Tp> tanh(const complex<_Tp>&); 102 //@} 103 104 105 // 26.2.2 Primary template class complex 106 /** 107 * Template to represent complex numbers. 108 * 109 * Specializations for float, double, and long double are part of the 110 * library. Results with any other type are not guaranteed. 111 * 112 * @param Tp Type of real and imaginary values. 113 */ 114 template<typename _Tp> 115 struct complex 116 { 117 /// Value typedef. 118 typedef _Tp value_type; 119 120 /// Default constructor. First parameter is x, second parameter is y. 121 /// Unspecified parameters default to 0. 122 complex(const _Tp& = _Tp(), const _Tp & = _Tp()); 123 124 // Lets the compiler synthesize the copy constructor 125 // complex (const complex<_Tp>&); 126 /// Copy constructor. 127 template<typename _Up> 128 complex(const complex<_Up>&); 129 130 /// Return real part of complex number. 131 _Tp& real(); 132 /// Return real part of complex number. 133 const _Tp& real() const; 134 /// Return imaginary part of complex number. 135 _Tp& imag(); 136 /// Return imaginary part of complex number. 137 const _Tp& imag() const; 138 139 /// Assign this complex number to scalar @a t. 140 complex<_Tp>& operator=(const _Tp&); 141 /// Add @a t to this complex number. 142 complex<_Tp>& operator+=(const _Tp&); 143 /// Subtract @a t from this complex number. 144 complex<_Tp>& operator-=(const _Tp&); 145 /// Multiply this complex number by @a t. 146 complex<_Tp>& operator*=(const _Tp&); 147 /// Divide this complex number by @a t. 148 complex<_Tp>& operator/=(const _Tp&); 149 150 // Lets the compiler synthesize the 151 // copy and assignment operator 152 // complex<_Tp>& operator= (const complex<_Tp>&); 153 /// Assign this complex number to complex @a z. 154 template<typename _Up> 155 complex<_Tp>& operator=(const complex<_Up>&); 156 /// Add @a z to this complex number. 157 template<typename _Up> 158 complex<_Tp>& operator+=(const complex<_Up>&); 159 /// Subtract @a z from this complex number. 160 template<typename _Up> 161 complex<_Tp>& operator-=(const complex<_Up>&); 162 /// Multiply this complex number by @a z. 163 template<typename _Up> 164 complex<_Tp>& operator*=(const complex<_Up>&); 165 /// Divide this complex number by @a z. 166 template<typename _Up> 167 complex<_Tp>& operator/=(const complex<_Up>&); 168 169 const complex& __rep() const; 170 171 private: 172 _Tp _M_real; 173 _Tp _M_imag; 174 }; 175 176 template<typename _Tp> 177 inline _Tp& 178 complex<_Tp>::real() { return _M_real; } 179 180 template<typename _Tp> 181 inline const _Tp& 182 complex<_Tp>::real() const { return _M_real; } 183 184 template<typename _Tp> 185 inline _Tp& 186 complex<_Tp>::imag() { return _M_imag; } 187 188 template<typename _Tp> 189 inline const _Tp& 190 complex<_Tp>::imag() const { return _M_imag; } 191 192 template<typename _Tp> 193 inline 194 complex<_Tp>::complex(const _Tp& __r, const _Tp& __i) 195 : _M_real(__r), _M_imag(__i) { } 196 197 template<typename _Tp> 198 template<typename _Up> 199 inline 200 complex<_Tp>::complex(const complex<_Up>& __z) 201 : _M_real(__z.real()), _M_imag(__z.imag()) { } 202 203 template<typename _Tp> 204 complex<_Tp>& 205 complex<_Tp>::operator=(const _Tp& __t) 206 { 207 _M_real = __t; 208 _M_imag = _Tp(); 209 return *this; 210 } 211 212 // 26.2.5/1 213 template<typename _Tp> 214 inline complex<_Tp>& 215 complex<_Tp>::operator+=(const _Tp& __t) 216 { 217 _M_real += __t; 218 return *this; 219 } 220 221 // 26.2.5/3 222 template<typename _Tp> 223 inline complex<_Tp>& 224 complex<_Tp>::operator-=(const _Tp& __t) 225 { 226 _M_real -= __t; 227 return *this; 228 } 229 230 // 26.2.5/5 231 template<typename _Tp> 232 complex<_Tp>& 233 complex<_Tp>::operator*=(const _Tp& __t) 234 { 235 _M_real *= __t; 236 _M_imag *= __t; 237 return *this; 238 } 239 240 // 26.2.5/7 241 template<typename _Tp> 242 complex<_Tp>& 243 complex<_Tp>::operator/=(const _Tp& __t) 244 { 245 _M_real /= __t; 246 _M_imag /= __t; 247 return *this; 248 } 249 250 template<typename _Tp> 251 template<typename _Up> 252 complex<_Tp>& 253 complex<_Tp>::operator=(const complex<_Up>& __z) 254 { 255 _M_real = __z.real(); 256 _M_imag = __z.imag(); 257 return *this; 258 } 259 260 // 26.2.5/9 261 template<typename _Tp> 262 template<typename _Up> 263 complex<_Tp>& 264 complex<_Tp>::operator+=(const complex<_Up>& __z) 265 { 266 _M_real += __z.real(); 267 _M_imag += __z.imag(); 268 return *this; 269 } 270 271 // 26.2.5/11 272 template<typename _Tp> 273 template<typename _Up> 274 complex<_Tp>& 275 complex<_Tp>::operator-=(const complex<_Up>& __z) 276 { 277 _M_real -= __z.real(); 278 _M_imag -= __z.imag(); 279 return *this; 280 } 281 282 // 26.2.5/13 283 // XXX: This is a grammar school implementation. 284 template<typename _Tp> 285 template<typename _Up> 286 complex<_Tp>& 287 complex<_Tp>::operator*=(const complex<_Up>& __z) 288 { 289 const _Tp __r = _M_real * __z.real() - _M_imag * __z.imag(); 290 _M_imag = _M_real * __z.imag() + _M_imag * __z.real(); 291 _M_real = __r; 292 return *this; 293 } 294 295 // 26.2.5/15 296 // XXX: This is a grammar school implementation. 297 template<typename _Tp> 298 template<typename _Up> 299 complex<_Tp>& 300 complex<_Tp>::operator/=(const complex<_Up>& __z) 301 { 302 const _Tp __r = _M_real * __z.real() + _M_imag * __z.imag(); 303 const _Tp __n = std::norm(__z); 304 _M_imag = (_M_imag * __z.real() - _M_real * __z.imag()) / __n; 305 _M_real = __r / __n; 306 return *this; 307 } 308 309 template<typename _Tp> 310 inline const complex<_Tp>& 311 complex<_Tp>::__rep() const { return *this; } 312 313 // Operators: 314 //@{ 315 /// Return new complex value @a x plus @a y. 316 template<typename _Tp> 317 inline complex<_Tp> 318 operator+(const complex<_Tp>& __x, const complex<_Tp>& __y) 319 { 320 complex<_Tp> __r = __x; 321 __r += __y; 322 return __r; 323 } 324 325 template<typename _Tp> 326 inline complex<_Tp> 327 operator+(const complex<_Tp>& __x, const _Tp& __y) 328 { 329 complex<_Tp> __r = __x; 330 __r.real() += __y; 331 return __r; 332 } 333 334 template<typename _Tp> 335 inline complex<_Tp> 336 operator+(const _Tp& __x, const complex<_Tp>& __y) 337 { 338 complex<_Tp> __r = __y; 339 __r.real() += __x; 340 return __r; 341 } 342 //@} 343 344 //@{ 345 /// Return new complex value @a x minus @a y. 346 template<typename _Tp> 347 inline complex<_Tp> 348 operator-(const complex<_Tp>& __x, const complex<_Tp>& __y) 349 { 350 complex<_Tp> __r = __x; 351 __r -= __y; 352 return __r; 353 } 354 355 template<typename _Tp> 356 inline complex<_Tp> 357 operator-(const complex<_Tp>& __x, const _Tp& __y) 358 { 359 complex<_Tp> __r = __x; 360 __r.real() -= __y; 361 return __r; 362 } 363 364 template<typename _Tp> 365 inline complex<_Tp> 366 operator-(const _Tp& __x, const complex<_Tp>& __y) 367 { 368 complex<_Tp> __r(__x, -__y.imag()); 369 __r.real() -= __y.real(); 370 return __r; 371 } 372 //@} 373 374 //@{ 375 /// Return new complex value @a x times @a y. 376 template<typename _Tp> 377 inline complex<_Tp> 378 operator*(const complex<_Tp>& __x, const complex<_Tp>& __y) 379 { 380 complex<_Tp> __r = __x; 381 __r *= __y; 382 return __r; 383 } 384 385 template<typename _Tp> 386 inline complex<_Tp> 387 operator*(const complex<_Tp>& __x, const _Tp& __y) 388 { 389 complex<_Tp> __r = __x; 390 __r *= __y; 391 return __r; 392 } 393 394 template<typename _Tp> 395 inline complex<_Tp> 396 operator*(const _Tp& __x, const complex<_Tp>& __y) 397 { 398 complex<_Tp> __r = __y; 399 __r *= __x; 400 return __r; 401 } 402 //@} 403 404 //@{ 405 /// Return new complex value @a x divided by @a y. 406 template<typename _Tp> 407 inline complex<_Tp> 408 operator/(const complex<_Tp>& __x, const complex<_Tp>& __y) 409 { 410 complex<_Tp> __r = __x; 411 __r /= __y; 412 return __r; 413 } 414 415 template<typename _Tp> 416 inline complex<_Tp> 417 operator/(const complex<_Tp>& __x, const _Tp& __y) 418 { 419 complex<_Tp> __r = __x; 420 __r /= __y; 421 return __r; 422 } 423 424 template<typename _Tp> 425 inline complex<_Tp> 426 operator/(const _Tp& __x, const complex<_Tp>& __y) 427 { 428 complex<_Tp> __r = __x; 429 __r /= __y; 430 return __r; 431 } 432 //@} 433 434 /// Return @a x. 435 template<typename _Tp> 436 inline complex<_Tp> 437 operator+(const complex<_Tp>& __x) 438 { return __x; } 439 440 /// Return complex negation of @a x. 441 template<typename _Tp> 442 inline complex<_Tp> 443 operator-(const complex<_Tp>& __x) 444 { return complex<_Tp>(-__x.real(), -__x.imag()); } 445 446 //@{ 447 /// Return true if @a x is equal to @a y. 448 template<typename _Tp> 449 inline bool 450 operator==(const complex<_Tp>& __x, const complex<_Tp>& __y) 451 { return __x.real() == __y.real() && __x.imag() == __y.imag(); } 452 453 template<typename _Tp> 454 inline bool 455 operator==(const complex<_Tp>& __x, const _Tp& __y) 456 { return __x.real() == __y && __x.imag() == _Tp(); } 457 458 template<typename _Tp> 459 inline bool 460 operator==(const _Tp& __x, const complex<_Tp>& __y) 461 { return __x == __y.real() && _Tp() == __y.imag(); } 462 //@} 463 464 //@{ 465 /// Return false if @a x is equal to @a y. 466 template<typename _Tp> 467 inline bool 468 operator!=(const complex<_Tp>& __x, const complex<_Tp>& __y) 469 { return __x.real() != __y.real() || __x.imag() != __y.imag(); } 470 471 template<typename _Tp> 472 inline bool 473 operator!=(const complex<_Tp>& __x, const _Tp& __y) 474 { return __x.real() != __y || __x.imag() != _Tp(); } 475 476 template<typename _Tp> 477 inline bool 478 operator!=(const _Tp& __x, const complex<_Tp>& __y) 479 { return __x != __y.real() || _Tp() != __y.imag(); } 480 //@} 481 482 /// Extraction operator for complex values. 483 template<typename _Tp, typename _CharT, class _Traits> 484 basic_istream<_CharT, _Traits>& 485 operator>>(basic_istream<_CharT, _Traits>& __is, complex<_Tp>& __x) 486 { 487 _Tp __re_x, __im_x; 488 _CharT __ch; 489 __is >> __ch; 490 if (__ch == '(') 491 { 492 __is >> __re_x >> __ch; 493 if (__ch == ',') 494 { 495 __is >> __im_x >> __ch; 496 if (__ch == ')') 497 __x = complex<_Tp>(__re_x, __im_x); 498 else 499 __is.setstate(ios_base::failbit); 500 } 501 else if (__ch == ')') 502 __x = __re_x; 503 else 504 __is.setstate(ios_base::failbit); 505 } 506 else 507 { 508 __is.putback(__ch); 509 __is >> __re_x; 510 __x = __re_x; 511 } 512 return __is; 513 } 514 515 /// Insertion operator for complex values. 516 template<typename _Tp, typename _CharT, class _Traits> 517 basic_ostream<_CharT, _Traits>& 518 operator<<(basic_ostream<_CharT, _Traits>& __os, const complex<_Tp>& __x) 519 { 520 basic_ostringstream<_CharT, _Traits> __s; 521 __s.flags(__os.flags()); 522 __s.imbue(__os.getloc()); 523 __s.precision(__os.precision()); 524 __s << '(' << __x.real() << ',' << __x.imag() << ')'; 525 return __os << __s.str(); 526 } 527 528 // Values 529 template<typename _Tp> 530 inline _Tp& 531 real(complex<_Tp>& __z) 532 { return __z.real(); } 533 534 template<typename _Tp> 535 inline const _Tp& 536 real(const complex<_Tp>& __z) 537 { return __z.real(); } 538 539 template<typename _Tp> 540 inline _Tp& 541 imag(complex<_Tp>& __z) 542 { return __z.imag(); } 543 544 template<typename _Tp> 545 inline const _Tp& 546 imag(const complex<_Tp>& __z) 547 { return __z.imag(); } 548 549 // 26.2.7/3 abs(__z): Returns the magnitude of __z. 550 template<typename _Tp> 551 inline _Tp 552 __complex_abs(const complex<_Tp>& __z) 553 { 554 _Tp __x = __z.real(); 555 _Tp __y = __z.imag(); 556 const _Tp __s = std::max(abs(__x), abs(__y)); 557 if (__s == _Tp()) // well ... 558 return __s; 559 __x /= __s; 560 __y /= __s; 561 return __s * sqrt(__x * __x + __y * __y); 562 } 563 564#if _GLIBCXX_USE_C99_COMPLEX 565 inline float 566 __complex_abs(__complex__ float __z) { return __builtin_cabsf(__z); } 567 568 inline double 569 __complex_abs(__complex__ double __z) { return __builtin_cabs(__z); } 570 571 inline long double 572 __complex_abs(const __complex__ long double& __z) 573 { return __builtin_cabsl(__z); } 574 575 template<typename _Tp> 576 inline _Tp 577 abs(const complex<_Tp>& __z) { return __complex_abs(__z.__rep()); } 578#else 579 template<typename _Tp> 580 inline _Tp 581 abs(const complex<_Tp>& __z) { return __complex_abs(__z); } 582#endif 583 584 585 // 26.2.7/4: arg(__z): Returns the phase angle of __z. 586 template<typename _Tp> 587 inline _Tp 588 __complex_arg(const complex<_Tp>& __z) 589 { return atan2(__z.imag(), __z.real()); } 590 591#if _GLIBCXX_USE_C99_COMPLEX 592 inline float 593 __complex_arg(__complex__ float __z) { return __builtin_cargf(__z); } 594 595 inline double 596 __complex_arg(__complex__ double __z) { return __builtin_carg(__z); } 597 598 inline long double 599 __complex_arg(const __complex__ long double& __z) 600 { return __builtin_cargl(__z); } 601 602 template<typename _Tp> 603 inline _Tp 604 arg(const complex<_Tp>& __z) { return __complex_arg(__z.__rep()); } 605#else 606 template<typename _Tp> 607 inline _Tp 608 arg(const complex<_Tp>& __z) { return __complex_arg(__z); } 609#endif 610 611 // 26.2.7/5: norm(__z) returns the squared magintude of __z. 612 // As defined, norm() is -not- a norm is the common mathematical 613 // sens used in numerics. The helper class _Norm_helper<> tries to 614 // distinguish between builtin floating point and the rest, so as 615 // to deliver an answer as close as possible to the real value. 616 template<bool> 617 struct _Norm_helper 618 { 619 template<typename _Tp> 620 static inline _Tp _S_do_it(const complex<_Tp>& __z) 621 { 622 const _Tp __x = __z.real(); 623 const _Tp __y = __z.imag(); 624 return __x * __x + __y * __y; 625 } 626 }; 627 628 template<> 629 struct _Norm_helper<true> 630 { 631 template<typename _Tp> 632 static inline _Tp _S_do_it(const complex<_Tp>& __z) 633 { 634 _Tp __res = std::abs(__z); 635 return __res * __res; 636 } 637 }; 638 639 template<typename _Tp> 640 inline _Tp 641 norm(const complex<_Tp>& __z) 642 { 643 return _Norm_helper<__is_floating<_Tp>::__value 644 && !_GLIBCXX_FAST_MATH>::_S_do_it(__z); 645 } 646 647 template<typename _Tp> 648 inline complex<_Tp> 649 polar(const _Tp& __rho, const _Tp& __theta) 650 { return complex<_Tp>(__rho * cos(__theta), __rho * sin(__theta)); } 651 652 template<typename _Tp> 653 inline complex<_Tp> 654 conj(const complex<_Tp>& __z) 655 { return complex<_Tp>(__z.real(), -__z.imag()); } 656 657 // Transcendentals 658 659 // 26.2.8/1 cos(__z): Returns the cosine of __z. 660 template<typename _Tp> 661 inline complex<_Tp> 662 __complex_cos(const complex<_Tp>& __z) 663 { 664 const _Tp __x = __z.real(); 665 const _Tp __y = __z.imag(); 666 return complex<_Tp>(cos(__x) * cosh(__y), -sin(__x) * sinh(__y)); 667 } 668 669#if _GLIBCXX_USE_C99_COMPLEX 670 inline __complex__ float 671 __complex_cos(__complex__ float __z) { return __builtin_ccosf(__z); } 672 673 inline __complex__ double 674 __complex_cos(__complex__ double __z) { return __builtin_ccos(__z); } 675 676 inline __complex__ long double 677 __complex_cos(const __complex__ long double& __z) 678 { return __builtin_ccosl(__z); } 679 680 template<typename _Tp> 681 inline complex<_Tp> 682 cos(const complex<_Tp>& __z) { return __complex_cos(__z.__rep()); } 683#else 684 template<typename _Tp> 685 inline complex<_Tp> 686 cos(const complex<_Tp>& __z) { return __complex_cos(__z); } 687#endif 688 689 // 26.2.8/2 cosh(__z): Returns the hyperbolic cosine of __z. 690 template<typename _Tp> 691 inline complex<_Tp> 692 __complex_cosh(const complex<_Tp>& __z) 693 { 694 const _Tp __x = __z.real(); 695 const _Tp __y = __z.imag(); 696 return complex<_Tp>(cosh(__x) * cos(__y), sinh(__x) * sin(__y)); 697 } 698 699#if _GLIBCXX_USE_C99_COMPLEX 700 inline __complex__ float 701 __complex_cosh(__complex__ float __z) { return __builtin_ccoshf(__z); } 702 703 inline __complex__ double 704 __complex_cosh(__complex__ double __z) { return __builtin_ccosh(__z); } 705 706 inline __complex__ long double 707 __complex_cosh(const __complex__ long double& __z) 708 { return __builtin_ccoshl(__z); } 709 710 template<typename _Tp> 711 inline complex<_Tp> 712 cosh(const complex<_Tp>& __z) { return __complex_cosh(__z.__rep()); } 713#else 714 template<typename _Tp> 715 inline complex<_Tp> 716 cosh(const complex<_Tp>& __z) { return __complex_cosh(__z); } 717#endif 718 719 // 26.2.8/3 exp(__z): Returns the complex base e exponential of x 720 template<typename _Tp> 721 inline complex<_Tp> 722 __complex_exp(const complex<_Tp>& __z) 723 { return std::polar(exp(__z.real()), __z.imag()); } 724 725#if _GLIBCXX_USE_C99_COMPLEX 726 inline __complex__ float 727 __complex_exp(__complex__ float __z) { return __builtin_cexpf(__z); } 728 729 inline __complex__ double 730 __complex_exp(__complex__ double __z) { return __builtin_cexp(__z); } 731 732 inline __complex__ long double 733 __complex_exp(const __complex__ long double& __z) 734 { return __builtin_cexpl(__z); } 735 736 template<typename _Tp> 737 inline complex<_Tp> 738 exp(const complex<_Tp>& __z) { return __complex_exp(__z.__rep()); } 739#else 740 template<typename _Tp> 741 inline complex<_Tp> 742 exp(const complex<_Tp>& __z) { return __complex_exp(__z); } 743#endif 744 745 // 26.2.8/5 log(__z): Reurns the natural complex logaritm of __z. 746 // The branch cut is along the negative axis. 747 template<typename _Tp> 748 inline complex<_Tp> 749 __complex_log(const complex<_Tp>& __z) 750 { return complex<_Tp>(log(std::abs(__z)), std::arg(__z)); } 751 752#if _GLIBCXX_USE_C99_COMPLEX 753 inline __complex__ float 754 __complex_log(__complex__ float __z) { return __builtin_clogf(__z); } 755 756 inline __complex__ double 757 __complex_log(__complex__ double __z) { return __builtin_clog(__z); } 758 759 inline __complex__ long double 760 __complex_log(const __complex__ long double& __z) 761 { return __builtin_clogl(__z); } 762 763 template<typename _Tp> 764 inline complex<_Tp> 765 log(const complex<_Tp>& __z) { return __complex_log(__z.__rep()); } 766#else 767 template<typename _Tp> 768 inline complex<_Tp> 769 log(const complex<_Tp>& __z) { return __complex_log(__z); } 770#endif 771 772 template<typename _Tp> 773 inline complex<_Tp> 774 log10(const complex<_Tp>& __z) 775 { return std::log(__z) / log(_Tp(10.0)); } 776 777 // 26.2.8/10 sin(__z): Returns the sine of __z. 778 template<typename _Tp> 779 inline complex<_Tp> 780 __complex_sin(const complex<_Tp>& __z) 781 { 782 const _Tp __x = __z.real(); 783 const _Tp __y = __z.imag(); 784 return complex<_Tp>(sin(__x) * cosh(__y), cos(__x) * sinh(__y)); 785 } 786 787#if _GLIBCXX_USE_C99_COMPLEX 788 inline __complex__ float 789 __complex_sin(__complex__ float __z) { return __builtin_csinf(__z); } 790 791 inline __complex__ double 792 __complex_sin(__complex__ double __z) { return __builtin_csin(__z); } 793 794 inline __complex__ long double 795 __complex_sin(const __complex__ long double& __z) 796 { return __builtin_csinl(__z); } 797 798 template<typename _Tp> 799 inline complex<_Tp> 800 sin(const complex<_Tp>& __z) { return __complex_sin(__z.__rep()); } 801#else 802 template<typename _Tp> 803 inline complex<_Tp> 804 sin(const complex<_Tp>& __z) { return __complex_sin(__z); } 805#endif 806 807 // 26.2.8/11 sinh(__z): Returns the hyperbolic sine of __z. 808 template<typename _Tp> 809 inline complex<_Tp> 810 __complex_sinh(const complex<_Tp>& __z) 811 { 812 const _Tp __x = __z.real(); 813 const _Tp __y = __z.imag(); 814 return complex<_Tp>(sinh(__x) * cos(__y), cosh(__x) * sin(__y)); 815 } 816 817#if _GLIBCXX_USE_C99_COMPLEX 818 inline __complex__ float 819 __complex_sinh(__complex__ float __z) { return __builtin_csinhf(__z); } 820 821 inline __complex__ double 822 __complex_sinh(__complex__ double __z) { return __builtin_csinh(__z); } 823 824 inline __complex__ long double 825 __complex_sinh(const __complex__ long double& __z) 826 { return __builtin_csinhl(__z); } 827 828 template<typename _Tp> 829 inline complex<_Tp> 830 sinh(const complex<_Tp>& __z) { return __complex_sinh(__z.__rep()); } 831#else 832 template<typename _Tp> 833 inline complex<_Tp> 834 sinh(const complex<_Tp>& __z) { return __complex_sinh(__z); } 835#endif 836 837 // 26.2.8/13 sqrt(__z): Returns the complex square root of __z. 838 // The branch cut is on the negative axis. 839 template<typename _Tp> 840 complex<_Tp> 841 __complex_sqrt(const complex<_Tp>& __z) 842 { 843 _Tp __x = __z.real(); 844 _Tp __y = __z.imag(); 845 846 if (__x == _Tp()) 847 { 848 _Tp __t = sqrt(abs(__y) / 2); 849 return complex<_Tp>(__t, __y < _Tp() ? -__t : __t); 850 } 851 else 852 { 853 _Tp __t = sqrt(2 * (std::abs(__z) + abs(__x))); 854 _Tp __u = __t / 2; 855 return __x > _Tp() 856 ? complex<_Tp>(__u, __y / __t) 857 : complex<_Tp>(abs(__y) / __t, __y < _Tp() ? -__u : __u); 858 } 859 } 860 861#if _GLIBCXX_USE_C99_COMPLEX 862 inline __complex__ float 863 __complex_sqrt(__complex__ float __z) { return __builtin_csqrtf(__z); } 864 865 inline __complex__ double 866 __complex_sqrt(__complex__ double __z) { return __builtin_csqrt(__z); } 867 868 inline __complex__ long double 869 __complex_sqrt(const __complex__ long double& __z) 870 { return __builtin_csqrtl(__z); } 871 872 template<typename _Tp> 873 inline complex<_Tp> 874 sqrt(const complex<_Tp>& __z) { return __complex_sqrt(__z.__rep()); } 875#else 876 template<typename _Tp> 877 inline complex<_Tp> 878 sqrt(const complex<_Tp>& __z) { return __complex_sqrt(__z); } 879#endif 880 881 // 26.2.8/14 tan(__z): Return the complex tangent of __z. 882 883 template<typename _Tp> 884 inline complex<_Tp> 885 __complex_tan(const complex<_Tp>& __z) 886 { return std::sin(__z) / std::cos(__z); } 887 888#if _GLIBCXX_USE_C99_COMPLEX 889 inline __complex__ float 890 __complex_tan(__complex__ float __z) { return __builtin_ctanf(__z); } 891 892 inline __complex__ double 893 __complex_tan(__complex__ double __z) { return __builtin_ctan(__z); } 894 895 inline __complex__ long double 896 __complex_tan(const __complex__ long double& __z) 897 { return __builtin_ctanl(__z); } 898 899 template<typename _Tp> 900 inline complex<_Tp> 901 tan(const complex<_Tp>& __z) { return __complex_tan(__z.__rep()); } 902#else 903 template<typename _Tp> 904 inline complex<_Tp> 905 tan(const complex<_Tp>& __z) { return __complex_tan(__z); } 906#endif 907 908 909 // 26.2.8/15 tanh(__z): Returns the hyperbolic tangent of __z. 910 911 template<typename _Tp> 912 inline complex<_Tp> 913 __complex_tanh(const complex<_Tp>& __z) 914 { return std::sinh(__z) / std::cosh(__z); } 915 916#if _GLIBCXX_USE_C99_COMPLEX 917 inline __complex__ float 918 __complex_tanh(__complex__ float __z) { return __builtin_ctanhf(__z); } 919 920 inline __complex__ double 921 __complex_tanh(__complex__ double __z) { return __builtin_ctanh(__z); } 922 923 inline __complex__ long double 924 __complex_tanh(const __complex__ long double& __z) 925 { return __builtin_ctanhl(__z); } 926 927 template<typename _Tp> 928 inline complex<_Tp> 929 tanh(const complex<_Tp>& __z) { return __complex_tanh(__z.__rep()); } 930#else 931 template<typename _Tp> 932 inline complex<_Tp> 933 tanh(const complex<_Tp>& __z) { return __complex_tanh(__z); } 934#endif 935 936 937 // 26.2.8/9 pow(__x, __y): Returns the complex power base of __x 938 // raised to the __y-th power. The branch 939 // cut is on the negative axis. 940 template<typename _Tp> 941 inline complex<_Tp> 942 pow(const complex<_Tp>& __z, int __n) 943 { return std::__pow_helper(__z, __n); } 944 945 template<typename _Tp> 946 complex<_Tp> 947 pow(const complex<_Tp>& __x, const _Tp& __y) 948 { 949#ifndef _GLIBCXX_USE_C99_COMPLEX 950 if (__x == _Tp()) 951 return _Tp(); 952#endif 953 if (__x.imag() == _Tp() && __x.real() > _Tp()) 954 return pow(__x.real(), __y); 955 956 complex<_Tp> __t = std::log(__x); 957 return std::polar(exp(__y * __t.real()), __y * __t.imag()); 958 } 959 960 template<typename _Tp> 961 inline complex<_Tp> 962 __complex_pow(const complex<_Tp>& __x, const complex<_Tp>& __y) 963 { return __x == _Tp() ? _Tp() : std::exp(__y * std::log(__x)); } 964 965#if _GLIBCXX_USE_C99_COMPLEX 966 inline __complex__ float 967 __complex_pow(__complex__ float __x, __complex__ float __y) 968 { return __builtin_cpowf(__x, __y); } 969 970 inline __complex__ double 971 __complex_pow(__complex__ double __x, __complex__ double __y) 972 { return __builtin_cpow(__x, __y); } 973 974 inline __complex__ long double 975 __complex_pow(const __complex__ long double& __x, 976 const __complex__ long double& __y) 977 { return __builtin_cpowl(__x, __y); } 978 979 template<typename _Tp> 980 inline complex<_Tp> 981 pow(const complex<_Tp>& __x, const complex<_Tp>& __y) 982 { return __complex_pow(__x.__rep(), __y.__rep()); } 983#else 984 template<typename _Tp> 985 inline complex<_Tp> 986 pow(const complex<_Tp>& __x, const complex<_Tp>& __y) 987 { return __complex_pow(__x, __y); } 988#endif 989 990 template<typename _Tp> 991 inline complex<_Tp> 992 pow(const _Tp& __x, const complex<_Tp>& __y) 993 { 994 return __x > _Tp() ? std::polar(pow(__x, __y.real()), 995 __y.imag() * log(__x)) 996 : std::pow(complex<_Tp>(__x, _Tp()), __y); 997 } 998 999 // 26.2.3 complex specializations 1000 // complex<float> specialization 1001 template<> 1002 struct complex<float> 1003 { 1004 typedef float value_type; 1005 typedef __complex__ float _ComplexT; 1006 1007 complex(_ComplexT __z) : _M_value(__z) { } 1008 1009 complex(float = 0.0f, float = 0.0f); 1010 1011 explicit complex(const complex<double>&); 1012 explicit complex(const complex<long double>&); 1013 1014 float& real(); 1015 const float& real() const; 1016 float& imag(); 1017 const float& imag() const; 1018 1019 complex<float>& operator=(float); 1020 complex<float>& operator+=(float); 1021 complex<float>& operator-=(float); 1022 complex<float>& operator*=(float); 1023 complex<float>& operator/=(float); 1024 1025 // Let's the compiler synthetize the copy and assignment 1026 // operator. It always does a pretty good job. 1027 // complex& operator= (const complex&); 1028 template<typename _Tp> 1029 complex<float>&operator=(const complex<_Tp>&); 1030 template<typename _Tp> 1031 complex<float>& operator+=(const complex<_Tp>&); 1032 template<class _Tp> 1033 complex<float>& operator-=(const complex<_Tp>&); 1034 template<class _Tp> 1035 complex<float>& operator*=(const complex<_Tp>&); 1036 template<class _Tp> 1037 complex<float>&operator/=(const complex<_Tp>&); 1038 1039 const _ComplexT& __rep() const { return _M_value; } 1040 1041 private: 1042 _ComplexT _M_value; 1043 }; 1044 1045 inline float& 1046 complex<float>::real() 1047 { return __real__ _M_value; } 1048 1049 inline const float& 1050 complex<float>::real() const 1051 { return __real__ _M_value; } 1052 1053 inline float& 1054 complex<float>::imag() 1055 { return __imag__ _M_value; } 1056 1057 inline const float& 1058 complex<float>::imag() const 1059 { return __imag__ _M_value; } 1060 1061 inline 1062 complex<float>::complex(float r, float i) 1063 { 1064 __real__ _M_value = r; 1065 __imag__ _M_value = i; 1066 } 1067 1068 inline complex<float>& 1069 complex<float>::operator=(float __f) 1070 { 1071 __real__ _M_value = __f; 1072 __imag__ _M_value = 0.0f; 1073 return *this; 1074 } 1075 1076 inline complex<float>& 1077 complex<float>::operator+=(float __f) 1078 { 1079 __real__ _M_value += __f; 1080 return *this; 1081 } 1082 1083 inline complex<float>& 1084 complex<float>::operator-=(float __f) 1085 { 1086 __real__ _M_value -= __f; 1087 return *this; 1088 } 1089 1090 inline complex<float>& 1091 complex<float>::operator*=(float __f) 1092 { 1093 _M_value *= __f; 1094 return *this; 1095 } 1096 1097 inline complex<float>& 1098 complex<float>::operator/=(float __f) 1099 { 1100 _M_value /= __f; 1101 return *this; 1102 } 1103 1104 template<typename _Tp> 1105 inline complex<float>& 1106 complex<float>::operator=(const complex<_Tp>& __z) 1107 { 1108 __real__ _M_value = __z.real(); 1109 __imag__ _M_value = __z.imag(); 1110 return *this; 1111 } 1112 1113 template<typename _Tp> 1114 inline complex<float>& 1115 complex<float>::operator+=(const complex<_Tp>& __z) 1116 { 1117 __real__ _M_value += __z.real(); 1118 __imag__ _M_value += __z.imag(); 1119 return *this; 1120 } 1121 1122 template<typename _Tp> 1123 inline complex<float>& 1124 complex<float>::operator-=(const complex<_Tp>& __z) 1125 { 1126 __real__ _M_value -= __z.real(); 1127 __imag__ _M_value -= __z.imag(); 1128 return *this; 1129 } 1130 1131 template<typename _Tp> 1132 inline complex<float>& 1133 complex<float>::operator*=(const complex<_Tp>& __z) 1134 { 1135 _ComplexT __t; 1136 __real__ __t = __z.real(); 1137 __imag__ __t = __z.imag(); 1138 _M_value *= __t; 1139 return *this; 1140 } 1141 1142 template<typename _Tp> 1143 inline complex<float>& 1144 complex<float>::operator/=(const complex<_Tp>& __z) 1145 { 1146 _ComplexT __t; 1147 __real__ __t = __z.real(); 1148 __imag__ __t = __z.imag(); 1149 _M_value /= __t; 1150 return *this; 1151 } 1152 1153 // 26.2.3 complex specializations 1154 // complex<double> specialization 1155 template<> 1156 struct complex<double> 1157 { 1158 typedef double value_type; 1159 typedef __complex__ double _ComplexT; 1160 1161 complex(_ComplexT __z) : _M_value(__z) { } 1162 1163 complex(double = 0.0, double = 0.0); 1164 1165 complex(const complex<float>&); 1166 explicit complex(const complex<long double>&); 1167 1168 double& real(); 1169 const double& real() const; 1170 double& imag(); 1171 const double& imag() const; 1172 1173 complex<double>& operator=(double); 1174 complex<double>& operator+=(double); 1175 complex<double>& operator-=(double); 1176 complex<double>& operator*=(double); 1177 complex<double>& operator/=(double); 1178 1179 // The compiler will synthetize this, efficiently. 1180 // complex& operator= (const complex&); 1181 template<typename _Tp> 1182 complex<double>& operator=(const complex<_Tp>&); 1183 template<typename _Tp> 1184 complex<double>& operator+=(const complex<_Tp>&); 1185 template<typename _Tp> 1186 complex<double>& operator-=(const complex<_Tp>&); 1187 template<typename _Tp> 1188 complex<double>& operator*=(const complex<_Tp>&); 1189 template<typename _Tp> 1190 complex<double>& operator/=(const complex<_Tp>&); 1191 1192 const _ComplexT& __rep() const { return _M_value; } 1193 1194 private: 1195 _ComplexT _M_value; 1196 }; 1197 1198 inline double& 1199 complex<double>::real() 1200 { return __real__ _M_value; } 1201 1202 inline const double& 1203 complex<double>::real() const 1204 { return __real__ _M_value; } 1205 1206 inline double& 1207 complex<double>::imag() 1208 { return __imag__ _M_value; } 1209 1210 inline const double& 1211 complex<double>::imag() const 1212 { return __imag__ _M_value; } 1213 1214 inline 1215 complex<double>::complex(double __r, double __i) 1216 { 1217 __real__ _M_value = __r; 1218 __imag__ _M_value = __i; 1219 } 1220 1221 inline complex<double>& 1222 complex<double>::operator=(double __d) 1223 { 1224 __real__ _M_value = __d; 1225 __imag__ _M_value = 0.0; 1226 return *this; 1227 } 1228 1229 inline complex<double>& 1230 complex<double>::operator+=(double __d) 1231 { 1232 __real__ _M_value += __d; 1233 return *this; 1234 } 1235 1236 inline complex<double>& 1237 complex<double>::operator-=(double __d) 1238 { 1239 __real__ _M_value -= __d; 1240 return *this; 1241 } 1242 1243 inline complex<double>& 1244 complex<double>::operator*=(double __d) 1245 { 1246 _M_value *= __d; 1247 return *this; 1248 } 1249 1250 inline complex<double>& 1251 complex<double>::operator/=(double __d) 1252 { 1253 _M_value /= __d; 1254 return *this; 1255 } 1256 1257 template<typename _Tp> 1258 inline complex<double>& 1259 complex<double>::operator=(const complex<_Tp>& __z) 1260 { 1261 __real__ _M_value = __z.real(); 1262 __imag__ _M_value = __z.imag(); 1263 return *this; 1264 } 1265 1266 template<typename _Tp> 1267 inline complex<double>& 1268 complex<double>::operator+=(const complex<_Tp>& __z) 1269 { 1270 __real__ _M_value += __z.real(); 1271 __imag__ _M_value += __z.imag(); 1272 return *this; 1273 } 1274 1275 template<typename _Tp> 1276 inline complex<double>& 1277 complex<double>::operator-=(const complex<_Tp>& __z) 1278 { 1279 __real__ _M_value -= __z.real(); 1280 __imag__ _M_value -= __z.imag(); 1281 return *this; 1282 } 1283 1284 template<typename _Tp> 1285 inline complex<double>& 1286 complex<double>::operator*=(const complex<_Tp>& __z) 1287 { 1288 _ComplexT __t; 1289 __real__ __t = __z.real(); 1290 __imag__ __t = __z.imag(); 1291 _M_value *= __t; 1292 return *this; 1293 } 1294 1295 template<typename _Tp> 1296 inline complex<double>& 1297 complex<double>::operator/=(const complex<_Tp>& __z) 1298 { 1299 _ComplexT __t; 1300 __real__ __t = __z.real(); 1301 __imag__ __t = __z.imag(); 1302 _M_value /= __t; 1303 return *this; 1304 } 1305 1306 // 26.2.3 complex specializations 1307 // complex<long double> specialization 1308 template<> 1309 struct complex<long double> 1310 { 1311 typedef long double value_type; 1312 typedef __complex__ long double _ComplexT; 1313 1314 complex(_ComplexT __z) : _M_value(__z) { } 1315 1316 complex(long double = 0.0L, long double = 0.0L); 1317 1318 complex(const complex<float>&); 1319 complex(const complex<double>&); 1320 1321 long double& real(); 1322 const long double& real() const; 1323 long double& imag(); 1324 const long double& imag() const; 1325 1326 complex<long double>& operator= (long double); 1327 complex<long double>& operator+= (long double); 1328 complex<long double>& operator-= (long double); 1329 complex<long double>& operator*= (long double); 1330 complex<long double>& operator/= (long double); 1331 1332 // The compiler knows how to do this efficiently 1333 // complex& operator= (const complex&); 1334 template<typename _Tp> 1335 complex<long double>& operator=(const complex<_Tp>&); 1336 template<typename _Tp> 1337 complex<long double>& operator+=(const complex<_Tp>&); 1338 template<typename _Tp> 1339 complex<long double>& operator-=(const complex<_Tp>&); 1340 template<typename _Tp> 1341 complex<long double>& operator*=(const complex<_Tp>&); 1342 template<typename _Tp> 1343 complex<long double>& operator/=(const complex<_Tp>&); 1344 1345 const _ComplexT& __rep() const { return _M_value; } 1346 1347 private: 1348 _ComplexT _M_value; 1349 }; 1350 1351 inline 1352 complex<long double>::complex(long double __r, long double __i) 1353 { 1354 __real__ _M_value = __r; 1355 __imag__ _M_value = __i; 1356 } 1357 1358 inline long double& 1359 complex<long double>::real() 1360 { return __real__ _M_value; } 1361 1362 inline const long double& 1363 complex<long double>::real() const 1364 { return __real__ _M_value; } 1365 1366 inline long double& 1367 complex<long double>::imag() 1368 { return __imag__ _M_value; } 1369 1370 inline const long double& 1371 complex<long double>::imag() const 1372 { return __imag__ _M_value; } 1373 1374 inline complex<long double>& 1375 complex<long double>::operator=(long double __r) 1376 { 1377 __real__ _M_value = __r; 1378 __imag__ _M_value = 0.0L; 1379 return *this; 1380 } 1381 1382 inline complex<long double>& 1383 complex<long double>::operator+=(long double __r) 1384 { 1385 __real__ _M_value += __r; 1386 return *this; 1387 } 1388 1389 inline complex<long double>& 1390 complex<long double>::operator-=(long double __r) 1391 { 1392 __real__ _M_value -= __r; 1393 return *this; 1394 } 1395 1396 inline complex<long double>& 1397 complex<long double>::operator*=(long double __r) 1398 { 1399 _M_value *= __r; 1400 return *this; 1401 } 1402 1403 inline complex<long double>& 1404 complex<long double>::operator/=(long double __r) 1405 { 1406 _M_value /= __r; 1407 return *this; 1408 } 1409 1410 template<typename _Tp> 1411 inline complex<long double>& 1412 complex<long double>::operator=(const complex<_Tp>& __z) 1413 { 1414 __real__ _M_value = __z.real(); 1415 __imag__ _M_value = __z.imag(); 1416 return *this; 1417 } 1418 1419 template<typename _Tp> 1420 inline complex<long double>& 1421 complex<long double>::operator+=(const complex<_Tp>& __z) 1422 { 1423 __real__ _M_value += __z.real(); 1424 __imag__ _M_value += __z.imag(); 1425 return *this; 1426 } 1427 1428 template<typename _Tp> 1429 inline complex<long double>& 1430 complex<long double>::operator-=(const complex<_Tp>& __z) 1431 { 1432 __real__ _M_value -= __z.real(); 1433 __imag__ _M_value -= __z.imag(); 1434 return *this; 1435 } 1436 1437 template<typename _Tp> 1438 inline complex<long double>& 1439 complex<long double>::operator*=(const complex<_Tp>& __z) 1440 { 1441 _ComplexT __t; 1442 __real__ __t = __z.real(); 1443 __imag__ __t = __z.imag(); 1444 _M_value *= __t; 1445 return *this; 1446 } 1447 1448 template<typename _Tp> 1449 inline complex<long double>& 1450 complex<long double>::operator/=(const complex<_Tp>& __z) 1451 { 1452 _ComplexT __t; 1453 __real__ __t = __z.real(); 1454 __imag__ __t = __z.imag(); 1455 _M_value /= __t; 1456 return *this; 1457 } 1458 1459 // These bits have to be at the end of this file, so that the 1460 // specializations have all been defined. 1461 // ??? No, they have to be there because of compiler limitation at 1462 // inlining. It suffices that class specializations be defined. 1463 inline 1464 complex<float>::complex(const complex<double>& __z) 1465 : _M_value(__z.__rep()) { } 1466 1467 inline 1468 complex<float>::complex(const complex<long double>& __z) 1469 : _M_value(__z.__rep()) { } 1470 1471 inline 1472 complex<double>::complex(const complex<float>& __z) 1473 : _M_value(__z.__rep()) { } 1474 1475 inline 1476 complex<double>::complex(const complex<long double>& __z) 1477 : _M_value(__z.__rep()) { } 1478 1479 inline 1480 complex<long double>::complex(const complex<float>& __z) 1481 : _M_value(__z.__rep()) { } 1482 1483 inline 1484 complex<long double>::complex(const complex<double>& __z) 1485 : _M_value(__z.__rep()) { } 1486 1487_GLIBCXX_END_NAMESPACE 1488 1489#endif /* _GLIBCXX_COMPLEX */ 1490