@@ -192,6 +192,17 @@ return( (*x)>=0 ?
192
192
}
193
193
194
194
195
+ #ifdef KR_headers
196
+ double floor ();
197
+ integer i_nint (x ) real * x ;
198
+ #else
199
+ #undef abs
200
+ integer i_nint (real * x )
201
+ #endif
202
+ {
203
+ return (integer )(* x >= 0 ? floor (* x + .5 ) : - floor (.5 - * x ));
204
+ }
205
+
195
206
#ifdef KR_headers
196
207
double pow ();
197
208
double pow_dd (ap , bp ) doublereal * ap , * bp ;
@@ -272,6 +283,79 @@ if(n != 0)
272
283
}
273
284
return (pow );
274
285
}
286
+
287
+ #ifdef KR_headers
288
+ VOID pow_zi (p , a , b ) /* p = a**b */
289
+ doublecomplex * p , * a ; integer * b ;
290
+ #else
291
+ extern void z_div (doublecomplex * , doublecomplex * , doublecomplex * );
292
+ void pow_zi (doublecomplex * p , doublecomplex * a , integer * b ) /* p = a**b */
293
+ #endif
294
+ {
295
+ integer n ;
296
+ unsigned long u ;
297
+ double t ;
298
+ doublecomplex q , x ;
299
+ static doublecomplex one = {1.0 , 0.0 };
300
+
301
+ n = * b ;
302
+ q .r = 1 ;
303
+ q .i = 0 ;
304
+
305
+ if (n == 0 )
306
+ goto done ;
307
+ if (n < 0 )
308
+ {
309
+ n = - n ;
310
+ z_div (& x , & one , a );
311
+ }
312
+ else
313
+ {
314
+ x .r = a -> r ;
315
+ x .i = a -> i ;
316
+ }
317
+
318
+ for (u = n ; ; )
319
+ {
320
+ if (u & 01 )
321
+ {
322
+ t = q .r * x .r - q .i * x .i ;
323
+ q .i = q .r * x .i + q .i * x .r ;
324
+ q .r = t ;
325
+ }
326
+ if (u >>= 1 )
327
+ {
328
+ t = x .r * x .r - x .i * x .i ;
329
+ x .i = 2 * x .r * x .i ;
330
+ x .r = t ;
331
+ }
332
+ else
333
+ break ;
334
+ }
335
+ done :
336
+ p -> i = q .i ;
337
+ p -> r = q .r ;
338
+ }
339
+
340
+ #ifdef KR_headers
341
+ VOID pow_ci (p , a , b ) /* p = a**b */
342
+ complex * p , * a ; integer * b ;
343
+ #else
344
+ extern void pow_zi (doublecomplex * , doublecomplex * , integer * );
345
+ void pow_ci (complex * p , complex * a , integer * b ) /* p = a**b */
346
+ #endif
347
+ {
348
+ doublecomplex p1 , a1 ;
349
+
350
+ a1 .r = a -> r ;
351
+ a1 .i = a -> i ;
352
+
353
+ pow_zi (& p1 , & a1 , b );
354
+
355
+ p -> r = p1 .r ;
356
+ p -> i = p1 .i ;
357
+ }
358
+
275
359
/* Unless compiled with -DNO_OVERWRITE, this variant of s_cat allows the
276
360
* target of a concatenation to appear on its right-hand side (contrary
277
361
* to the Fortran 77 Standard, but in accordance with Fortran 90).
0 commit comments