First, your explanation at the beginning isn't right. The centrifugal potential is counted in the hydrogen atom calculation but the energy of the state only turns out to be a function of $n$, not $l$ or $m$.
The preference for $l=0$ states only applies to multi-electron atoms and it boils down to electron-electron interactions. Electrons with lower values $l$ are closer to the nucleus so they feel the full electric charge of it which gives these bound states a higher, more negative potential energy. On the other hand, electrons with higher $l$ are further and the nuclei charges are partially "screened" or "shielded" by other electrons, so the "gain" for the binding energy is smaller.
More detailed questions about the ordering of the electron states and their spin in multi-electron atoms are answered by Hund's rules.
The energy difference between $s$ and $p$ electron states is substantial, comparable to other energy differences between major levels, i.e. nearly an electronvolt (values of course depend on the exact atoms and states), and may only be suppressed for higher $Z$ (the relative difference i.e. ratio is suppressed, the absolute energy difference is large, too).
The energy of all levels of an atom may always be measured from the spectra. Each atom only has one truly stable state, the ground state with the lowest energy, and there exist transitions from all higher states that emit photons of measurable frequencies.
The Lamb shift is a very subtle effect in the hydrogen atom implied by QED. It only makes qualitative sense to calculate it for the hydrogen atom. For multi-electron atoms, the degeneracy (originally equal energies of levels) is lifted (removed) by much stronger effects than the Lamb shift. So although the Lamb shift is present even for multi-electron atoms, it almost never influences the qualitative ordering of levels.
For non-relativistic hydrogen atom, as I said, the energy depends on the quantum number $n$ via $-1/n^2$ only. It doesn't depend on $l,m,s$. The Dirac equation adds some small, relativistic correction dependence on $j$ which comes from $\vec L+\vec S$ but keeps some accidental degeneracy for different values (pairs) of $l$. This degeneracy is removed in QED by the Lamb shift. Those things are calculable, fine, and predictable for the hydrogen atom only. More complicated atoms already lift the degeneracy at the non-relativistic level, before we even consider the Dirac equation or QED, and Dirac equation or QED only make small corrections to a "maximally split" spectrum.