CODEWITHSUNDEEP
c++precision
Jian Gao
Jian Gao

Reputation: 11

calculating the mean of an array c++

I encountered a problem when I tried to calculate the mean of an array in two ways. Below is the code: 

float sum1, sum2, tmp, mean1, mean2;
double sum1_double, sum2_double, tmp_double;
int i, j;
int Nt=29040000;  //array size
int piecesize=32;
int Npiece=Nt/piecesize;
float* img;
float* d_img;
double* img_double;
img_double = (double*)calloc(Nt, sizeof(double));
cudaHostAlloc((void**)&img, sizeof(float)*Nt, cudaHostAllocDefault);
cudaMalloc((void**)&d_img, sizeof(float)*Nt);
...
//Some calculation is done in GPU and the results are stored in d_img;
...    
cudaMemcpy(img, d_img, Nt*sizeof(float), cudaMemcpyDeviceToHost);
for (i=0;i<Nt;i++) img_double[i]=(double)img[i];

//Method 1
sum1=0;
for (i=0;i<Nt;i++) 
{ sum1 += img[i]; }

sum1_double=0;
for (i=0;i<Nt;i++) 
{ sum1_double += img_double[i]; }

//Method 2
sum2=0;
for (i=0;i<Npiece;i++)
{   tmp=0; 
      for (j=0;j<piecesize;j++)
        { tmp += img[i*piecesize+j];}
    sum2 += tmp;
}

sum2_double=0;
for (i=0;i<Npiece;i++)
{   tmp_double=0; 
      for (j=0;j<piecesize;j++)
        { tmp_double += img_double[i*piecesize+j];}
    sum2_double += tmp_double;
}

mean1=sum1/(float)Nt;
mean2=sum2/(float)Nt;
mean1_double=sum1_double/(double)Nt;
mean2_double=sum2_double/(double)Nt;

cout<<setprecision(15)<<mean1<<endl;
cout<<setprecision(15)<<mean2<<endl;
cout<<setprecision(15)<<mean1_double<<endl;
cout<<setprecision(15)<<mean2_double<<endl;

Output:

132.221862792969
129.565872192383
129.565938340543
129.565938340543

The results obtained from the two methods, mean1=129.6, mean2=132.2, are significantly different. May I know why?

Thanks a lot in advance!

Upvotes: 0

Views: 124

Answers (1)

geza
geza

Reputation: 29962

The reason is that floating point arithmetic is not precise. When you accumulate integers, float becomes imprecise when abs(value) becomes larger than 224 (I'm supposing IEEE-754 32-bit here). For example, float is incapable to store 16777217 precisely (it will become 16777216 or 16777218, depending on the rounding mode).

Supposedly your second calculation is the more precise one, as less precision is lost, because of the separate tmp accumulation.

Change your sum1, sum2, tmp variables to long long int, and hopefully you'll get the same result for both calculations.

Note: I've supposed that your img stores integer data. If it stores floats, then there is no easy way to fix this perfectly. One way is to use double instead of float for sum1, sum2 and tmp. The difference will be there, but it will be much smaller. And there are techniques how to accumuluate floats more precisely than simple summing. Like Kahan Summation.

Upvotes: 4

Related Questions