/* Copyright (C) 1996 by Thomas Glen Smith. All Rights Reserved. */ /* powerpx APL2 V1.0.0 ************************************************ * Called from powerp when *left is negative, and the denominator of the* * fraction representing *rite (the exponent) is even. The result can * * be expressed by ((|*left)**rite)X-1**rite. * ***********************************************************************/ #define INCLUDES MATH #include "includes.h" #define ABS(v) ((v) > 0e0 ? (v) : -(v)) void powerpx(left,rite,pret,ret,num,den) double *left, /* Real number to be taken to the power in *rite. */ *rite, /* Exponent, or power, to which *left is taken. */ *pret, /* Place to store the answer. */ *ret; /* *ret == (|*left)**rite, *(ret+1) == 0e0. */ double num,den; /* Numerator and denominator of exponent. */ { Timesx; extern double fuzz; double angle,wrka[2],wrkb[2]; for (;;) { /* Lets me use break. */ if (fmod(num,2.0) < fuzz) break; /* Num. even, -1**rite == 1. */ num = 1e0; /* Numerator might as well be 1. */ den /= 2e0; /* Make it for (-1*.5)*1%den. */ angle = PI / 2 * den; *wrka = cos(angle); *(wrka+1) = sin(angle); *wrkb = *ret; *(wrkb+1) = *(ret+1); timesx(wrka,wrkb,ret); break; /* Final break from for(;;) */ } *pret = *ret; *(pret+1) = *(ret+1); }