Intel® Fortran Compiler
Build applications that can scale for the future with optimized code designed for Intel® Xeon® and compatible processors.
28945 Discussions

Floating-Point Precision Issue: Exponentiation vs. Multiplication in Fortran

Haku
Novice
353 Views

During runtime, I calculate the value of a variable ust, which is of double precision type.

I then compute its second and third powers using two different methods:

  1. The exponentiation operator (ust**2, ust**3)
  2. Repeated multiplication (ust * ust, ust * ust * ust)

Next, I transfer the results of these calculations to their bit representations for comparison. As shown in the attached output, there is a slight difference in the bit patterns between the results obtained from the exponentiation operator(last digit is 1) and those from direct multiplication(last digit is 0).

The issue becomes more pronounced when I divide these results by another number. The output values differ slightly, which seems to propagate the discrepancy.

I’d like to understand:ust * ust * ust

  1. Why the bit patterns differ between the two methods (** vs. *).
  2. How this impacts subsequent calculations, such as division.
  3. Any recommended practices to minimize or resolve these differences.

 

real*8 :: ust
integer*8 :: ust_i, ust2_i, ust2_i2, ust3_i, ust3_i2
ust_i = transfer(ust,ust_i)
ust2_i = transfer(ust,ust2_i)
ust2_i2 = transfer(ust,ust2_i2)
ust3_i = transfer(ust,ust3_i)
ust3_i2 = transfer(ust,ust3_i2)

print *, "ust=", ust
print *, "ust_i=", ust_i
print *, "ust2_i=", ust2_i
print *, "ust2_i2=", ust2_i2
print *, "ust3_i=", ust3_i
print *, "ust3_i2=", ust3_i2
print '(5(A,E40.30,/))', "ust**3", ust**3, "ust*ust*ust", ust * ust * ust, "ust**3/karman", ust * ust * ust / karman

 

My software and hardware version is below:
IFX: 2024.2.1

Visual Studio 2022: 17.12.3

Windows: 11, 24H2

Processor: Intel(R) Xeon(R) Gold 6240 CPU

 

My project related property is below(Floating Point Section and Code Generation Section):

Haku_2-1737093534818.png

Haku_0-1737094718042.png

 

The output of the code:

Haku_1-1737095004052.png

 

 

Labels (1)
0 Kudos
3 Replies
Arjen_Markus
Honored Contributor I
304 Views

In general it helps us to help you if you post a complete program. In this particular case it does not really matter, at least as far as general considerations are likely to be the answer you are looking for.

The differences you see are caused by the finite precision and the rounding that takes place. DIfferent expressions that from a mathematical point of view are identical therefore can produce slightly different answers. It is inevitable, unfortunately. If you encounter a situation where these differences matter to the extent that you get practical differences, then the only real solution is to do the calculations in higher precision or to seek a different algorithm.

For instance:

1/3 + 1/3 + 1/3 = 1

But if you calculate it with a finite number of decimals, say 4 decimals, then:

0.3333 + 0.3333 + 0.3333 = 0.9999

This sort of differences occur with any arithmetic operation.

Haku
Novice
250 Views

Thank you for your reply, let me confirm that even I didn't use finite number of decimals explicitly, may the different expression x**3 and x * x * x generate slightly different results?

0 Kudos
JohnNichols
Valued Contributor III
211 Views

The following code 

!****************************************************************************

    program Console1

    implicit none

    ! Variables

    ! Body of Console1
    integer n 
    integer i
    
    real*8  m, p
    
    m = 1.0000
    n = 1000001
    p = 0.000001
    open(12,file="a.csv")
    do i = 1,n
    
    write(12,100)i,m,m*m,m**2,(m*m - m**2)
100 Format(i7, "  , ", F15.13, " ,  ", F15.13, "  , ", F12.10, "  , ", F12.10)
    m = m + p

end do

    end program Console1

A sample of the results are present in the XLSX file.  

If I take the remainder of m*m - m**2 multiply the results by a large 1000000000 (something like this)  I get a histogram of the results that looks like this 

 

Screenshot 2025-01-18 105240.png

Entirely not random or Gaussian and must be a function of the algorithm, I am sure the experts can explain same. 

 

 

0 Kudos
Reply