C pow()

[jj2007]

I'm working on a better version of ExpXY(a, b):

Code Select Expand

  For_ fct=1.7 To 2.3 Step 0.1
Print Str$("\n-1.23 ** %3f", fct), Str$("\tMB: %4f  ", ExpXY2(-1.23, fct)v)
invoke crt_pow, FP8(-1.23), fct
Print Str$("\t C pow: %4f", ST(0)v)
  Next

Output:

Code Select Expand

-1.23 ** 1.70   MB: -1.422       C pow: 0.0
-1.23 ** 1.80   MB: -1.452       C pow: 0.0
-1.23 ** 1.90   MB: -1.482       C pow: 0.0
-1.23 ** 2.00   MB: -1.513       C pow: 1.513
-1.23 ** 2.10   MB: -1.545       C pow: 0.0
-1.23 ** 2.20   MB: -1.577       C pow: 0.0
-1.23 ** 2.30   MB: -1.610       C pow: 0.0

The Windows calculator has no problems with the negative input, and comes to the same results as MB. The docs on pow() are remarkably silent on that point, e.g. https://www.tutorialspoint.com/c_standard_library/c_function_pow.htm

I vaguely remember an old math rule that a negative number squared yields a positive result. But it looks rather inconsistent, and the Windows calculator confuses me even more:
-1.23^2 = -1.5129 using x^y
-1.23^2 = +1.5129 using x²

[Biterider]

Hi JJ
On my Calc.exe (10.1906.53.0 - Win10) the result using x^y is positive, which is mathematically correct.

Biterider

[TimoVJL]

ISO/IEC 9899:TC2 Committee Draft — May 6, 2005 WG14/N1124

7.12.7.4 The pow functions
Synopsis
1 #include <math.h>
double pow(double x, double y);
float powf(float x, float y);
long double powl(long double x, long double y);
Description
2 The pow functions compute x raised to the power y. A domain error occurs if x is finite
and negative and y is finite and not an integer value. A range error may occur. A domain
error may occur if x is zero and y is zero. A domain error or range error may occur if x
is zero and y is less than zero.
Returns
3 The pow functions return x y

[jj2007]

On my Calc.exe (10.1906.53.0 - Win10) the result using x^y is positive, which is mathematically correct.

You mean this, I suppose: -1.23^2 = +1.5129 using x^y?
What is your result for -1.23^2.00000000001 then?

@Timo: Thanks. Which means pow() cannot handle -1.23^2.0000001 but -1.23^2 is ok.

I found an article on SOF that makes it crystal clear. Well, almost

And here is a test page by Wolfram Alpha:

-1.23^2.00000000000 ≈ 1.51290
-1.23^2.00000000001 ≈ 1.51290 + 0. 10-10*i
-1.23^2.50000000000 ≈ 1.67789*i

Same for Windows Calc.exe:
-1.23^2.00000000001 = -1.5129000000031319173686187117352
-1.23^2.50000000000 = -1.6778872680546807224467865795239

Isn't it fun?

[Biterider]

Hi JJ
Here is what I get:

-1.23^2 =  1.5129
-1.23^2.1 => invalid input
-1.23^2.00000000001 => invalid input

Biterider

[jj2007]

So Win7 and Win10 behave differently. See above results for the Wolfram Alpha calculator. I think I'll make negative inputs optional, and clarify what it does in the documentation.

[FORTRANS]

Hi,

   The results should be complex numbers,
that much I remember.  Here are the results
from my HP calculator.  Since I have forgotten
how to it with paper and pencil, I'll have to
believe them.

Code Select Expand

-1.23 ** 1.70   0.83571 -1.15026i
-1.23 ** 1.80   1.17432 -0.85319i
-1.23 ** 1.90   1.40937 -0.45793i
-1.23 ** 2.00   1.51290
-1.23 ** 2.10   1.46895  0.47729i
-1.23 ** 2.20   1.27570  0.92685i
-1.23 ** 2.30   0.94624  1.30239i

Cheers,

Steve N.

[jj2007]

Steve,

Your HP is in line with Wolfram Alpha's results

