We have plenty of good reasons to expect that standard quantum electrodynamics accurately predicts the electronic structure of all neutral atoms in their ground states, including the Madelung rule and its exceptions, but I don't know how completely this has actually been checked. First-principles calculations are difficult, to say the least, and approximations that make them tractable can also introduce large errors.
The rule and its exceptions were originally established experimentally [1]. Are they understood theoretically? We should distinguish between two different questions:
First question: Have they been derived, either analytically or numerically, from some adequate approximation to quantum electrodynamics?
Second question: Are they understood intuitively?
The distinction is important. We have a proof that Fermat's last theorem follows from the laws of arithmetic (first question), but I doubt anybody understands it intuitively (second question). The question being asked here is the first one.
Here's an excerpt from 2019 [2]:
the Madelung rule has not yet been derived from quantum mechanics or other fundamental physical principles. In 1969, on the 100th anniversary of the periodic table, chemist Per-Olov Löwdin declared this derivation to be one of chemistry’s major theoretical challenges. It still is, 50 years on.
That excerpt seems to say that the answer is "no," period. Not understood intuitively, and not understood numerically. However, I'm not sure that's what the author meant. That 2019 commentary cites a 2009 paper [3] which says this:
no general convincing theoretical derivation exists [of the $(n+\ell,n)$ rule]. That was already noted by Löwdin a long time ago... [It] has been derived, however, in the framework of the Thomas-Fermi approximation [refs], for a very special case [namely for parts of rows 4 and 5 in the periodic table]. ... One can create an atomic operator by intelligent design that just yields the intended Madelung eigenvalue sequence [ref], or design a symmetry group that just fits to the Madelung AO [atomic orbital] pattern [ref]. Another attempt to create the Madelung sequence refers to the Sturm–Liouville node-energy theorem [ref]. The number of nodes of an $n\ell$ AO in the domain $0<r<\infty$, $0<\text{angle}<2\pi$ is $n + \ell$, indeed. However, the Sturm–Liouville theorem holds only for one dimension; already for two dimensions, the eigen-energies do not strictly increase with the number of nodes of the eigenfunctions.
This excerpt also says that the answer is "no," but in contrast to the previous one, this excerpt seems to be merely saying that the rule is not understood intuitively, because it mentions simple things like the Thomas-Fermi approximation but doesn't mention more sophisticated approximations that might give good results numerically even if they don't convey much intuition.
Here's an example of what looks like a positive numerical result. Table 1 in a 1967 report [4] shows total energies of atoms calculated using the Hartree-Fock approximation. For Cerium ($Z=58$, one of the exceptions), the table shows energies of two different electronic configurations, and the rule-breaking configuration for Ce is just barely lower in energy than the rule-respecting configuration. This seems to correctly predict that Ce is an exception to the Madelung rule — and this is from 1967. That was a long time ago, and surely we can do even better today, but how big are the errors in these approximations? I didn't notice a clear error analysis in the report. (To be fair, maybe I just missed it. I only skimmed through it.) Also, Hartree-Fock calculations use a variational approach in which the set of candidate orbitals must be assumed, with some adjustable parameters that the computer can vary to minimize the energy. That might be why the 1967 report showed energies for two different electronic configurations for some of the atoms (including Cerium): to see which of the two assumed configurations gave the lowest energy. The strength of these results depends partly on how much the computer is allowed to vary the orbitals. All of these things need to be considered before we can assess how compelling the results are.
By the way, the text that precedes the table says this:
To date the calculations have been carried out for all the elements to atomic number 103 in their ground state configuration, as well as various other alternate configurations. ... The configurations of the first 103 elements listed are believed to be ground-state configurations. Since the ground-state configuration of the heaviest elements (97 to 102) have not been determined by analysis of the optical spectra, calculations have been made for likely alternate configurations of these elements. Cerium and terbium have also been included with alternate configurations...
This seems to say that, at least in 1967, the ground-state configurations printed on periodic-table posters were a mix of empirically-determined (for $Z\leq 96$) and theoretically-determined (for $Z\geq 97$) configurations.
So... are the Madelung rule and all of its exceptions understood theoretically, in the sense that we have clear confirmation that quantum electrodynamics correctly predicts them? I'm not sure, but the answer seems to be something like... "partly, but not yet."
References:
[1] Page 157 in "What and how physics contributes to understanding the periodic law" (http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.23.5192&rep=rep1&type=pdf)
[2] "Can quantum ideas explain chemistry’s greatest icon?" (https://media.nature.com/original/magazine-assets/d41586-019-00286-8/d41586-019-00286-8.pdf)
[3] "Some solved problems of the periodic system of chemical elements" (https://onlinelibrary.wiley.com/doi/abs/10.1002/qua.22277)
[4] "Atomic Structure Calculations I. Hartree-Fock Energy Results for the Elements Hydrogen to Lawrencium" (https://www.osti.gov/servlets/purl/4297173)