home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 7 / 07.iso / c / c082_144 / 2.ddi / CLIBSRC3.ZIP / CPLX2.CPP < prev    next >
Encoding:
C/C++ Source or Header  |  1992-06-10  |  5.6 KB  |  220 lines
  1. /*-----------------------------------------------------------------------*
  2.  * filename - cplx2.cpp
  3.  * C++ Complex Library Routines
  4.  *-----------------------------------------------------------------------*/
  5.  
  6. /*
  7.  *      C/C++ Run Time Library - Version 5.0
  8.  *
  9.  *      Copyright (c) 1990, 1992 by Borland International
  10.  *      All Rights Reserved.
  11.  *
  12.  */
  13.  
  14.  
  15. #include <complex.h>
  16.  
  17. static const complex _Cdecl complex_zero(0, 0);
  18. inline long double sqr(double x) { return x*x;}
  19.  
  20. #pragma warn -lvc
  21.  
  22. complex _FAR & _Cdecl complex::operator/=(complex _FAR & z2)
  23. {
  24.     long double sum_sqrs = norm(z2);
  25.     *this *= conj(z2);
  26.     re /= sum_sqrs;
  27.     im /= sum_sqrs;
  28.     return *this;
  29. }
  30.  
  31. complex _Cdecl _FARFUNC operator/(complex _FAR & z1, complex _FAR & z2)
  32. {
  33.     return z1 * conj(z2) / norm(z2);
  34. }
  35.  
  36. complex _Cdecl _FARFUNC operator/(double re_val1, complex _FAR & z2)
  37. {
  38.     return re_val1 * conj(z2) / norm(z2);
  39. }
  40.  
  41. complex _Cdecl _FARFUNC acos(complex _FAR & z)
  42. {
  43.     // complex i(0,1);
  44.     // return -i*log(z+i*sqrt(1-z*z));
  45.     complex temp1(1 - sqr(z.re) + sqr(z.im), -2*z.re*z.im);
  46.     double phi = arg(temp1)/2;
  47.     double rp = sqrt(abs(temp1));
  48.     complex temp2(z.re - rp*sin(phi), z.im + rp*cos(phi));
  49.     return complex(arg(temp2), -log(abs(temp2)));
  50. }
  51.  
  52. double _Cdecl _FARFUNC arg(complex _FAR & z)
  53. {
  54.     return z == complex_zero ? 0.0 : atan2(z.im, z.re);
  55. }
  56.  
  57. complex _Cdecl _FARFUNC asin(complex _FAR & z)
  58. {
  59.     // complex i(0,1);
  60.     // return -i*log(i*z+sqrt(1-z*z));
  61.  
  62.     complex temp1(1 - sqr(z.re) + sqr(z.im), -2*z.re*z.im);
  63.     double phi = arg(temp1) / 2;
  64.     double rp = sqrt(abs(temp1));
  65.     complex temp2(-z.im + rp*cos(phi), z.re + rp*sin(phi));
  66.     return complex(arg(temp2), -log(abs(temp2)));
  67. }
  68.  
  69. complex _Cdecl _FARFUNC atan(complex _FAR & z)
  70. {
  71.     // complex i(0, 1);
  72.     // return -0.5*i * log((1+i*z)/(1-i*z));
  73.     // if z=i then a floating-point exception may occur
  74.  
  75.     double opb = 1 + z.im;
  76.     double a2 = sqr(z.re);
  77.     double den = opb*opb + a2;
  78.     complex temp(((1-z.im)*opb - a2)/den, 2*z.re/den);
  79.     return complex(arg(temp)/2, -log(abs(temp))/2);
  80. }
  81.  
  82. complex _Cdecl _FARFUNC cos(complex _FAR & z)
  83. {
  84.     // complex i(0, 1);
  85.     // return (exp(i*z) + exp(-i*z))/2;
  86.  
  87.     long double eb  = exp(z.im);
  88.     long double emb = 1 / eb;
  89.     return complex(cos(z.re)*(emb+eb)/2, sin(z.re)*(emb-eb)/2);
  90. }
  91.  
  92. complex _Cdecl _FARFUNC cosh(complex _FAR & z)
  93. {
  94.     // return (exp(z) + exp(-z))/2;
  95.     long double ea = exp(z.re);
  96.     long double eainv = 1 / ea;
  97.     return complex(cos(z.im)*(ea + eainv)/2, sin(z.im)*(ea - eainv)/2);
  98. }
  99.  
  100. complex _Cdecl _FARFUNC exp(complex _FAR & z)
  101. {
  102.     return polar(exp(z.re),z.im);
  103. }
  104.  
  105. complex _Cdecl _FARFUNC log(complex _FAR & z)
  106. {
  107.     // exception if z = 0
  108.     return complex(log(abs(z)), arg(z));
  109. }
  110.  
  111. complex _Cdecl _FARFUNC log10(complex _FAR & z)
  112. {
  113.     // return log10e*log(z);
  114.     return log(z) * M_LOG10E;
  115. }
  116.  
  117. complex _Cdecl _FARFUNC pow(complex _FAR & base, double expon)
  118. {
  119.     // return exp(expon*log(base));
  120.  
  121.     if (base == complex_zero && expon > 0) return complex_zero;
  122.     return polar(pow(abs(base), expon), expon*arg(base));
  123. }
  124.  
  125. complex _Cdecl _FARFUNC pow(double base, complex _FAR & expon)
  126. {
  127.     // return exp(expon*log(base));
  128.     if (base == 0 && real(expon) > 0) return complex_zero;
  129.     double lnx = log(fabs(base));
  130.     if (base > 0.0)
  131.         return exp(expon * lnx);
  132.     return exp(expon * complex(lnx,M_PI));
  133. }
  134.  
  135. complex _Cdecl _FARFUNC pow(complex _FAR & base, complex _FAR & expon)
  136. {
  137.     if (base == complex_zero && real(expon) > 0) return complex_zero;
  138.     return exp(expon * log(base));
  139. }
  140.  
  141. complex _Cdecl _FARFUNC sin(complex _FAR & z)
  142. {
  143.     // complex i(0,1);
  144.     // return (exp(i*z) - exp(-i*z))/(2*i);
  145.  
  146.     long double eb  = exp(z.im);
  147.     long double emb = 1 / eb;
  148.     return complex(sin(z.re)*(emb+eb)/2, -0.5*cos(z.re)*(emb-eb));
  149. }
  150.  
  151. complex _Cdecl _FARFUNC sinh(complex _FAR & z)
  152. {
  153.     // return (exp(z) - exp(-z))/2
  154.     long double ea = exp(z.re);
  155.     long double eainv = 1 / ea;
  156.     return complex(cos(z.im)*(ea - eainv)/2, sin(z.im)*(ea + eainv)/2);
  157. }
  158.  
  159. complex _Cdecl _FARFUNC sqrt(complex _FAR & z)
  160. {
  161.     // return pow(z, 0.5);
  162.     return polar(sqrt(abs(z)), arg(z) / 2);
  163. }
  164.  
  165. complex _Cdecl _FARFUNC tan(complex _FAR & z)
  166. {
  167.     // return sin(z)/cos(z)
  168.     double sina = sin(z.re);
  169.     double cosa = cos(z.re);
  170.     long double eb = exp(z.im);
  171.     long double emb = 1 / eb;
  172.     double emin = emb - eb;
  173.     double eplus = emb + eb;
  174.     return complex(4*sina*cosa, -emin*eplus)
  175.         / (sqr(cosa*eplus) + sqr(sina*emin));
  176. }
  177.  
  178. complex _Cdecl _FARFUNC tanh(complex _FAR & z)
  179. {
  180.     // return sinh(z)/cosh(z);
  181.     double sinb = sin(z.im);
  182.     double cosb = cos(z.im);
  183.     long double ea = exp(z.re);
  184.     long double eainv = 1 / ea;
  185.     double eamin = ea - eainv;
  186.     double eapls = ea + eainv;
  187.     return complex(eamin*eapls, 4*sinb*cosb)
  188.         / (sqr(cosb*eapls) + sqr(sinb*eamin));
  189. }
  190.  
  191.  
  192. // Stream I/O function definitions
  193.  
  194. ostream _FAR & _Cdecl _FARFUNC operator<<(ostream _FAR & s, complex _FAR & z)
  195. {
  196.     return s << "(" << real(z) << ", " << imag(z) << ")";
  197. }
  198.  
  199. istream _FAR & _Cdecl _FARFUNC operator>>(istream _FAR & s, complex _FAR & z)
  200. {
  201.     double re_val = 0.0, im_val = 0.0;
  202.     char c = 0;
  203.  
  204.     s >> c;
  205.     if (c == '(') {
  206.         s >> re_val >> c;
  207.          if (c == ',') s >> im_val >> c;
  208.          if (c != ')') s.clear(ios::badbit | s.rdstate());
  209.     }
  210.     else {
  211.         s.putback(c);
  212.         s >> re_val;
  213.     }
  214.     if (s) z = complex(re_val, im_val);
  215.     return s;
  216. }
  217.  
  218. #pragma warn .lvc
  219.  
  220.