Bizarre math bug

[Daniel Corbier]

I ran across a strange problem in a DLL compiled with PB.
Step by step I removed code and modified it, until I simplified
it to the small Exe program below, so others could easily try it.
For some reason instead of correctly displaying 0, it displays
2.4426950… I don’t know if this is a PB problem, or if it
happens only on certain computers (I’m using a PIII). I hope
others can check this out to let me know their findings.

$Compile Exe

Function WinMain (ByVal hCurInstance  As Long, _
                  ByVal hPrevInstance As Long, _
                  lpszCmdLine         As Asciiz Ptr, _
                  ByVal nCmdShow      As Long) Export As Long

   Dim a As Long, b As Ext, c As Ext, d As Ext

   MsgBox Str$(-1 + Abs((-1) ^ 2) / 1)
                                                    
   a = c
   
End Function

[Tim Wisseman]

I am seeing the very same problem on my P200.

If I add the line:
#Register None

The code works correctly and returns a zero. Very strange.



------------------

[Florent Heyworth]

I can confirm the problem tested on my box (NT 4.0 PII). It
seems to be a parsing error since the following give the
correct answers

FUNCTION WINMAIN (BYVAL hCurInstance  AS LONG, _
                  BYVAL hPrevInstance AS LONG, _
                  lpszCmdLine         AS ASCIIZ PTR, _
                  BYVAL nCmdShow      AS LONG) EXPORT AS LONG

   DIM a AS LONG, b AS EXT, c AS EXT, d AS EXT
   
   'PARSING ERROR - NOT OKAY
   MSGBOX STR$(-1 + ABS((-1)^2) / 1)

   'RESULTS ARE FINE
   MSGBOX STR$(-1 + ABS(-1^2) / 1)
   MSGBOX STR$(-1 + ABS(SQR(-1)) / 1)

END FUNCTION

Cheers

Florent

[UPDATE]
Wrote this while Tim posted his message. Tested with #REGISTER NONE
and got a correct answer. Hmmm...


------------------


[This message has been edited by Florent Heyworth (edited August 12, 2000).]

[Florent Heyworth]

Note that this problem only occurs when declaring
the unused variables and not using #REGISTER NONE

DIM a AS LONG, b AS EXT, c AS EXT, d AS EXT
or
DIM a AS INTEGER, b AS EXT, c AS EXT, d AS EXT
etc..

Doing
DIM a AS EXT, b AS EXT, c AS EXT, d AS EXT

yields the correct answer with or without #REGISTER NONE

Cheers

Florent


------------------

[Bob Houle]

daniel,

this appears to be related to the problem i detected last june, and never did get a satisfactory answer. see

in this case unused variables were causing the program to 'speed up'.

this might help you,
--bob

------------------

[Tim Wisseman]

This code also fails
Now we know the STR$() function and the ABS() function is
not the problem.

#COMPILE EXE
DECLARE SUB TEST()
FUNCTION WINMAIN (BYVAL hCurInstance  AS LONG, _
                  BYVAL hPrevInstance AS LONG, _
                  lpszCmdLine         AS ASCIIZ PTR, _
                  BYVAL nCmdShow      AS LONG) EXPORT AS LONG

    test
END FUNCTION

SUB Test()
   DIM a AS LONG, b AS EXT, c AS EXT, d AS EXT
   DIM z AS SINGLE
   z = -1 + ((-1) ^ 2) / 1
   MSGBOX STR$(z)
   a = c
END SUB

The (-1) value can be replaced with any negative value
and it will still fail.

[This message has been edited by Tim Wisseman (edited August 13, 2000).]

[Tim Wisseman]

This code also fails.

#COMPILE EXE
FUNCTION WINMAIN (BYVAL hCurInstance  AS LONG, _
                  BYVAL hPrevInstance AS LONG, _
                  lpszCmdLine         AS ASCIIZ PTR, _
                  BYVAL nCmdShow      AS LONG) EXPORT AS LONG
                  
   DIM b AS EXT, c AS EXT, d AS EXT
   DIM z AS SINGLE
   DIM r AS SINGLE
   DIM r2 AS SINGLE
   r = -1
   r2 = 1
   z = -1 + (r ^ 2) / r2
   MSGBOX STR$(z)
   b = 0
   c = 0
   d = 0
   b = c + d
END FUNCTION

------------------

[Tim Wisseman]

This code also fails:

z = (((-1)^ 2) * 1) - 1

------------------

[Charles Dietz]

I wonder if the problem has something to do with setting 32 bit registers
with an 80 bit value. It seems that if the the second declared variable
is EXT the problem occurs, but if it is declared as DOUBLE, the right answer
is obtained. Maybe the 80 bits is somehow spilling over into other registers.

