1#include "f2c.h"
2
3#ifdef KR_headers
4VOID pow_zi(p, a, b) 	/* p = a**b  */
5 doublecomplex *p, *a; integer *b;
6#else
7extern void z_div(doublecomplex*, doublecomplex*, doublecomplex*);
8void pow_zi(doublecomplex *p, doublecomplex *a, integer *b) 	/* p = a**b  */
9#endif
10{
11	integer n;
12	unsigned long u;
13	double t;
14	doublecomplex q, x;
15	static doublecomplex one = {1.0, 0.0};
16
17	n = *b;
18	q.r = 1;
19	q.i = 0;
20
21	if(n == 0)
22		goto done;
23	if(n < 0)
24		{
25		n = -n;
26		z_div(&x, &one, a);
27		}
28	else
29		{
30		x.r = a->r;
31		x.i = a->i;
32		}
33
34	for(u = n; ; )
35		{
36		if(u & 01)
37			{
38			t = q.r * x.r - q.i * x.i;
39			q.i = q.r * x.i + q.i * x.r;
40			q.r = t;
41			}
42		if(u >>= 1)
43			{
44			t = x.r * x.r - x.i * x.i;
45			x.i = 2 * x.r * x.i;
46			x.r = t;
47			}
48		else
49			break;
50		}
51 done:
52	p->i = q.i;
53	p->r = q.r;
54	}
55