xref: /original-bsd/usr.bin/f77/libF77/pow_zi.c (revision 7921151c)
1 /*-
2  * Copyright (c) 1980 The Regents of the University of California.
3  * All rights reserved.
4  *
5  * %sccs.include.proprietary.c%
6  */
7 
8 #ifndef lint
9 static char sccsid[] = "@(#)pow_zi.c	5.3 (Berkeley) 04/12/91";
10 #endif /* not lint */
11 
12 #include "complex"
13 
14 #define	Z_MULEQ(A,B)	\
15 	t = (A).dreal * (B).dreal - (A).dimag * (B).dimag,\
16 	(A).dimag = (A).dreal * (B).dimag + (A).dimag * (B).dreal,\
17 	(A).dreal = t	/* A *= B */
18 
19 void
20 pow_zi(p, a, b) 	/* p = a**b  */
21 	dcomplex *p, *a;
22 	long int *b;
23 {
24 	register long n = *b;
25 	double t;
26 	dcomplex x;
27 
28 	x = *a;
29 	p->dreal = (double)1, p->dimag = (double)0;
30 	if (!n)
31 		return;
32 	if (n < 0) {
33 		z_div(&x, p, a);
34 		n = -n;
35 	}
36 	while (!(n&1)) {
37 		Z_MULEQ(x, x);
38 		n >>= 1;
39 	}
40 	for (*p = x; --n > 0; Z_MULEQ(*p, x))
41 		while (!(n&1)) {
42 			Z_MULEQ(x, x);
43 			n >>= 1;
44 		}
45 }
46