Why do electrons in singly occupied states have parallel spin? I've been reading about Hund's rule and I realized that the Pauli exclusion principle only says that the electrons singly occupy states first, and nothing about their spins needing to be parallel. Why is this? I suspect that it's due to the symmetry of spin product states (maximum and minimum spin being symmetric and alternating with antisymmetric in between, or that there might be some part of the spin-statistics theorem that causes this. I'm mostly asking because I would expect spins to antialign due to magnetic dipole interaction. So another way to phrase my question is: why do singly occupying electrons have parallel spins instead of the seemingly more energetically favorable mixed direction spins?
Thank you very much for answering, I'm doing a paper on ferromagnetism and this will help a lot.
Edit: A clarifying example
In a 2p orbital with 2 electrons, if the electrons had the same m, then they'd need to have antiparallel spins [up down][ _ _ ][ _ _ ]. They generally don't have the same m because of reduced screening of the nucleus. In this case, why is the "maximum spin" or "maximum multiplicity" state [up _ ][up _  ][ _ _ ] preferred as opposed to [up _][down _ ][ _ _ ].
 A: After digging around for quite some time, I finally found the answer, Fermi heaps and holes!. The wiki page lists two references by Dan Dill, a webpage with animations and a PDF explaining this qualitatively, and the webpage references yet another PDF which explains it more quantitatively.
The long and the short of it is that when the electrons have the same spin (a spin symmetric state), the spatial function needs to be antisymmetric. The non-obvious consequence of this is that the probability density of the electrons being at the same location is zero. Conversely, when the electrons are in a spin antisymmetric state (antiparallel spins for two electrons) the spatial function has a maximum where the locations of electrons are the same.
To make this clear, when the spatial function is antisymmetric, it's zero  wherever electrons are in the same location, a "Fermi hole." When the spatial function is symmetric, it is amplified wherever electrons are in the same location, a "Fermi heap".
The Coulomb repulsion is much stronger generally than the magnetic dipole interaction (this is classical E&M that is covered in a first semester E&M course). As a result, this means that the overall energy is lower when the spatial function is antisymmetric, and this only happens when the spin function is symmetric. When the electron spins are added like so, you can see that the maximum total spin states are the only symmetric spin states (although this might not be obvious). This happens even though the electrons are in different orbitals.
The upshot of this is that electron Coulomb repulsion is reduced when the spins of electrons in different orbitals are parallel because they, as Dill says "stay away from one another", and Levine says "keep out of each other’s way" (in Quantum Chemistry 2013 pg. 311).
It's frustrating to me that Levine says "(recall the idea of Fermi holes)", but he never discusses them in his book. Also, Dr. Dill works through 2 electron examples, but I find that his approach could be readily extended for multiple electrons.
I hope that others (and my future self) will find this useful, instead of just taking Hund's rule for face value. Good day.
Edit: I made a plot for my paper showing the dipole and coulomb energies calculated by classical E&M, and I'll share it for those who don't know how to do the calculation.

The dipole dies off at about a $\frac{1}{|r|^3}$ rate, whereas coulomb dies off like $\frac{1}{r}$, but the coefficients matter to know where the regime transition occurs. As I claimed, the dipole interaction only matters at a very small range, less than 1% of the Hydrogen atom radius.
