libf2c/libF77/pow_zi.c

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