Hey,
in several places in the code, indices are computed like
short aStartingGridPoint = gridDimensionA_ * aCoord;
with some precautions to ensure 0 <= aCoord < 1 so that aStartingGridPoint < gridDimensionA_ should always evaluate to true.
If aStartingGridPoint >= gridDimensionA_, the code will run into an out-of-bounds access of an array.
Unfortunately, despite the precautions taken, this happens on extremely rare occasions. I do not completely understand why.
My best guess is some kind of numerical error originating in the way short * Real is rounded on some platforms. However, I couldn't reproduce this problem in a minimal example. My code (hundreds of thousands of charges and coordinates) consistently runs into this issue though.
I could confirm this is happening by inserting if (aStartingGridPoint == gridDimensionA_)-checks into the code.
Something like this fixes the problem, though a little dirtily:
inline short computeClampedIndex(unsigned int dim, Real coord) {
short startingGridPoint = dim * coord;
if (startingGridPoint == dim) {
startingGridPoint -= 1;
}
return startingGridPoint;
}
...
short aStartingGridPoint = computeClampedIndex(gridDimensionA_, aCoord);
Any ideas why this is happening, despite the precautions which are already taken in the code to prevent this very problem?
Best,
Eik
Hey,
in several places in the code, indices are computed like
short aStartingGridPoint = gridDimensionA_ * aCoord;with some precautions to ensure
0 <= aCoord < 1so thataStartingGridPoint < gridDimensionA_should always evaluate to true.If
aStartingGridPoint >= gridDimensionA_, the code will run into an out-of-bounds access of an array.Unfortunately, despite the precautions taken, this happens on extremely rare occasions. I do not completely understand why.
My best guess is some kind of numerical error originating in the way
short * Realis rounded on some platforms. However, I couldn't reproduce this problem in a minimal example. My code (hundreds of thousands of charges and coordinates) consistently runs into this issue though.I could confirm this is happening by inserting
if (aStartingGridPoint == gridDimensionA_)-checks into the code.Something like this fixes the problem, though a little dirtily:
Any ideas why this is happening, despite the precautions which are already taken in the code to prevent this very problem?
Best,
Eik