47 lines
760 B
C++
47 lines
760 B
C++
|
// power function for double precision floating point
|
||
|
|
|
||
|
|
#include "stdafx.h"
|
||
|
|
#include "defs.h"
|
||
|
|
|
||
|
|
void
|
||
|
|
dpow(void)
|
||
|
|
{
|
||
|
|
double a, b, base, expo, result, theta;
|
||
|
|
|
||
|
|
expo = pop_double();
|
||
|
|
base = pop_double();
|
||
|
|
|
||
|
|
// divide by zero?
|
||
|
|
|
||
|
|
if (base == 0.0 && expo < 0.0)
|
||
|
|
stop("divide by zero");
|
||
|
|
|
||
|
|
// nonnegative base or integer power?
|
||
|
|
|
||
|
|
if (base >= 0.0 || fmod(expo, 1.0) == 0.0) {
|
||
|
|
result = pow(base, expo);
|
||
|
|
push_double(result);
|
||
|
|
return;
|
||
|
|
}
|
||
|
|
|
||
|
|
result = pow(fabs(base), expo);
|
||
|
|
|
||
|
|
theta = M_PI * expo;
|
||
|
|
|
||
|
|
// this ensures the real part is 0.0 instead of a tiny fraction
|
||
|
|
|
||
|
|
if (fmod(expo, 0.5) == 0.0) {
|
||
|
|
a = 0.0;
|
||
|
|
b = sin(theta);
|
||
|
|
} else {
|
||
|
|
a = cos(theta);
|
||
|
|
b = sin(theta);
|
||
|
|
}
|
||
|
|
|
||
|
|
push_double(a * result);
|
||
|
|
push_double(b * result);
|
||
|
|
push(imaginaryunit);
|
||
|
|
multiply();
|
||
|
|
add();
|
||
|
|
}
|