Using 1D Global Array Index and rounding errors

Accelerated Computing CUDA CUDA Programming and Performance
1

I am using a 1D global index constructed from the blockidx and threadidx:

indx = blockDim%x*(blockIdx%x-1) + threadIdx%x

I then deconstruct this 1D array (e.g. 1,2,3…1000000) into a 7D index, i.e. array(1,3,2,4,1,2,3), through some maths, e.g.:

Iterator1 = MOD(indx-1, 20)+1
Iterator2 = MOD((indx-1)/20, 20)+1
Iterator3 = MOD((indx-1)/(20*20), 10)+1

and so on allowing me to reference the array, array(iterator1,iterator2,…)

I then print out what it is supposed to be (i.e. me referencing the 7D array element) and then deconstructing the value back:

WRITE(5,'(10I3, F20.6)') Iterator1, Iterator2,..., Array(Iterator1,Iterator2,...)

I get the following output:

First line is input, second line is output

Various deconstructed values      Overall 1D index
12 -1 -2  1 18 17 11  1  1        33444738.000000
12 -1 -2  1 18 17 11  1  1        33444738

12 -1 -2  1 19 17 11  1  1        33444740.000000
12 -1 -2  1 20 17 11  1  1        33444740

12 -1 -2  1 20 17 11  1  1        33444740.000000
12 -1 -2  1 20 17 11  1  1        33444740

12 -1 -2  1  1 18 11  1  1        33444740.000000
12 -1 -2  1 20 17 11  1  1        33444740

12 -1 -2  1  2 18 11  1  1        33444742.000000
12 -1 -2  1  2 18 11  1  1        33444742

As can be seen I am getting rounding errors, the first set is 33444738, it should be followed by 33444739, but instead is followed by 33444740, as a result I am incorrectly calculating the values.

This doesnt occur at low values for the 1D index, only at larger values (>1m). What could be causing this?
The larger the index the more errors occur (at around ~150m there are 15 values merged).

2

one possibility: overflow of exact representation of a 32-bit floating point quantity (it has about 23 mantissa bits, allowing for exact representation of integers up to about 16 million). That’s really just a guess, however.

3

Hi, thanks for the reply, is there away for me to exceed this limit?

I tried defining indx as a interger*8, but the inaccuracy still occurs.

Does the imprecision (error) occur when defining indx?

4

it might help if you showed a more complete example. For example, not showing how indx and IteratorN are defined leaves me to guess.

5

I define them as:

attributes(global) subroutine geArray(AnArray)
implicit none
integer*8 :: indx, maxIDX
integer :: Iterator1, Iterator2
integer, parameter :: CAPItr1 = 20, CAPItr2 = 20

REAL, intent(inout) :: AnArray(1:CAPItr, 1:CAPItr1, ...)

maxIDX = Iterator1*Iterator2*...

indx = blockDim%x*(blockIdx%x-1) + threadIdx%x

IF(indx<=maxIdx) THEN
   Iterator1 = MOD(indx-1,CAPItr1)+1
   Iterator2 = MOD((indx-1)/CAPItr1, CAPItr2)+1

   AnArray(Iterator1,Iterator2,...) = indx
ENDIF

They are defined the same way on both the device code and host code.

6

The problem was that my array, AnArray was defined as a real4 not real8, thats why the rounding was occuring