#COMPILE EXE
FUNCTION PBMAIN
   'DIM a AS LONG, b AS EXT, c AS EXT, d AS EXT    'a = 2
   DIM a AS LONG, b AS DOUBLE, c AS EXT, d AS EXT  'a = 0  
   a = -1 + ABS((-1) ^ 2) / 1
   MSGBOX STR$(a)
END FUNCTION

------------------

[Tim Wisseman]

The working code:
  z = (((3)^ 2) * 4) - 1
  Compiles into:

004010aa   fld1
004010ac   fild      word ptr [004022b0]
004010b2   fild      word ptr [004022ac]
004010b8   fild      word ptr [004022a8]
004010be   fxch      st(1)
004010c0   call      004016f5
004010c5   fmulp     st(1),st
004010c7   fsubrp    st(1),st
004010c9   fstp      dword ptr [ebp-000000a0]
004010cf   fld       dword ptr [ebp-000000a0]

  This translates into:

;FPU Stack setup
004010aa   fld1   ;Push 1 onto FPU Stack
004010ac   fild      word ptr [004022b0] ;Push the value 4
004010b2   fild      word ptr [004022ac] ;Push the value 2
004010b8   fild      word ptr [004022a8] ;Push the value 3 
004010be   fxch      st(1)  ;Exchange ST(1) with ST(0)

   ;Here is our FPU Stack
 ST(0) = +2.00000000000000000e+0000
 ST(1) = +3.00000000000000000e+0000
 ST(2) = +4.00000000000000000e+0000
 ST(3) = +1.00000000000000000e+0000
 ST(4) = +0.00000000000000000e+0000
 ST(5) = +0.00000000000000000e+0000
 ST(6) = +0.00000000000000000e+0000
 ST(7) = +0.00000000000000000e+0000

004010c0   call      004016f5 ;<<< Raise ST(1) to the ST(0) Power! 
                              ;(3) ^ 2

   ;The stack after the call
 ST(0) = +9.00000000000000145e+0000  <--There is the value we want
 ST(1) = +4.00000000000000000e+0000
 ST(2) = +1.00000000000000000e+0000
 ST(3) = +0.00000000000000000e+0000
 ST(4) = +0.00000000000000000e+0000
 ST(5) = 1#IND                     
 ST(6) = +3.00000000000000000e+0000
 ST(7) = +9.00000000000000145e+0000

004010c5   fmulp     st(1),st  ;ST(0) * ST(1), POP the stack

004010c7   fsubrp    st(1),st  ;ST(0) - ST(1), POP the stack

   ;FPU stack after the * 4 and - 1 opps
 ST(0) = +3.50000000000000058e+0001  ;<<  35! Correct!
 ST(1) = +0.00000000000000000e+0000
 ST(2) = +0.00000000000000000e+0000
 ST(3) = 1#IND                     
 ST(4) = +3.00000000000000000e+0000
 ST(5) = +9.00000000000000145e+0000
 ST(6) = +9.00000000000000145e+0000
 ST(7) = +3.60000000000000058e+0001

004010c9   fstp      dword ptr [ebp-000000a0] ;Store ST(0) at z and pop FPU Stack

 z = +3.50000000000000058e+0001 '(CORRECT!!!)

004010cf   fld       dword ptr [ebp-000000a0]

/////////////////////////////////
/////////////////////////////////


The broken code:
  z = (((-3)^ 2) * 4) - 1

  Produces: 

;FPU Stack setup
004010aa   fld1   ;Push 1 onto FPU Stack
004010ac   fild      word ptr [004022b0] ;Push the value 4
004010b2   fild      word ptr [004022ac] ;Push the value 2
004010b8   fild      word ptr [004022a8] ;Push the value 3 
004010be   fchs             ;Change the sign of ST(0)
004010c0   fxch      st(1)  ;Exchange ST(1) with ST(0)

 FPU STACK
 ST(0) = +2.00000000000000000e+0000
 ST(1) = -3.00000000000000000e+0000
 ST(2) = +4.00000000000000000e+0000
 ST(3) = +1.00000000000000000e+0000
 ST(4) = +0.00000000000000000e+0000
 ST(5) = +0.00000000000000000e+0000
 ST(6) = +0.00000000000000000e+0000
 ST(7) = +0.00000000000000000e+0000

004010c2   call      004016f5  ;<<<Call Raise ST(1) to the ST(0) Power! 
                               ; (-3) ^ 2

 FPU STACK after the call
 ST(0) = +1.00000000000000000e+0000 << ERROR!!!
 ST(1) = +6.93147180559945309e-0001 << ERROR!!!
 ST(2) = +4.00000000000000000e+0000 << WRONG PLACE ON STACK!
 ST(3) = +1.00000000000000000e+0000 << WRONG PLACE ON STACK!
 ST(4) = +0.00000000000000000e+0000
 ST(5) = +0.00000000000000000e+0000
 ST(6) = 1#IND                     
 ST(7) = +1.44269504088896340e+0000

 The returning FPU stack not only has the wrong value on top, there is an extra value on
 the stack. This messes up the next to FPU math opperations that expect a 4 at ST(1) and a 1 at ST(2).

