Rounding numbers (double- to single- precision)

[haminin]

Hi everyone

Function D2R is to round a double precision value to asingle one but sometimes does not work. For instance when A1 = .004D0; A2 = .031D0 the programs echoes 4.0000002E-03 3.0999999E-02 instead of 4.0000000E-03 3.1000000E-02.

Could someone tell me 1- whats wrong and 2- if in general there is a function which rounds anumber to a given decimal place or not.

Hamid

DOUBLE PRECISION A1,A2

REAL*4 D2R

A1 = .004D0; A2 = .031D0; PRINT*, D2R(A1), D2R(A2)

END PROGRAM

REAL*4 FUNCTION D2R(A)

DOUBLE PRECISION A

D2R = A

RETURN

END FUNCTION D2R

[TimP]

Nothing wrong. It looks like you are finding the nearest representable single precision binary numbers to the double precision values you set. Those will usually (not always) match to 7 decimal digits when range, not often to 8 decimal digits, except in the cases where the conversion is exact.

You can set a format to display only 6 decimal digits, which should give results more to your expectation.

[haminin]

Many thanks Tim

I need to compare their values not displaying them, so do you mean I can use internal writing to write only 6 decimal digits? Is it the best solution to round numbers?

Hamid

[TimP]

Your purpose isn't clear. Do you want to find out whether 2 numbers would round (usually, or always) to the same number in decimal format? You could use internal write to make character string representations, and compare them. Or, a more usual goal, to check whether they are within a specified relative tolerance?

If this is homework, the "point" may be to encourage you to understand something about binary floating point.

[rase]

You should never compare two floating point numbers directly, neither for single nor for double precision. The following code snippet can be used to do the trick:

real :: a, b
real, parameter :: eps=1.E-5 ! set to the precision needed
if ( abs(a-b) .lt. eps) then ! if absolute difference is less than eps, than the term becomes true

[haminin]

Many thanks again for your comment.

I wanted to round a number to a given decimal place. I wrote the function using internal writings. If J is not updated it clearly echoes the digits after the decimal point without rounding afterward digits.

The function set to change double- to single precision (N = 7).

REAL*8 FUNCTION D2R(A)

IMPLICIT NONE

DOUBLE PRECISION A, B, C

INTEGER E, J, N, I

CHARACTER*23 TEMP

N = 7 ! The number of places after the decimal point

TEMP = ''

WRITE(TEMP, FMT = '(E23.16E2)') A

! Number without exponent

READ(UNIT = TEMP, FMT = "(E.)") B

! The (N+1)-th digit after decimal point

READ(UNIT = TEMP, FMT = "(<3+N>X,I1)") J

! The number of digits in the exponent

READ(UNIT = TEMP, FMT = "(20X,I3)") E

!_____________________

IF (J .GE. 5) B = B + SIGN(10.D0 ** -N,B) ! Update last digit

IF (J .EQ. 4) THEN

LOOP_A : DO I = 1, 8 ! Check next digit

READ(UNIT = TEMP, FMT = "(<3+N+I>X,I1)") J

IF (J .LT. 4) EXIT

IF (J .GT. 4) THEN

B = B + SIGN(10.D0 ** -N,B) ! Update last digit

EXIT LOOP_A

END IF

! If J = 4 --> continue checking next numbers

END DO LOOP_A

END IF

!_____________________

D2R = B * 10.D0 ** E

RETURN

END FUNCTION D2R

[haminin]

Many thanks Rase. I couldnt understand you well, but my problem is now solved by function D2R.

[Steven_L_Intel1]

I'd suggest instead something like this:

[plain]REAL(4) FUNCTION D2R (A,PLACES)
DOUBLE PRECISION A
INTEGER PLACES

D2R = ANINT(A*(10.D0**PLACES))/(10.D0**PLACES)
RETURN
END[/plain]

This relies on the fact that floating point values can exactly represent integers (within a range based on the precision)

Do keep in mind that single precision is good to only about six significant digits so that it can't exactly represent all possible rounded values of more digits.

[haminin]

Its interesting. I didnt know about ANINT function. It should be useful for this purpose, but I think your D2R function needs to be revised, for instance when A = -.031123354445444443D0 * 10.D0 ** 10 it always gives -3.1123354E+08. Anyway thank you very much for this.

Hamid

[Steven_L_Intel1]

That's why there's a divide after the ANINT. The ANINT rounds and the divide reestablishes the magnitude. But again you're limited by the number of bits available in single precision. Don't expect reliable results when rounding past 6 significant digits (including any to the left of the decimal.