Is it possible to compute only the adjacencies of Voronoi cells in the Voronoi diagram efficiently?

Yes, it can be done very efficiently by linear programming (LP), and very importantly this can be done for very large scale problems, with practically no bounds on the size with an efficient LP solver.

The method is simple. The lifting technique we described in 3.2 immediately gives the idea. Recall that the Voronoi diagram of a set $S$ of $n$ points in $R^d$ is the projection of the following $(d+1)$-polyhedron to $R^d$ space of the first $d$ components.

$$P = \{ x\in R^{d+1} \; \vert \; \sum_{j=1}^d p_j^2 - \sum_{j=1}^d 2 p_j x_j + x_{d+1} \ge 0 \quad \forall p \in S \}.$$

For simplicity, denote it as

$$P = \{ x\in R^{d+1} \; \vert \; b - A x \ge 0 \},$$

where $A$ is a given $n \times (d+1)$ matrix and $b$ is a $n$-vector. Now for each $i=1,\ldots, n$, consider the $i$th facet $F_i$ of $P$:

$$F_i = \{ x\in R^{d+1} \; \vert \; b - A\; x \ge 0 b_i - A_i \; x \le 0 \},$$ (7)

Two facets $F_i$ and $F_j$ are called adjacent if the intersection $F_i \cap F_j$ is a facet of both, i.e. has dimension $d-2$. An equivalent definition is: they are adjacent if (*) the facet $F_i$ becomes larger once the facet $F_j$ is removed from the polyhedron, i.e. the $j$th inequality is removed from the system $b - A \; x \ge 0$.

It is easy to see that two Voronoi cells $vo(p^i)$ and $vo(p^j)$ are adjacent if and only if the corresponding facets $F_i$ and $F_j$ are adjacent in the polyhedron $P$. Now, we formulate the following LP for any distinct $i$, $j = 1,2, \ldots,n$:

$$\begin{align*}\begin{array}{lllllll} &\text{minimize} & \quad f(x) := & b_j & - ... ... & x \; \ge \; 0 \\ & & & b_i & - & A_i & x \; \le \; 0, \end{array}\end{align*}$$ (8)

where $b'$ is equal to $b$ except for $j$th component $b'_j = b_j + 1$. The new inequality system $b' - A \; x \ge 0$ is simply a modification of the original system obtained by relaxing the $j$th inequality a little bit. An important remark is, by definition (*), $F_j$ and $F_i$ are adjacent if and only if the objective value $f(x)$ is negative at an optimum solution. Thus we formulated the Voronoi adjacency computation as an LP problem.

How much do we gain by using LP for the adjacency computation, instead of computing the whole Voronoi diagram? A lot. It is hard to exaggerate this, because the LP (8) (in fact any LP) is solvable in polynomial time, whereas the associated Voronoi computation is exponential in $d$ and $n$. Using the standard simplex method, the time complexity of solving an LP is not polynomial, but the practical complexity is roughly $O(n d^3)$.