誘導されて来ました、
プログラムでこの部分がかなりの時間をとっているので改善したいと考えています。

#define M 300000
#define N 50000

int XY[M], b[N];
double x[N], y[N], z[N];

for(int i = 0; i < 40000; i++)
{
　for( int j = i+1; j < 40000; j++)
　{
　　r = sqrt((x[i]-x[j])*(x[i]-x[j]) +
　　　(y[i]-y[j])*(y[i]-y[j]) +
　　　(z[i]-z[j])*(z[i]-z[j]))+0.5;
　　XY[r] += b[i]*b[j];
　}
}

・x,y,zの有効数字は整数桁まで
・rの最大値は200000程度
・インラインアセンブラやMMXを使用せずに速度を上げたい

インデックスへのアクセスとsqrtの計算の部分がボトルネックのようです。
よろしくお願いします