1 /* 2 * Copyright (c) 1980 Regents of the University of California. 3 * All rights reserved. The Berkeley software License Agreement 4 * specifies the terms and conditions for redistribution. 5 * 6 * @(#)pow_zi.c 5.2 01/24/88 7 */ 8 9 #include "complex" 10 11 #define Z_MULEQ(A,B) \ 12 t = (A).dreal * (B).dreal - (A).dimag * (B).dimag,\ 13 (A).dimag = (A).dreal * (B).dimag + (A).dimag * (B).dreal,\ 14 (A).dreal = t /* A *= B */ 15 16 void 17 pow_zi(p, a, b) /* p = a**b */ 18 dcomplex *p, *a; 19 long int *b; 20 { 21 register long n = *b; 22 double t; 23 dcomplex x; 24 25 x = *a; 26 p->dreal = (double)1, p->dimag = (double)0; 27 if (!n) 28 return; 29 if (n < 0) { 30 z_div(&x, p, a); 31 n = -n; 32 } 33 while (!(n&1)) { 34 Z_MULEQ(x, x); 35 n >>= 1; 36 } 37 for (*p = x; --n > 0; Z_MULEQ(*p, x)) 38 while (!(n&1)) { 39 Z_MULEQ(x, x); 40 n >>= 1; 41 } 42 } 43