Is the Mendeleev table explained in quantum mechanics?

I am not offended by the suggestion that physicists should follow the standards of mathematical proof, but I think this suggestion and the phrasing of the question demonstrate a lack of understanding of how physicists think about such things and more importantly why they put such little emphasis on axioms.

In my view it is rarely useful to think of physics as an axiomatic system, and I think this question reflects the difficulty with thinking of it as such. A different question, which is much more in tune with a physicist's point of view, would be to ask what physical description is required to explain various features of the structure of atoms as reflected in the periodic table at a prescribed level of accuracy. Until you specify what features you want to understand, and at what level of accuracy, you don't even know what the correct starting point should be. If you want just the crudest structure of the periodic table, then indeed non-relativistic quantum mechanics along with the Pauli exclusion principle will give you the rough structure as described in any standard QM textbook. If you want to understand the detailed quantum numbers of large atoms then you have to start including relativistic effects. Spin-orbit coupling is one of the most important and its effects are often summarized by a set of Hund's rules which are described in many QM textbooks or physical chemistry textbooks. If you want very accurate numerical values for ionization energies or the detailed structure of wave functions then one must do hard numerical work which probably becomes impossibly difficult for large atoms. As you ask for greater and greater precision you should eventually use a fully relativistic description. This is even harder. The Dirac equation is not sufficient, one cannot restrict to a Hilbert space with a finite number of particles in a relativistic quantum theory, and bound state problems in Quantum Field Theory are notoriously difficult. So as one asks more detailed and more precise questions, one has to keep changing the mathematical framework used to formulate the theory. Of course this process could end and there could be an axiomatic formulation of some ultimate theory of physics, but even if this were the case this would undoubtedly not be the most useful formulation for most problems of practical interest.

There is some rigorous work by Goddard and Friesecke on this, see

http://www.ma.hw.ac.uk/~chris/icms/GeomAnal/friesecke.pdf

My understanding is that even getting accurate numerics for the Schrodinger equation becomes very difficult once one has more than 10 or so electrons in play. The one regime where we do seem to have good asymptotics is when the atomic number is large but the number of electrons are small (i.e. extremely highly ionized heavy atoms).

At any rate, the foundations of the periodic table are pretty much uncontested (i.e. N-body fermionic Schrodinger equation with semi-classical Coulomb interactions as the only significant force). The main difficulty is being able to solve the resulting equations mathematically (or even numerically).

I doubt any answer will be satisfactory. My opinion is that we are still very far from a mathematical justification. If we accept the mathematical foundations of quantum mechanics, and if we make the approximation that the nucleus of the atom is just one heavy thing with $N$ positive charges, then the motion of the $N$ electrons is governed by a linear equation (Schrödinger) in ${\mathbb R}^{3N}$. The unknown is a function $\psi(r^1,\ldots,r^N,t)$ with the property (Pauli exclusion) that it has full skew-symmetry. For instance, $$\psi(r^2,r^1,\ldots,r^N,t)=-\psi(r^1,r^2,\ldots,r^N,t).$$ In practice, we look for steady states $e^{i\omega t}\phi(r^1,r^2,\ldots,r^N)$. Then $\omega$ is the energy level.

Because of the very large space dimension, one cannot perform reliable calculations on computer, when $N$ is larger than a few units. One attempt to simplify the problem has been to postulate that $\phi$ is a Slatter determinant, which means that $$\phi(r^1,r^2,\ldots,r^N)=\|a_i(r^j)\|_{1\le i,j\le N}.$$ The unknown is then an $N$-tuple of functions $a_i$ over ${\mathbb R}^3$. Of course, we do not expect that steady states be really Slater determinants; after all, the Schrödinger equation does not preserve the class of Slater determinants. Thus there is a price to pay, which is to replace the Schrödinger equation by an other one, obtained by an averaging process (Hartree--Fock model). The drawback is that the new equation is non-linear. Such approximate states have been studied by P.-L. Lions & I. Catto in the 90's.

Update. Suppose $N=2$ only. If we think to $\phi$ as a finite-dimensional object instead of an $L^2$-function, then it is nothing but a skew-symmetric matrix $A$. Approximation à la Slater consists in writing $A\sim XY^T-YX^T$, where $X$ and $Y$ are vectors. In other words, one approximate $A$ by a rank-two skew-symmetric matrix. The approximation must be in terms of the Hilbert-Schmidt norm (also named Frobenius, Schur): this norm is natural because of the requirement $\|\phi\|_{L^2}=N$. If $\pm a_1,\ldots,\pm a_m$ are the pairs of eigenvalues of $A$, with $0\le a_1\le\ldots\le a_m$, then the best Slater approximation $B$ satisfies $\|B\|^2=2a_m^2$, $\|A-B\|^2=2(a_1^2+\cdots+a_{m-1}^2)$. Not that good. Imagine how much worse it can be if $N$ is larger than $2$.