A largest subset of a cubic lattice with unique distances between its points Let $n$ be a positive integer. Consider a set $S_n=([1,n]\cap\mathbb N)^3$ having $n^3$ points in the Euclidean space $\mathbb R^3$ arranged in a cubic lattice. Let $a_n$ be the size of a largest subset of $S_n$ such that for each pair of points in the subset there are no other pair of points in the subset having the same distance between them (i.e. all distances between points in the subset are unique).

Can you suggest an efficient way to compute $a_n$, or find a general formula for this sequence? Or, at least, find $a_n$ for a few small positive integers?

Update: Lower or upper bounds, or asymptotic behavior are welcome. Is the sequence $\{a_n\}$ strictly increasing?

Update: The corresponding OEIS sequence: $A275672$ (only terms $n=0..6$ are currently known).
 A: Promoting the comments to an answer, as the upper bound seems to be surprisingly good early on.
We want to avoid repetitions of squared distances. The squared distance $d$ between two points in $S_n$ is of the form
$$
d=a^2+b^2+c^2,
$$
where $a,b,c$ are integers in the range $[0,n-1]$. Let us denote by $D_n$
the set of such sums $>0$. The exact size of the set $D_n$ is difficult to write down, but it is obviously $<3(n-1)^2$, and does grow roughly as a quadratic polynomial function of $n$.
Call a subset $S\subseteq S_n$ a 3-dimensional Golomb ruler, 3DGR for short, if there are no repetitions among the squared distances between its points. If $S$ is a 3DGR with $m$ elements, then the squared distances between its points form a subset of $\binom m2$ integers of the set $D_n$. Hence we have the inequality
$$
\binom m2\le |D_n|.\qquad(*)
$$
So if $a_n$ is the maximum size of a 3DGR $S\subseteq S_n$, then
we have an upper bound $a_n\le m$, where $m$ is the largest integer satisfying $(*)$.

Example. Consider the case $n=3$. We see easily that the set $D_3$ consists of integers
$$
D_3=\{1,2,3,4,5,6,8,9,12\},
$$
a total of nine elements. Because $\binom 42<9<\binom 52$ we can deduce that a 3DGR $\subseteq S_3$ can have at most four elements. It is easy to check that
$$
S=\{(1,1,1),(1,1,3),(2,3,3),(3,3,3)\}
$$
is a 3DGR (yielding the squared distances $1,4,5,8,9,12$). Therefore the bound is tight in this case, and we can conclude that $a_3=4$.

Here's a short table (thanks, Mathematica) of the upper bounds to $a_n$:
$$
\begin{array}{c|c|ccc|c|c}
n&|D_n|&a_n\le&&n&|D_n|&a_n\le\\
\hline
1&0&1&&11&178&19\\
2&3&3&&12&211&21\\
3&9&4&&13&259&23\\
4&18&6&&14&299&24\\
5&31&8&&15&346&26\\
6&44&9&&16&401&28\\
7&66&12&&17&463&30\\
8&87&13&&18&516&32\\
9&115&15&&19&591&34\\
10&144&17&&20&648&36
\end{array}
$$
So we see that the upper bound grows linearly as a function of $n$ (follows immediately from the estimates on the size of $D_n$). I dare not guess, how the true values of $a_n$ behave. The case of 1-dimensional Golomb rulers is difficult enough, and huge amounts of computing resources have been put towards improving their construction.
A: For each $r \in [1, \ldots, 3 n^2]$ let $E_r$ be the set of unordered pairs $(p,q)$ of points in $S_n$ such that $\|p-q\|^2 = r$.  We want a largest possible subset $X$ of $S_n$ such that $X$ contains at most one of each $E_r$.
It can be considered as a quadratically constrained integer programming problem:
maximize $\sum_{p \in S_n} x_p$ 
subject to $\sum_{(p,q) \in S_r} x_p x_q \le 1$ for each $r$
all $x_p \in \{0,1\}$
Not that there are efficient methods to solve this, but CPLEX ought to be able to do this for small $n$.
EDIT:  It appears that a careful integer linear programming model (with binary variables for both edges and vertices) works better. The results for small $n$ are 
$$ \matrix{n & a_n\cr
           1 & 1\cr
           2 & 3\cr
           3 & 4\cr
           4 & 6\cr}$$
I have some faint hopes of getting $a_5$ (so far CPLEX has found a solution with $7$ vertices and an upper bound of $19.51$)
EDIT:  CPLEX crashed before obtaining an optimal solution for $n=5$.
An example to show $a_5 \ge 7$ is 
$$\pmatrix{     3  &   4  &   2\cr
     1  &   2  &   2\cr
     4  &   4  &   2\cr
     5 &   5   &  2\cr
     5  &   1  &   4\cr
     5  &   5  &   4\cr
     2  &   1  &   5\cr}$$
EDIT:
For $n=6$, Jyrki's estimate shows $a_6 \le 9$, and (using tabu search) I was able to find an example to show $a_6 = 9$:
$$ \pmatrix{   4  &   3 &    1\cr
     1  &   5  &   4\cr
     1  &   2  &   6\cr
     6  &   1  &   1\cr
     1  &   6  &   4\cr
     2  &   5  &   3\cr
     2  &   2  &   2\cr
     6  &   2  &   6\cr
     1  &   6  &   6\cr}$$
EDIT: I was able to find an example to show $a_7 \ge 10$:
$$\pmatrix{
    1  &   7  &   7\cr
     5 &    7  &   1\cr
     7 &    2  &   2\cr
     7 &    5  &   5\cr
     1 &    2  &   2\cr
     1 &    4  &   7\cr
     2 &    7  &   7\cr
     7 &    6  &   4\cr
     1 &    7  &   5\cr
     6 &    3  &   1\cr}$$
Still trying for $11$.
