I'm not sure that I understand what this debate is about, and it's hard to read the entire thread. Since
@Charles Link has asked in a PM, let me try to put in another point of view to the debate.
First of all, it's of course clear that anything to understand about the "matter part of electromagnetism" from first principles has to do with quantum theory. For everyday matter in many cases you can use non-relativistic many-body theory, which simplifies some things, particularly if it comes to spin, because spin separates from position and momentum observables in the sense that in non-relativistic quantum physics spin is an additional angular momentum commuting with position and momentum observables. Thus you have a pretty simple concept of total angular momentum as the sum of orbital and spin angular momentum (the only technical problem being the understanding of this sum in terms of Clebsch-Gordan coefficients and all that, which is a somewhat more complicated topic of the QM 1 or maybe QM 2 lecture).
The closest intuitive picture of spin in classical physics are twofold: (a) the spin, i.e., intrinsic rotation of a "rigid body" as part of it's total angular momentum and (b) the magnetic moment of current distributions in a compact region of space seen from a far distance, so that the extension of the region can be negelected compared to the length scale where the magnetic field varies significantly. The latter picture is the one I think is most intuitive for introducing spin in the QM 1 lecture (this I've just done before Christmas in my theoretical-physics lecture for high school-teacher students, and I think it has been quite well understood).
Of course, one must emphasize that these are pictures, and the only correct description we have today is quantum mechanics. And all I could do in my QM 1 lecture is to thoroughly discuss the Stern-Gerlach experiment using quantum theory and only quantum theory in an idealized approximation for the interactin of an external static magnetic field with an uncharged "particle" with an magnetic moment (of course the historical example with a silver atom is appropriate, but one can also use neutrons as an example of an even more particle-like object too). This is a very interesting example for a quite simple analytically solvable dynamical (!) problem for a measurement process, which to my surprise seems not to be present in the textbook literature. I'm about to write a paper on it for AJP or EuJP as soon as I find the time after the semester.
To understand a permanent magnet, however, is a paradigmatic example for a collective many-body effect (something I cannot teach in a first lecture on quantum mechanics due to lack of time), i.e., you have a collection of spins interacting via their associated magnetic moments, and the most simple models are the often discussed Heisenberg and Ising models. It's an entire plethora to discuss fundamental many-body physics, including the "exchange forces" due to the indistinguishability of (fermionic) particles (microscopic level) and the concept of effective quantum field theory, using the mean-field approximation and spontaneous symmetry breaking of (global) rotational symmetry, leading to an example for Ginzburg-Landau theory etc. etc.
Of course, one can also treat permanent magnets on a completely classical level. Then of course the treatment of the matter part of electromagnetism completely reduces to the use of phenomenological parameters. In the usual E&M lecture one restricts oneself to the "linear-response approximation", which however goes pretty far towards a realistic description of the electromagnetism of macroscopic matter (see the excellent treatment of the subject in Vol. 8 of Landau&Lifshitz). The reason, of course, is that usually the electromagnetic fields imposed on the matter in everyday life are small compared to the inneratomar fields which hold the matter together, over which one "coarse grains" to the relevant macroscopic degrees of freedom.
For a permanent magnet in the most simple model treatment you simply assume some plausible distributions of magnetization (model of a "hard ferromagnet"). Then the task is to solve Maxwell's equations with the appropriate boundary conditions, in the static case for ##\vec{B}## and ##\vec{H}## with the constituent equation ##\vec{B}=\mu_0 (\vec{H}+\vec{M})## (SI units). Since then there are no free (exernal) currents, the Maxwell equations for the static case reduce to
$$\vec{\nabla} \times \vec{H}=0, \quad \vec{\nabla} \cdot \vec{B}=0\; \Rightarrow \; \vec{\nabla} \cdot \vec{H}=-\vec{\nabla} \cdot \vec{M}.$$
One can in this case thus use a scalar potential for ##\vec{H}##, ##\vec{H}=-\vec{\nabla} \phi_m## and then solve for the Poisson equation in the usual way:
$$\phi_m(\vec{x})=\int_{\mathbb{R}^3} \mathrm{d}^3 \vec{x}' \frac{-\vec{\nabla}' \cdot \vec{M}(\vec{x}')}{4 \pi |\vec{x}-\vec{x}'|}.$$
Using Gauss's integral theorem and the fact that at infinity ##\vec{M}## vanishes, one can transform this into the more useful form
$$\phi_m(\vec{x}) = \int_{\mathbb{R}^3} \frac{\vec{M}(\vec{x}') \cdot (\vec{x}-\vec{x}')}{4 \pi |\vec{x}-\vec{x}'|^3} = -\vec{\nabla} \cdot \int_{\mathbb{R}^3} \mathrm{d}^3 \vec{x}' \frac{\vec{M}(\vec{x}')}{|\vec{x}-\vec{x}'|}.$$
Another, equivalent way is to use the vector potential for the magnetic field, $$\vec{B}=\vec{\nabla} \times \vec{A}$$ and
$$\vec{\nabla} \times \vec{H}=\frac{1}{\mu_0} \vec{\nabla} \times \vec{B}-\vec{\nabla} \times \vec{M}=0.$$
Then, imposing the Coulomb-gauge condition ##\vec{\nabla} \cdot \vec{A}=0## one finds
$$\Delta \vec{A}=-\mu_0 \vec{\nabla} \times \vec{M},$$
i.e., the magnetization is equivalent with a current density,
$$\vec{j}_m = \vec{\nabla} \times \vec{M}$$
as far as the magnetic field ##\vec{B}## is concerned.
If you use homogeneous magnetization as a simple model, this effective current reduces to the surface current, corresponding to the jump of the magnetization at the surface of the permanent magnet. Which technique to use to calculate the fields, of course, depends on the problem. Usually the scalar-potential method is more straight forward, because one has not so much problems to struggle with the singulartities of the magnetization and the evaluation of the effective surface current.
You find the example of a homogeneously magnetized sphere, where both methods can be used to evaluate the magnetic field analytically in my lecture notes on E&M for high school-teacher students (in German, but I think with a high enough "equation density" to be understandable also by non-German speakers):
https://th.physik.uni-frankfurt.de/~hees/publ/theo2-l3.pdf (Sect. 3.3.3).
## \\ ## Much of my working years were spent doing electro-optic experiments, but had I been doing all kinds of experiments with magnets and devices employing magnets, I think think I would have found it equally exciting.
