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

custom INT to handle rounding errors of operations

Leos_P_
Beginner
324 Views

Problem of INT() is that it is not necessarily predictable in the sense e.g. INT(5.D0/1.D0) is not necessarily 5 but cal also happen to be 4. NINT() does not really help as NINT(4.8) = 5. In other words, I would need exact arithmetics. So I tried to write a simple function:

FUNCTION f_int(x) RESULT(i_x)
REAL(KIND=8) x, DEPS
INTEGER i_x
PARAMETER ( DEPS = 5.D-15 )
IF ( ABS( NINT(x) - x ) .GT. DEPS) THEN
	i_x = INT(x)
ELSE
	i_x = INT(NINT(x))
ENF IF
END FUNCTION f_int

My question is, does this make sense for the purpose I need and what can I be missing (while I realise that if I wanted to do things like e.g. 5.99999999978/1.999999999999999 on purpose, it would probably fail but I can live with this)?

0 Kudos
3 Replies
jimdempseyatthecove
Honored Contributor III
324 Views

>>INT(5.D0/1.D0) is not necessarily 5 but cal also happen to be 4.

Incorrect, the above will produce 5.

However note, INT(SomeEpressionWithErrorThatMathematicallyShouldResultIn5.0D0 / 1.0D0) may not necessarily be 5.

Sketch of questionable code

ValueWithError = 0.1D0 ! infinitely repeating fraction (IOW contains error in LSB)
...
Result = DoubleFunctionThatProduces~5.0D0using(ValueWithError)
I_X = INT(Result / 1.0D0) ! not necessarily 5.

Consider  using:

nBitsPrecision = ... ! you specify here

if(Result .NE. 0.0D0) Result = Result + (NEAREST(Result,Result) * (2.0D0**nBitsPrecision))
I_X = INT(Result)

You can also use EPSILON to similar effect.

Jim Dempsey

0 Kudos
mecej4
Honored Contributor III
324 Views

My question is, does this make sense for the purpose...

You have not yet declared what your purpose is.

In other words, I would need exact arithmetics.

Not available with machine floating point.

 

 

0 Kudos
FortranFan
Honored Contributor II
324 Views

@Leos P.,

It would appear you're looking for the ANINT standard intrinsic function in Fortran: https://software.intel.com/en-us/node/679004

If not, you would need to furnish more details.

 

0 Kudos
Reply