# Developers Planet

Cnoobplusplus February 2016

### Eigen matrix/array compute distances efficiently

I need to create an Eigen array with all distances between each source and recorder position.

I have three Eigen::Array sx, sy and sz denoting the source positions and three Eigen::Array rx, ry and rz denoting the recorder positions. The length of source and recording position arrays are not necessary the same.

Using for-loops I can calculate the distance matrix as follow

``````Eigen::ArrayXf  sx, sy, sz;
Eigen::ArrayXf  rx, ry, rz;
Eigen::ArrayXXf dist;

for (int s=0; s<nrSources; s++ ) {

for (int r=0; r<nrReceivers; r++) {

dist(h,g) = sqrt(pow(rx(h)-sx(g),2.)+
pow(ry(h)-sy(g),2.)+
pow(sz(h)-sz(g),2.));
}
}
``````

I need to compute dist array like 500 x 1000 times for the experiment. Using the for-loop obvious works but it is probably not the most efficient method since it does not make use of the vectorization.

Rewriting the sx, sx, sz and hx, hy and hz into sxyz and hxyz arrays should be possible.

Is it possible to write the equation more efficiently!?

Avi Ginsburg February 2016

You probably want to discard the inner loop in favor of an expression. By doing so, we allow calculations to be lumped together and hopefully enable vectorization. Assuming the `r` and `s` variables were declared as `Eigen::ArrayX3f sxyz(nrSources,3), rxyz(nrReceivers,3);`, you could then write the outer loop as:

``````for (int s = 0; s < nrSources; s++)
{
dist.col(s) =
(rxyz.rowwise() - sxyz.row(s)).matrix().rowwise().norm();
}
``````

Here we use `rxyz.rowwise() - sxyz.row(s)` to subtract the s'th `s` from all `r`s. The `matrix()` is needed to give access to `norm()`.

Benchmark and compare to your implementation.