Community
cancel
Showing results for
Did you mean:
Highlighted
Beginner
22 Views

DO WHILE and LESS THAN in REAL numbers

I have two codes to reproduce the situation:

Code1

PROGRAM test
REAL(KIND=8) dt, t, tmax
INTEGER i, imax
PARAMETER (tmax = 10.2D0, dt = 0.2D0, imax = 10)
t = 0
i = 0
DO WHILE (t .LT. tmax)
PRINT *, "before", t, i
t = t + dt
i = i + 1
PRINT *, "after", t, i
END DO

END PROGRAM test

Code 2

PROGRAM test
REAL(KIND=8) dt, t, tmax
INTEGER i, imax
PARAMETER (tmax = 10.2D0, dt = 0.2D0, imax = 10)
t = 0
i = 0
DO WHILE (i .LT. imax)
PRINT *, "before", t, i
t = t + dt
i = i + 1
PRINT *, "after", t, i
END DO

END PROGRAM test

If I look at the output, the last 3 terms code 1 prints:

after   10.2000000000000               51
before   10.2000000000000               51
after   10.4000000000000               52

and the code 2

after   1.80000000000000                9
before   1.80000000000000                9
after   2.00000000000000               10

The code 2 is doing what it is supposed to, code 1 is doing one more loop it is not supposed to do. I believe the reason is the machine representation of numbers. But what do I do to get around it, to get a real "less than" condition?

Tags (1)
4 Replies
Highlighted
Black Belt
22 Views

You could determine the loop count before you start and use an integer counter.

E.G. nloops = NINT((tmax-tmin)/dt)

assuming tmin is set.

Testing reals for equality is always going to lead to grief unless you allow for the possibility of round off error.

Highlighted
Beginner
22 Views

Thank you for that. But then, there is an interesting (and for me embarassing). How do you determine the loop cound (I think  you should have INT( (tmax-tmin)/dt ) +1 in your post, but anyway, this is valid only for the loop of the type where I wanted DO WHILE (t .LT. tmax)  but this no longer works for the loop type DO WHILE (t .LE. tmax) and there is a difference from the above number of loops by 1 and this depends on dt by I have not been able to figure out how. It is always one less or the same as INT( (tmax-tmin)/dt ) +1

To clarify, I want the same number of loops as would have happened in the case t .LE. tmax in the case we were able to do precise arithmentics (no round off errors etc).

Highlighted
Valued Contributor II
22 Views

The phenomenon you are seeing is precisely the reason why real numbers are tricky to use as a loop variable. One improvement:

Use t = tmin + dt * (i-1)

instead of accumulating the variable t. Accumulation means that also the rounding error is accumulated.

If you seek a lower-than loop:

n = int( (tmax - tmin - 0.5d0*dt) / dt)

ought to do it

For a lower-equal:

n = int( (tmax - tmin + 0.5d0*dt) / dt)

Caveat: I have not tested these formulae in any formal sense, so I may be wrong :).

Highlighted
Black Belt
22 Views

To avoid the accumulation of roundoff errors...
you avoid accumulation altogether.

PROGRAM test
REAL(KIND=8) dt, t, tmax
INTEGER i, imax
PARAMETER (tmax = 10.2D0, dt = 0.2D0, imax = 10)
t = 0
i = 0
DO WHILE (i*dt .LT. tmax)
PRINT *, "before", t, i
t = t + dt
i = i + 1
PRINT *, "after", t, i
END DO

END PROGRAM test

Jim Dempsey