#include "f2c.h" extern void z_div (doublecomplex *, doublecomplex *, doublecomplex *); void pow_zi (doublecomplex * p, doublecomplex * a, integer * b) /* p = a**b */ { integer n; unsigned long u; double t; doublecomplex q, x; static doublecomplex one = { 1.0, 0.0 }; n = *b; q.r = 1; q.i = 0; if (n == 0) goto done; if (n < 0) { n = -n; z_div (&x, &one, a); } else { x.r = a->r; x.i = a->i; } for (u = n;;) { if (u & 01) { t = q.r * x.r - q.i * x.i; q.i = q.r * x.i + q.i * x.r; q.r = t; } if (u >>= 1) { t = x.r * x.r - x.i * x.i; x.i = 2 * x.r * x.i; x.r = t; } else break; } done: p->i = q.i; p->r = q.r; }