Incorrect Results - Floating Point Issue?

[FortranBeginner]

Hi

I am having an issue with a very simple calculation that ends up having "incorrect" results.

VOLUME, LENGTH, WIDTH and THICK are declared as REAL

VOLUME = (LENGTH / 1000.) * (WIDTH / 1000.) * (THICK / 1000.)

When I use the following values:
length = 130 ; width=95 ; thick=14
The debugger shows 1.7289999E-04. Yet I am expecting 1.729E-04 as a result.

I tried declaring VOLUME as REAL*8, yet this did not help. So I tried to change the FloatingPoint settings (eventhough I really don't understand any of their meanings). So far I did not manage to get the "correct" result.

Does anyone know which FloatingPoint Settings (or other settings which I haven't thought about) I have to use to get the "correct" results?
Or are those "incorrect" result normal with Fortran?

Thanks for your Help

[mecej4]

Yet I am expecting 1.729E-04 as a result.

That expectation is based on performing calculations using decimal floating point calculation. Real*4 variables are stored in the form ±a ×2^b, where a and b are integers represented in 24 and 8 bits, respectively^*, and floating point arithmetic is done in binary IEEE. The number in question (1.729E-04) cannot be represented exactly in the number of bits available. In fact, even the simpler number 0.1 does not have an exact representation as a REAL*4 constant.

It will be very useful for you to read enough about floating point arithmetic that findings such as this no longer perturb you. Among many articles,  you may read http://docs.oracle.com/cd/E19957-01/800-7895/800-7895.pdf.

Here is an instructive example:

[fortran]

program f3
real x,y,z
data x,y,z/4.0,3.0,1.0/
print *,z-(x/y-1)*y
end[/fortran]

Compile and run the example, and compare the result with what you get using (i) manual calculation and (ii) one or more calculators.

I tried declaring VOLUME as REAL*8, yet this did not help.

It should make no difference, because the expression on the right hand side is of type REAL*4. Try declaring WIDTH, LENGTH and THICK as REAL*8.

___________________

^* The exponent is adjusted with a bias such that the mantissa lies between 1 and 2, and its leading bit is implicitly = 1, but these details are less important.

[Andrew_Smith]

Your example showed assignments of integers to single precision floating point variables. The integers will get converted in the assignment so when you change the declaration of the variables to *8,  the assigned values should be as precise as you might expect. However if you assign using:

length = 130.0 ; width=95.0; thick=14.0

You will then only have single precsion accuracy constants assigned to double precsion variables. Here is a better way to assign constants to double precision variables:

length = 130.0_8 ; width=95.0_8 ; thick=14.0_8

or the old fashioned:

length = 130.0D0; width=95.0D0 ; thick=14.0D0

Other languages have their own notations for precision of floating point constants.

[Les_Neilson]

Steve Lionel (as Dr Fortran) has published much useful information covering this and other points which you may find interesting.

see http://software.intel.com/en-us/forums/topic/275071 

in particular the "Perils of Real Numbers" parts 1 2 and 3

Les

[FortranBeginner]

Thanks both for your help.

Declaring all variables as REAL*8 solved the problem.
I now also know how to assign constants properly.

[jimdempseyatthecove]

FortranBeginner,

>>Declaring all variables as REAL*8 solved the problem.
It may have hid the problem.

The problem is some nice clean decimal fractions cannot be stored in a finite number of bits as a binary fraction. 0.1 is a good example:

[fortran]

real(4) :: r4
real(8) :: r8fromr4, r8
r4 = 0.1
r8fromr4 = r4
r8 = 0.1_8

Watch Window
R8       0.100000000000000 REAL(8)
R8FROMR4 0.100000001490116 REAL(8)
R4       0.1000000         REAL(4)

Looks good eh.

Now choose Hexadecimal display

R8       #3FB999999999999A REAL(8)
R8FROMR4 #3FB99999A0000000 REAL(8)
R4       #3DCCCCCD         REAL(4)
[/fortran]

0.1 produces a repeating binary fraction (similar to 1/3 producing a repeating decimal fraction)

The decimal display is showing you the rounded number
The Hexadecimal display discloses two things:
a) the ...9999A is indicitive of rounding a repeating fraction (same with ...CCCCD)
b) the R8FROMR4 shows the truncation effect of an assignment of a rounded real(4) number.

Jim Dempsey

[Les_Neilson]

I tried to post this earlier but it hasn't arrived yet.
You may find Steve's Dr Fortran articles helpful
see http://software.intel.com/en-us/forums/topic/275071
especially "The Perils of Real Numbers" parts 1,2, and 3
Les

[SergeyKostrov]

>>...The debugger shows 1.7289999E-04. Yet I am expecting 1.729E-04 as a result. There is some UI issue with several Visual Studios, like 2005, 2008, etc. I see such "mis-roundings" every day when debugging codes and verifying results of some computations. However, if actual binary results are different from what you're expecting then there is some problem and, of course, it needs to be investigated. I will also take a look at your test case.

[mecej4]

Sergey Kostrov wrote:
I see such "mis-roundings" every day

Why do you think that this is an instance of "mis-rounding"? Looking at the OP's example in the VS debugger, using REAL*4 for all variables, and asking the debugger to show values in hexadecimal, for the volume I see 0x39354C7D, which is closer to the exact result, 1.729E-4, than the next closest machine number, 0x39354C7E. The number of decimal digits shown by the debugger is consistent with the 24-bit precision of REAL*4 -- it shows eight significant decimal digits; 24 bits corresponds to 24 log₁₀2 = 7.22 decimal digits.

For more advanced applications, finer control over the conversion of the internal representation to decimal output is afforded by the format control edit descriptors RU, RD, RZ, RN, RC, and RP.

[Steven_L_Intel1]

Just a point of order - while I reposted the "Perils of Real Numbers" articles in this forum, they were written by compiler developer Dave Eklund for the original Digital Visual Fortran newsletter.

[SergeyKostrov]

FortranBeginner, Let me go back to your questions from the Initial post. >>...The debugger shows 1.7289999E-04. Yet I am expecting 1.729E-04 as a result... This is Not from the Debugger because Visual Studio has an Expression Evaluation software subsystem ( UI ) and it is used to display everything, like ( A + B - C ), in the Watch Window. >>...Does anyone know which FloatingPoint Settings (or other settings which I haven't thought about) I have to use to get the "correct" results? There are no such Unknown Floating-Point Unit settings. Once again, the Expression Evaluation software subsystem ( 1st ) and Floating-Point Unit settings ( 2nd and it is controled from a command line of the compiler or at runtime ) are Not related directly. The 1st one simply displays any values,depending on a data type, or expression(s), and it could have its own Rounding limitations. >>...Or are those "incorrect" result normal with Fortran? No. I would say the following: your program should produce correct results at runtime in Release Configuration and verify it first of all. Also, take a look at two screenshots.

[SergeyKostrov]

Screenshot #1: Note: Even if it is for VS 2005 you see that everything is correct.

[SergeyKostrov]

Screenshot #2: Note: At runtime it also shows a correct result.

[SergeyKostrov]

>>...Just a point of order - while I reposted the "Perils of Real Numbers" articles in this forum, they were written by compiler >>developer Dave Eklund for the original Digital Visual Fortran newsletter... Could you make a Sticky Post with links to these articles, and another as well, please? Thanks.

[Steven_L_Intel1]

See the "Useful Information and Frequently Asked Questions" link in the forum header.  Links are in that article. I didn't want to clutter the forum with lots of sticky posts.