I'm working on a program that solves for the roots of a quadratic equation. I was able to get the first root in the root1 subroutine. However, when I try to solve for the second root in root2, the "/2a" part of the quadratic formula keeps yielding NaN.
Here's the code:
INCLUDE Irvine32.inc
INCLUDE macros.inc
.data
a real8 ?
b real8 ?
cc real8 ?
a2 real8 ?
b2 real8 ?
cc2 real8 ?
two real8 2.0
four real8 4.0
two2 real8 2.0
four2 real8 4.0
ten real8 1000.0
num real8 10.0
.code
main PROC
finit
mWrite "Enter coefficient (a): "
call ReadFloat
fst a
fstp a2
mWrite "Enter coefficient (b): "
call ReadFloat
fst b
fstp b2
mWrite "Enter coefficient (c): "
call ReadFloat
fst cc
fstp cc2
mWrite "Roots: "
call root1
call Crlf
call root2
;call showfpustack
exit
main ENDP
root1 PROC
; b^2
fld b
fmul b
fchs ; flip sign
fst b
; 4 * a * c
fld four
fmul a
fmul cc
fchs
fsub b
fsqrt
fst four
fld b
fchs
fsqrt
fchs
fadd four
fst b
fld two
fmul a
fst two
fld b
fdiv two
call WriteFloat
call showfpustack
ret
root1 endp
root2 PROC
fld b2
fmul b2
fchs
fst b2
fld four2
fmul a2
fmul cc2
fchs
fsub b2
fsqrt
fst four2
fld b2
fchs
fsqrt
fchs
fsub four2
fst b2
call Crlf
call WriteFloat
fld two2
fmul a2
fst two2
fld b2
fdiv two2
call showfpustack
ret
root2 endp
end main
I was able to verify the results of previous calculations. It's just this part that I'm having trouble with.
Instead of speculating about why you get a NAN, I'll walk you through the code that you need to write in order to calculate the roots of the quadratic equation.
You see that you don't need any constants in memory nor do you need multiple copies of the inputs a, b, and c.
Reset the FPU environment
fninit
Store an interesting constant that we will re-use muliple times:
fld a ; (st0) a
fadd st0 ; (st0) a + a == 2a
Calculating D = b^2 - 4ac
fld b ; (st0) b (st1) 2a
fmul st0 ; (st0) b * b == b^2 (st1) 2a
fld c ; (st0) c (st1) b^2 (st2) 2a
fmul st2 ; (st0) c * 2a == 2ac (st1) b^2 (st2) 2a
fadd st0 ; (st0) 2ac + 2ac == 4ac (st1) b^2 (st2) 2a
fsubp ; (st0) b^2 - 4ac == D (st1) 2a
Here we need to test if D is not negative because no roots exist for D<0.
ftst ; This compares st0 to 0.0 and sets flags C3, C2, and C0
fnstsw ax ; Copies those flags to AX
sahf ; Copies those flags to EFLAGS
jp IsUnordered ; This should not happen
jc IsNegative ; This is very possible
We can take the square root only if D is not negative.
fsqrt ; (st0) sqrt(D) (st1) 2a
Going for the 1st root R1:
fld b ; (st0) b (st1) sqrt(D) (st2) 2a
fchs ; (st0) -b (st1) sqrt(D) (st2) 2a
fsub st1 ; (st0) -b - sqrt(D) (st1) sqrt(D) (st2) 2a
fdiv st2 ; (st0) (-b - sqrt(D)) / 2a == R1 (st1) sqrt(D) (st2) 2a
Going for the 2nd root R2:
fld b ; (st0) b (st1) R1 (st2) sqrt(D) (st3) 2a
fchs ; (st0) -b (st1) R1 (st2) sqrt(D) (st3) 2a
fadd st2 ; (st0) -b - sqrt(D) (st1) R1 (st2) sqrt(D) (st3) 2a
fdiv st3 ; (st0) (-b - sqrt(D)) / 2a == R2 (st1) R1 (st2) sqrt(D) (st3) 2a
At this point the FPU register stack has only 4 entries:
st3 = 2a
st2 = sqrt(D)
st1 = R1
st0 = R2
Please note that nowhere in this code did I have to store anything back in memory. Working with the FPU is a matter of careful planning. Re-arranging an expression is often advantageous and re-using previously calculated values is always a win.