------------------

[Lance Edmonds]

Thanks guys... R&D tell me that this problem has been identified and will be sorted out in the next update to the compiler.


------------------
Lance
PowerBASIC Support
mailto:[email protected][email protected]</A>

[Daniel Corbier]

while on this subject, there's a message that was posted in the
faq section (at http://www.powerbasic.com/support/pb...hread.php?t=45
about the requirement of some assembly code so that extended
precision would work properly. i was never able to duplicate that
problem on my computer. so i haven't used the assembly code.
lance, can you give an example of code where that problem would
show up without the fix? also, will this be addressed in the
next version? the archives in these forums do not seem to go
prior to 1999. i may be wrong, but somehow, i seem to remember
reading about that asm fix well before version 6.0 was out.

also, although i couldn't verify the problem described in the
faq, there's another extended-precision problem that a few
customers who use my dll did notice. at least i think it's
different from the faq problem, because it doesn't do
calculations using 54-bit float. instead, on certain
computers it just crashes, with the error "invalid
floating-point operation". after searching around, i found
that apparently certain versions of particular popular
programs turn off the 80-bit mode. the remedy for c++ users
of my dll who encountered this problem was to add these
lines:

#include <float.h>
_fpreset();
_control87(mcw_em, mcw_em);

i don't know if the asm code from the faq solves this or not.
in any event, it would be nice for the next version to handle
things so that no asm code or external function call would be
required just to make extended precision work properly.

------------------
daniel corbier
ucalc fast math parser
http://www.ucalc.com

[Lance Edmonds]

I'm not sure if I have any example code to test if the FPU is running in 54-bit precision mode, which I must stress is a behavior of Windows not PowerBASIC! I'll check when I get back to my DEV PC and post some if I still have any in my archives.

Regardless, R&D tell me that the next update to PB/CC and PB/DLL will automatically switch the FPU to 80-bit precision for you.


------------------
Lance
PowerBASIC Support
mailto:[email protected][email protected]</A>

[John Hutchins]

Is the problem not with raising a negative number to a power?
This is how i would have done it.
MSGBOX STR$ (-1 +((ABS(-1) ^2)/1))
so that the negation is removed before raising to the power of 2.
Answer 0

------------------

[Hanns Ackermann]

Daniel,

I'm very interested in that ASM fix - Is the code available
somewhere? Unfortunately, the URL given above yields HTTP
Error 404(Not found) ...

Thanks & Greetings
Hanns.

------------------


[This message has been edited by Hanns Ackermann (edited August 21, 2000).]

[Lance Edmonds]

it is right there in the faq forum:
http://www.powerbasic.com/support/pb...hread.php?t=45

------------------
lance
powerbasic support
mailto:[email protected][email protected]</a>

[Tim Wisseman]

>>
Is the problem not with raising a negative number to a power?
This is how i would have done it.
MSGBOX STR$ (-1 +((ABS(-1) ^2)/1))
so that the negation is removed before raising to the power of 2.
Answer 0
>>

My tests show that the problem only hits when you raise
a negative number to a power. Using ABS() is a good work
around to the problem. Your code will work fine.

Tim


------------------

[Hanns Ackermann]

Thanks Lance,

seems that the right bracket ")" is erronously included in Daniel's URL.

Cheers, Hanns.

------------------

[John Hutchins]

On all the calculators I have ever used, an error or overflow is always
returned when you try to raise a negative number to a power.
If I remember correctly, this involves complex numbers.

------------------

[Hanns Ackermann]

John,

because of a=exp(log(a)), a=x^n=exp(n*log(x)) yields an error
if x<=0. I think the calculators are using something like that,
and of course, f.ex., (-2)^3=exp(3*log(-2)) fails. That is why
the formula must work not only for integers, but for real numbers,
too. Of course, f.ex., (-2)^3=(-2)*(-2)*(-2)=-8. That's not the
problem with complex numbers, but you need complex numbers if n
is not an integer.

In PB6, I tried >> msgbox str$( (-2)^3.45 ) << and, without an
error, got the result: -8, (2)^3.45=10.92 is ok. Evidently, PB6
uses the next integer for computation, perhaps a good (?) idea.

Daniel,

in your original post you computed (-1)^2, in terms of the above
formula this reads as exp(2*log(-1)), and of course log(-1) fails.
Substituting (-1)^2 with exp(2*log(-1))in your program yields the
same wrong result (without an error!), inserting (-1)*(-1)
correctly yields zero. Did John found the key to a solution?

Regards, Hanns.


------------------




[This message has been edited by Hanns Ackermann (edited August 23, 2000).]