What sort of indexing method can I use to store the distances between X^2 vectors in an array without redundancy?

Question

I'm working on a demo that requires a lot of vector math, and in profiling, I've found that it spends the most time finding the distances between given vectors.

Right now, it loops through an array of X^2 vectors, and finds the distance between each one, meaning it runs the distance function X^4 times, even though (I think) there are only (X^2)/2 unique distances.

It works something like this: (pseudo c)

#define MATRIX_WIDTH 8

typedef float vec2_t[2];
vec2_t matrix[MATRIX_WIDTH * MATRIX_WIDTH];

...

for(int i = 0; i < MATRIX_WIDTH; i++)
{
    for(int j = 0; j < MATRIX_WIDTH; j++)
    {
        float xd, yd;
        float distance;

        for(int k = 0; k < MATRIX_WIDTH; k++)
        {
            for(int l = 0; l < MATRIX_WIDTH; l++)
            {
                int index_a = (i * MATRIX_LENGTH) + j;
                int index_b = (k * MATRIX_LENGTH) + l;

                xd = matrix[index_a][0] - matrix[index_b][0];
                yd = matrix[index_a][1] - matrix[index_b][1];

                distance = sqrtf(powf(xd, 2) + powf(yd, 2));
            }
        }

        // More code that uses the distances between each vector
    }
}

What I'd like to do is create and populate an array of (X^2) / 2 distances without redundancy, then reference that array when I finally need it. However, I'm drawing a blank on how to index this array in a way that would work. A hash table would do it, but I think it's much too complicated and slow for a problem that seems like it could be solved by a clever indexing method.

EDIT: This is for a flocking simulation.

Answer 1

performance ideas: a) if possible work with the squared distance, to avoid root calculation b) never use pow for constant, integer powers - instead use xd*xd

I would consider changing your algorithm - O(n^4) is really bad. When dealing with interactions in physics (also O(n^4) for distances in 2d field) one would implement b-trees etc and neglect particle interactions with a low impact. But it will depend on what "more code that uses the distance..." really does.

just did some considerations: the number of unique distances is 0.5*n*n(+1) with n = w*h. If you write down when unique distances occur, you will see that both inner loops can be reduced, by starting at i and j.

Additionally if you only need to access those distances via the matrix index, you can set up a 4D-distance matrix.

If memory is limited we can save up nearly 50%, as mentioned above, with a lookup function that will access a triangluar matrix, as Code-Guru said. We would probably precalculate the line index to avoid summing up on access

float distanceArray[(H*W+1)*H*W/2];
int lineIndices[H];

searchDistance(int i, int j)
{
    return i<j?distanceArray[i+lineIndices[j]]:distanceArray[j+lineIndices[i]];
}