Windows 7 calc.exe produces something completely different, and it's in line with what the fpu spits out. And it seems they are both wrong, so I'll better drop that idea

[aw27]

We can do it easily in C99 (1999) standard C (no need for C++), so it is easily converted to True Masm MASM™ (left as an exercise).

Code Select Expand


#include <stdio.h>   
#include <complex.h>

int main()
{
_Dcomplex base = { -1.23, 0 };
_Dcomplex result;

double f;

for (f = 1.7; f < 2.31; f += 0.1)
{
_Dcomplex exponent = { f,0 };
result = cpow(base, exponent);
//printf("-1.23^%.2f = %f %fi\n", f, creal(result), cimag(result));
                printf("-1.23^%.2f = %f %fi\n", f, result._Val[0], result._Val[1]);
}
return 0;
}


Output:
-1.23^1.70 = 0.835713 -1.150261i
-1.23^1.80 = 1.174321 -0.853194i
-1.23^1.90 = 1.409373 -0.457933i
-1.23^2.00 = 1.512900 -0.000000i
-1.23^2.10 = 1.468950 0.477291i
-1.23^2.20 = 1.275701 0.926851i
-1.23^2.30 = 0.946238 1.302385i

[TimoVJL]

@Jose, that example was for msvc, not quite C99.
msvc have troubles with C99 and double complex ?

Code Select Expand

double complex cpow(double complex, double complex);
double creal(double complex z);
double cimag(double complex z);
int printf(const char * restrict format, ...);

int main(void)
{
double complex base = -1.23;

for (double f = 1.7; f < 2.31; f += 0.1)
{
double complex exponent = f;
double complex result = cpow(base, exponent);
printf("-1.23^%.2f = %f %fi\n", f, creal(result), cimag(result));
//printf("-1.23^%.2f = %f %fi\n", f, result._Val[0], result._Val[1]);
}
return 0;
}
how intel C compiler works with this code?
clang compiles this code, but it don't have it's own cpow(),creal(),cimag().

[aw27]

@Timo,

You are correct, it is Microsoft specific but appears to be more comprehensive than C99 because it allows to power complex bases to complex exponents.  (Nonsense).
https://docs.microsoft.com/en-us/cpp/c-runtime-library/complex-math-support?view=vs-2019
"The Microsoft implementation of the complex.h header defines these types as equivalents for the C99 standard native complex types"
The need of Microsoft to be original 

The code I posted builds with LLVM from within VS without change. Sure it borrowed the Microsoft stuff.

[mikeburr]

i daresay you already know this .as its elementary maths..
AW is correct but hasn't explained the following x^y does grow  it has a real component and an imaginary component often shown plotted as on an x-y axis
and there are 3 values to consider
the real part .. the imaginary part .. and the resulting  vector length [on  diagram rhs .. the blue line with arrow]
https://en.wikipedia.org/wiki/Complex_number
so for AW results
the vector length is calculated by squaring the the real and imaginary components .. and then taking the square root of the sum  - simple pythagorean geomerty
this is sometimes known as the absolute value
so
1.23^1.90 = 1.409373 -0.457933i the vector length = root [1.409373^2 + 0.457933^2 ]= root[2.196 ] ===>   1.4819
1.23^2.0  =  1.512900 + 0i               the vector length =root[ 1.512900 ^2 + 0.0^2] = root[2.2889] ===> 1.5129
1.23^2.10 = 1.468950 0.477291i   the vector length = root [1.468950^2 + 0.477291^2 ] = root[2.3856] ====>  1.5445
and so on
regards  mikeb

[UniverseIsASimulation]

My web-app claims that the assembly for "pow(x,y)" is:

Code Select Expand

finit
fld dword [x]
fld dword [y]
fxch
fld1
fxch
fyl2x
fldl2e
fdivp st1,st0
fmulp st1,st0
fldl2e
fmulp st1,st0
fld1
fscale
fxch
fld1
fxch
fprem
f2xm1
faddp st1,st0
fmulp st1,st0
fstp dword [result]
Apparently, it returns NaN for negative bases.