First, the idealized clean non-relativistic problem with the $1/r$ Coulomb potential is integrable: the exact wave functions and energy eigenvalues may be precisely written down using elementary functions. It is exactly true that the eigenvalues of the bound states are $-13.6\,{\rm eV}/n^2$ etc.
Now, the first simple correction is the proton motion. This is easily incorporated by switching from the electron mass to the reduced mass $m_e m_p / (m_e+m_p)$. The corrections are comparable to 0.1%.
The Dirac equation gives relativistic corrections. Because the speed of the electron in the atom is comparable to the fine-structure constant $\alpha\sim 1/137.036$ times the speed of light, the relativistic corrections start as relative corrections of order $\alpha^2\sim 10^{-4}$ or so. The Dirac equation also splits some levels so it's no longer true that there's no energy dependence on $\ell$. Instead, a small dependence on $j$ emerges etc.
Aside from these, a simple correction is the so-called hyperfine structure which results from the interaction between the magnetic moment of the electron and the nucleus.
One also needs to include even higher-order corrections resulting from the fact that the proton motion and the relativistic corrections (assuming a static proton) aren't quite independent of each other – the system isn't quite linear. Those require complex calculations in Quantum Electrodynamics. We're talking about the $10^{-6}$ precision here and better.
Quantum Electrodynamics corrections also produce loops (in the sense of Feynman diagrams). The largest effect of this type is the Lamb shift, an effect of a virtual photon that gets emitted and reabsorbed by the atom in some state. This splits (by a very small amount) some levels that were previously degenerate even according to the Dirac equation. In principle, Quantum Electrodynamics gives a systematic formula with an infinite sequence of corrections that include relativity and quantum loops.
We have still assumed that the nucleus is point-like. The leading correction from the "complex nuclear physics" arises due to the extended geometry of the nucleus. It means that the potential felt by the electron isn't $1/r$ all the way up to $r=0$. Instead, at $r=r_N$, the nuclear radius (radius of proton etc.), the previously divergent $1/r$ basically develops a plateau. It's a nice exercise (we had it on an exam) to calculate how this "nonzero size of nucleus" effect affects the energy eigenvalues.
It may be calculated by the perturbation theory and it primarily affects $\ell=0$ states because for those states, the electron is likely to be "at the nucleus", while for a higher $\ell$, the wave function decreases as $r^\ell$ near $\ell=0$. The size-of-proton correction is the only one by which the nuclear effects influence the atomic spectrum measurably.
In principle, one should calculate the exact dynamics of the quarks and gluons in the proton etc. and incorporate all possible (very incomplete) entanglement/correlations between the electrons and properties of the quarks. But the tiny effects would be unmeasurably tiny: it's enough to assume that the proton and electron are basically independent and the energy of the proton is whatever it is without the electron.
In this approximation which is enough for the highest-precision measurements we can make, again, the nucleus only affects the energies by its motion (via the reduced mass), via the interaction of the magnetic moments (the hyperfine structure), and via the truncation of the Coulomb's potential (nonzero radius of the proton).
