5 /*complex divide, defensive against overflow from
6 * * and /, but not from + and -
7 * assumes underflow yields 0.0
9 * (a + bi)/(c + di) = ((a + bd/c) + (b - ad/c)i)/(c + dd/c)
10 * (a + bi)/(c + di) = (b - ai)/(d - ci)
13 cdiv(double a, double b, double c, double d, double *u, double *v)
27 cmul(double c1, double c2, double d1, double d2, double *e1, double *e2)
34 csq(double c1, double c2, double *e1, double *e2)
40 /* complex square root
41 * assumes underflow yields 0.0
42 * uses these identities:
43 * sqrt(x+_iy) = sqrt(r(cos(t)+_isin(t))
44 * = sqrt(r)(cos(t/2)+_isin(t/2))
45 * cos(t/2) = sin(t)/2sin(t/2) = sqrt((1+cos(t)/2)
46 * sin(t/2) = sin(t)/2cos(t/2) = sqrt((1-cos(t)/2)
49 csqrt(double c1, double c2, double *e1, double *e2)
79 void cpow(double c1, double c2, double *d1, double *d2, double pwr)
81 double theta = pwr*atan2(c2,c1);
82 double r = pow(hypot(c1,c2), pwr);