Turn on suggestions

Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.

Showing results for

- Intel Community
- Software Development Tools (Compilers, Debuggers, Profilers & Analyzers)
- Intel® Fortran Compiler
- Round-off error while using double precision variables

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Mute
- Printer Friendly Page

Highlighted
##

Nagarada_Gadde__Srin

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-06-2018
01:52 PM

53 Views

Round-off error while using double precision variables

I have a parallelised MPI code which solves non-linear PDEs (basically LES simulations), written in Fortran 90. In a particular section of the code I have to multiply double precision variables with an other double precision variable.

u(iwrap_p(:),j,k) = alpha(:) * u(iwrap_p(:),j,k) + beta(:) * u(iwrap_p(:),j,k)

where alpha and beta are double precision variables whose values vary between 0. and 1.0, u is the velocity,and the sum of alpha and beta at any position is 1.0. The output of this operation should be equal to u(iwrap_p(:), j, k). But, performing this operation i.e. multiplication of the velocity with alpha and beta introduces round-off errors which keep on growing with time and finally the solution becomes non-physical.

Alpha and beta are calculated by the following,

beta = 0.5_rprec * ( 1.0_rprec - cos (pi * real (i - istart, rprec) & / real(iend - istart, rprec)) )

alpha = 1._rprec - beta

where 'rprec' is the user defined double precision datatype.

I am using the following compiler arguments

FFLAGS = -O3 -fp-model source

I have tried compiling with and without optimisations, but nothing seem to work.

The error disappears when alpha and beta are declared as single precision variables. What might be the cause of this error?

6 Replies

Highlighted
##

jimdempseyatthecove

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-07-2018
08:17 AM

53 Views

Try:

-fp-model precise -fimf-precision=high

Also assure that pi is correct for double precision.

Jim Dempsey

Highlighted
##

Sorry, something got lost in the narrative. If alpha(:) + beta(:) = 1, the calculation of u can be skipped completely, because it accomplishes nothing (other than to put in some round-off error).

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-07-2018
08:42 AM

53 Views

Highlighted
##

Nagarada_Gadde__Srin

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-07-2018
08:58 AM

53 Views

@mecej4 I am using the operation just as a sanity check. The actual code is the following:

u(iwrap_p(:),j,k) = alpha(:) * u(iwrap_p(:),j,k) + beta(:) * u_previous_domain(iwrap_p(:),j,k)

where u_previous_domain is the velocity from the precursor domain which is being copied to the next computational domain. It is basically the addition of a fraction of velocity from the precursor domain to the consumer domain which utilises the velocity.

@jimdempseyatthecove I will try it and let you know if it works out.

Highlighted
##

jimdempseyatthecove

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-07-2018
09:56 AM

53 Views

Your statement in #4 appears to be an Adams integration adjustment of an array (indexed by an index array) where u was previously advanced using Euler integration. If so, why didn't you write u once using Adams integration? Sketch:

velocity(:) = velocity(:) + (alpha(:) * acceleration(:) + beta(:) * prior_acceleration(:)) * delta_t

Typically alpha and beta do not vary.

Jim Dempsey

Highlighted
##

Nagarada_Gadde__Srin

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-07-2018
11:14 AM

53 Views

@jimdempseyatthecove Your suggestion in #2 did not work. The error still keeps growing.

Sorry about the lack of information in my question, perhaps I over-simplified it. I am not doing any integration. The simulation has two computational domains, one produces atmospheric boundary layer flow by solving Navier-Stokes equation (precursor domain), the second domain has a wind turbine (consumer domain); velocity from the first domain has to be given as input to the turbine domain. In the absence of the turbine, velocities in both the domain have to be identical, even after copying. The following operation does the copying:

u(iwrap_p(:),j,k) = alpha(:) * u(iwrap_p(:),j,k) + beta(:) * u_precursor_domain(iwrap_p(:),j,k)

Fractions alpha and beta correspond to the coefficients which slowly set the velocity in the consumer domain to the velocity in the precursor domain. Alpha varies from 1 to 0 and beta varies from 0 to 1 making the velocity at the exit of the consumer domain equal to the velocity at the exit of the precursor domain. For this test case as there is no turbine, u(iwrap_p(:), j, k) and u_precursor_domain(iwrap_p(:), j, k) are identical.

I might be expecting too much after using alpha and beta which obviously introduces round-off errors. But, I don't see the error when I use a single precision or quadruple precision for alpha, beta and the calculation of those coefficients. The error renders the velocity non-physical. I get perfect, physically meaningful solutions for both the precursor and consumer domains if I use single or quadruple precision for the variables.

Here are the typical values of alpha and beta:

Alpha Beta 0.933012701892219 6.698729810778065E-002 0.750000000000000 0.250000000000000 0.500000000000000 0.500000000000000 0.250000000000000 0.750000000000000 6.698729810778059E-002 0.933012701892219 0.000000000000000E+000 1.00000000000000 0.000000000000000E+000 1.00000000000000 0.000000000000000E+000 1.00000000000000 0.000000000000000E+000 1.00000000000000

Output from one of the operations where the coefficients are defined with double precision.

u : Before 22.3470033697708

u: After 22.3470033697707

Highlighted
##

I haven't run the following in Fortran so I cannot tell from your typical values of alpha and beta are correct. Incorporating those numbers into Excell shows line 06 above as producing a sum that is 1 lsb less than 1.0. If this tends to be the case (error rounds down), then I would check the code that generates and then force the condition that alpha(n) + beta(n) == 1.0. This may be facilitated by manipulating alpha(n) and/or beta(n) in a manner that statically is neutral in balancing the factors. (IOW tweak using the +/- EPSILON(alpha(n)) and/or +/- EPSILON(beta(n)) to produce 1.0 exactly.

jimdempseyatthecove

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

05-07-2018
11:33 AM

53 Views

For more complete information about compiler optimizations, see our Optimization Notice.