Hello, I originally posted this earlier today as part of an earlier question submitted to the community, but I think maybe it should be a “new” post. So I am cutting and pasting the most important parts of the message here, hope you don’t mind.
I have a 3d matrix as part of a code, initialized to zeros. This is for finite difference derivatives, and the faces of the 3d matrix are the guard channels. While recently debugging the code, I noticed that after
populating the interior of the matrix, sum(entire matrix) was slightly different than sum(the interior parts which are nonzero).
Is there some issue with Julia doing the sum of the full matrix (including zeros) vs just doing the sum over the interior points? I have written this example, running it a few times shows the differences are typically zero, but about 1/10 of the time the differences are (curiously) +/- 5.82076609e-11 on the Linux machine I am on (with Julia 1.8.2)
########################################################################
function matrRand2(nx,nG)
iLen = nx-nG*2
matrXYZ = zeros(Float64,iLen,iLen,iLen)
matrXYZ .= rand(iLen,iLen,iLen)
return matrXYZ
end
########################################################################
function printdiff(M0,MX,nx,nG)
iaRange = nG-1:nx-nG
str = zeros(1,2);
str1 = [sum(M0) sum(M0)]
str2 = [sum(MX) sum(M0[iaRange,iaRange,iaRange])]
str .= str1 .- str2
@printf(“The final answer = %15.8e %15.8e \n”,str[1],str[2])
return nothing
end
########################################################################
########################################################################
nx = 91;
nG = 2;
println(“VERS X”)
matr = zeros(Float64,nx,nx,nx)
matrX = zeros(Float64,nx-2nG,nx-2nG,nx-2*nG)
matrX = matrRand2(nx,nG)
matr[nG+1:nx-nG,nG+1:nx-nG,nG+1:nx-nG] = matrX;
printdiff(matr,matrX,nx,nG)
#########################
println(“VERS Y”)
matr = zeros(Float64,nx,nx,nx)
matrY = zeros(Float64,nx-2nG,nx-2nG,nx-2*nG)
matrY .= matrRand2(nx,nG)
matr[nG+1:nx-nG,nG+1:nx-nG,nG+1:nx-nG] = matrY;
printdiff(matr,matrY,nx,nG)