If the graph $G=(V,E)$ has a lot of symmetries, then using **spectral realizations** might give you nice drawings that reflects these symmetries. The success of this method (e.g. whether the drawing is planar, whether all vertices are on a sphere) depends on a lot of factors, some of which are not completely clear to me. However, what I can tell you is that it works for the graphs of all uniform 3-polytopes (so, e.g., the dodecahedron).

I explain the most straight forward way to do it, some tweaks might be neccessary for the general case:

**Costruction.** Let $\theta$ be an eigenvalue of (the adjacency matrix of) $G$, and $v_1,v_2,v_3\in\Bbb R^n$ three ortho-normal eigenvectors to $\theta$. Construct the matrix $M:=(v_1,v_2,v_3)\in\Bbb R^{n\times 3}$ with the $v_i$ as columns. The rows of that matrix are a 3-dimensional embedding of the vertices of $G$.

Usually, you should take $\theta_2$, i.e., the second-largest eigenvalue of the adjacency matrix of $G$. Surprisingly, this eigenvalue has multiplicity three for most symmetric graphs that come from 3-polytopes (exceptions are, as far as I know, only prisms). This means, you cannot do anything wrong by choosing just any orthonormal basis of eigenvectors.

Here is code for Mathematica to automatically find a nice drawing of the dodecahedral graph:

```
G = GraphData["DodecahedralGraph"];
A = AdjacencyMatrix[G];
n = VertexCount[G];
eval = Eigenvalues[A // N];
th2 = RankedMax[eval, 2];
evec = NullSpace[A - th2*IdentityMatrix[n]];
GraphPlot3D[G,
VertexCoordinateRules -> Table[i -> evec[[{1,2,3}, i]], {i, 1, n}]
]
```

Output:

If $\theta_2$ has multiplicity $<3$, you can add eigenvectors of other eigenvalues until you have three, preferably from the next largest eigenvalues. Just do not use the largest eigenvalue.

For your example in the comments (the capped cube), we have the problem that $\theta_2$ has multiplicity one. However, as explained, we could use eigenvectors of $\theta_3$ (which has multiplicity two) to complete to a set of three vectors. Use this:

```
th2 = RankedMax[eval, 2];
th3 = RankedMax[eval, 3];
evec = NullSpace[A - th2*IdentityMatrix[n]] ~Join~ NullSpace[A - th3*IdentityMatrix[n]];
```

Output:

For general graphs one should certainly follow a more dynamic approach, e.g. sorting the eigenvectors by their eigenvalues and taking the largest three (but not the largest one). Something like this:

```
vec = SortBy[Transpose[Eigensystem[A // N]], First][[{-2, -3, -4}, 2]];
GraphPlot3D[G,
VertexCoordinateRules -> Table[i -> vec[[;; , i]], {i, 1, n}]
]
```

17th Internat Symp Experimental Algorithms. 2018. PDF download $\endgroup$Proceedings of the International Conference on Advanced Visual Interfaces, pp. 1-9. 2020. $\endgroup$