#include "f2c.h" #ifdef __cplusplus extern "C" { #endif #ifdef KR_headers VOID pow_zi(p, a, b) /* p = a**b */ doublecomplex *p, *a; integer *b; #else extern void z_div(doublecomplex*, doublecomplex*, doublecomplex*); void pow_zi(doublecomplex *p, doublecomplex *a, integer *b) /* p = a**b */ #endif { 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; } #ifdef __cplusplus } #endif