I will try to give a concise summary of what I wrote below. I understand that it is very long and apologize if I am wasting your time.
I used the Liénard-Wiechert potential and the Lorentz force formula to derive equations of motion for charged particles that do not involve the $\mathbf{E}$ and $\mathbf{B}$ fields. I use this in a two-body problem setup with two charges, and want to see whether the resulting equations result in spiraling orbits.
More specifically, I plug in the formulas for $\mathbf{E}_2$ and $\mathbf{B}_2$ given by Liénard-Wiechert into the Lorentz force formula $\mathbf{F}_{12}(\mathbf{r}_1(t), t) = q_1(\mathbf{E}_2(\mathbf{r}_1(t), t) + \mathbf{\dot{r}}_1(t)\times\mathbf{B}_2(\mathbf{r}_1(t), t))$ and use the classical equation of motion $m_1\mathbf{\ddot{r}}_1(t) = \mathbf{F}_{12}(\mathbf{r}_1(t), t)$. A second equation is obtained for particle 2 using the expressions for $\mathbf{E}_1$ and $\mathbf{B}_1$ instead.
The full equations are complicated so I take a limit as the mass of the nucleus goes to infinity, and (assuming my calculations are correct), I obtain a simpler equation of motion formally identical to that of planetary motion around the sun. This predicts elliptical trajectories for most initial conditions.
The concern is that the resulting motion is stable, and this seems to contradict what the usual accounts say about classical electromagnetism failing to explain atomic orbit stability.
If you want more detail, you can also read what I wrote below.
Suppose for a moment that we work classically and use the planetary model for the atom (hydrogen, for simplicity), as a positively charged nucleus with a negatively charged electron orbiting it in a similar way that planets orbit around the Sun.
The usual story goes that due to the electron orbiting around the nucleus, it undergoes acceleration which will cause the electron to radiate electromagnetic waves. The energy of this radiation is considered to be taken away over time from the electron's total energy via the Larmor formula, and therefore the model predicts a collapse of the electron's orbit over a short period of time because the radius of the orbit must decrease to compensate for the diminished energy of the electron.
At the risk of sounding ridiculous to more knowledgeable people, I would like to challenge this assumption with the following considerations (only as a way to clarify my own understanding of the problem). It appears to me that this problem is arising only because we consider the electromagnetic field to have an existence independent of the charges generating it, and carrying energy and momentum on its own.
But it seems possible to describe the fundamentals of (classical) electromagnetism without resorting to the concept of electromagnetic fields, by the use of a combination of the Lorentz force law and the Liénard-Wiechert potential. In particular, one can substitute the explicit expressions for the $\mathbf{E}$ and $\mathbf{B}$ fields obtained from the Liénard-Wiechert formulas into the Lorentz force formula to derive the force between two charged particles moving on arbitrary paths in space. One can then derive a classical equation of motion for the particles using Newtonian mechanics, or its special relativistic correction.
Explicitly, we obtain this system of two ODEs, where $m_1, m_2, q_1, q_2$ are the masses and charges of the two particles, $\mathbf{r}_1(t), \mathbf{r}_2(t)$ are the paths, and other quantities are defined in the Wikipedia page for the Liénard-Wiechert potential:
$$m_1{\mathbf{\ddot{r}}_1}(t) = {\frac {\mu_0c^2q_1q_2}{4\pi}}\left(1 + \left[\mathbf{\dot{r}}_1(t)\times\left[\frac{\mathbf{n}_2(t_r)}{c}\times\right]\right]\right)$$
$${\left[\frac{\mathbf {n}_2(t_r) -{\boldsymbol {\beta}_2(t_r)}}{\gamma_2 ^{2}(t_r)(1-\mathbf {n}_2(t_r) \cdot {\boldsymbol {\beta }_2(t_r)})^{3}|\mathbf {r}_1(t) -\mathbf {r}_2(t_r)|^{2}}+\frac{\mathbf {n}_2(t_r) \times {\big (}(\mathbf {n}_2(t_r) -{\boldsymbol {\beta }_2(t_r)})\times {{\boldsymbol {\dot{\beta} }_2(t_r)}}{\big )}}{c(1-\mathbf {n}_2(t_r) \cdot {\boldsymbol {\beta }_2(t_r)})^{3}|\mathbf {r}_1(t) -\mathbf {r}_2(t_r)|}\right]}$$
$$m_2{\mathbf{\ddot{r}}_2}(t) = {\frac {\mu_0c^2q_1q_2}{4\pi}}\left(1 + \left[\mathbf{\dot{r}}_2(t)\times\left[\frac{\mathbf{n}_1(t_r)}{c}\times\right]\right]\right)$$
$${\left[\frac{\mathbf {n}_1(t_r) -{\boldsymbol {\beta}_1(t_r)}}{\gamma_1 ^{2}(t_r)(1-\mathbf {n}_1(t_r) \cdot {\boldsymbol {\beta }_1(t_r)})^{3}|\mathbf {r}_2(t) -\mathbf {r}_1(t_r)|^{2}}+\frac{\mathbf {n}_1(t_r) \times {\big (}(\mathbf {n}_1(t_r) -{\boldsymbol {\beta }_1(t_r)})\times {{\boldsymbol {\dot{\beta} }_1(t_r)}}{\big )}}{c(1-\mathbf {n}_1(t_r) \cdot {\boldsymbol {\beta }_1(t_r)})^{3}|\mathbf {r}_2(t) -\mathbf {r}_1(t_r)|}\right]}$$
where $\mathbf {n}_2(t_r) = \frac{\mathbf{r}_1(t)-\mathbf{r}_2(t_r)}{|\mathbf{r}_1(t)-\mathbf{r}_2(t_r)|}$, $\boldsymbol {\beta}_2(t_r) = \frac{\mathbf{\dot{r}}_2(t_r)}{c}$ and $\gamma_2(t_r) = \frac{1}{\sqrt{1-|\boldsymbol {\beta}_2(t_r)|^2}}$, and similarly for $\mathbf {n}_1(t_r)$, $\boldsymbol {\beta}_1(t_r)$ and $\gamma_1(t_r)$. The retarded time is defined implicitly by the equation $t_r = t - \frac{1}{c}|\mathbf{r}_1(t)-\mathbf{r}_2(t_r)|$. I slightly abused notation by "factoring out" the cross product terms to avoid duplicating things, hopefully this is clear enough.
This can be generalized to a system of $N$ charged particles in a similar way. I have not performed the calculations myself due to the apparent complexity of the resulting equation of motion, but in principle one could check whether the solutions correspond to what is predicted with the electron spiraling down into the nucleus, or whether it leads to something close to elliptical orbits.
My intuition tells me that we won't observe the kind of spiralling down predicted by assuming energy is stored in the electric fields throughout all of space in this setup. Instead we consider the $\mathbf{E}$ and $\mathbf{B}$ fields as useful mathematical abstractions to simplify the expression given above into more manageable components. The interpretation of the Poynting vector would be as energy flux density which would exist only if there exist other charges that would get accelerated by the fields. In particular, the atomic stability problem would require additional charges being present near the hydrogen atom, which would obviously affect the electron's orbit through additional force terms as a multi-body problem. In that scenario, there would be an energy transfer between the electron and other charges nearby. But even then it is not clear that the electron automatically loses energy, because the nearby particles will in turn radiate and couple with the electron's motion.
We can simplify the equations above by assuming the mass of particle 2 to be very large, so that it remains effectively stationary in some inertial frame of reference, and located at the origin. Then the equations above simplify greatly and we have the following equation of motion for particle 1 (with $\mathbf{r}_2(t) = \mathbf{0}$):
$$m_1{\mathbf{\ddot{r}}_1}(t) = {\frac {\mu_0c^2q_1q_2}{4\pi}}\frac{\mathbf {r}_1(t)}{|\mathbf {r}_1(t)|^{3}}$$
which is just the equation of motion of a charged particle moving in an electrostatic potential with the Coulomb force. This limiting model is formally identical to the model of a planet orbiting a massive object under Newtonian gravitation, and we clearly have elliptical orbits. Therefore it would be incorrect to claim that classical electromagnetism predicts instabilities of the atom (in this limiting case with very massive nucleus, at least) if we don't consider the energy supposedly stored in the EM fields. Also it deals away with the self-energy problem of a charged particle (integrating the "electrical field energy density" over all spaces gives an infinite result, which would be simply a meaningless calculation since there is no actual energy stored in such a field).
I hope what I said above was sufficiently clear, and I would be curious whether solutions to the equations of motion I have described have been calculated (or approximated somewhat) to predict classical orbits of an electron around the nucleus.
Note also that I am not questioning the validity of quantum mechanics and more detailed theories of matter. I am simply wondering whether the specific problem of atomic instability supposedly predicted by classical electromagnetics arises only due to the assumed existence of electromagnetic fields carrying energy, or whether a field-free formulation of electromagnetics using the equations of motion above is also subject to this problem. I am sure there are other problems that this model cannot resolve, such as the existence of discrete atomic emission and absorption spectra. But the important observation I wanted to make is that starting from the classical Maxwell equations and Lorentz force, we can derive the Liénard-Wiechert potential, and then derive the explicit equations of motion above, and finally forget about the existence of the $\mathbf{E}$ and $\mathbf{B}$ fields. This leads to a classical model of a two-body atom with stable orbits.