The notion of primitive idempotents isn't that common in the literature. You have to be a bit of a conosewer to run into it. I'm sure there's references in the literature as it's kind of obvious, if you spend enough time playing with the things. My website, www.densitymatrix.com has a connection to Frank Porter's (Cal Tech) class notes on quantum mechanics. If you click on that link, you will see a link labeled "Physics 125c Course Notes Density Matrix Formalism" or similar. Page 11 of those class notes will verify what I've said about "Hermitian primitive idempotents" being a characterization of pure density matrices. Non Hermitian density matrices are even more rarely discussed than the primitive idempotents and pure density matrices. I seem to recall seeing some reference at a conference on non Hermiticity.
Re connection between pure density matrices and MUBs. The powerful thing about MUBs is that they define a set of quantum states that interact. The usual way one must define a set of mathematical objects (spinors) to represent the fermions, and then define a set of stuff that allows them to interact. It is in defining these TWO things that it is inevitable that you have to choose a symmetry. With states chosen from a complete MUB set, the interactions are already built into them.
The interactions between MUB basis states are trivial in that they all have the same amplitude, but the phases are not completely trivial and this gives them a structure which can be analyzed. As an example of the non triviality of the phases of pure density matrices chosen from the simplest (Pauli algebra) MUB, consider the following product of pure density matrices:
(1+x)/2 (1+y)/2 (1+z)/2 (1+x)/2
where "x", "y", and "z" stand for the three Pauli spin matrices. You will get a complex multiple k of the pure density matrix (1+x)/2. The complex phase that gets picked up in that product is a Berry-Pancharatnam or quantum phase, not one of the arbitrary complex phases that confuse spinor calculations. The magnitude of k is 1/8 because there are three transitions between quantum states with transition probabilities of 1/2. The phase of k is -pi/4 or +pi/4, I forget.
And what does the above product represent physically in quantum mechanics? Because these are operators, and because they are written in geometric language i.e. x, y, z, they have more meaning than the equivalent spinor calculation which is:
( +x | +y ) ( +y | +z ) ( +z | +x ). I hope that parsed.
The meaning of the product "(1+x)/2 (1+y)/2 (1+z)/2 (1+x)/2" is:
(a) It is a product of projection operator for a quantum state; this is a history of a particle. See the interpretation of quantum mechanics known as "Consistent Histories" for more.
(b) It represents a sequence of measurements of a quantum state. See Schwinger's "Measurement Algebra" for more on this. His papers are linked at another one of my websites, www.measurementalgebra.com.
(c) It represents a sequence of four Stern-Gerlach filters set up so as to allow only spins oriented in the x, z, y, and finally x directions to pass. That is, this is a representation of a sequence of four polarizing filters. The phase is the phase that a particle passing all four filters would possess. From a field theory point of view, each Stern-Gerlach filter is a source of gauge bosons that interact with a beam of fermions. The point here is that with a Stern-Gerlach filter the physicist is not thinking of the effect of the beam on the filters, but only the reverse effect, how does the filter act upon the beam.
(d) The individual terms like (1+x)/2 represents the field that is present in the Stern-Gerlach field. That is, a Stern-Gerlach filter oriented in the +x direction has a magnetic field that is inhomogeneous in the +x direction. This is a hint that there is something going on here that has something to do with gauge bosons and the geometry of space-time.
(e) Finally, most important to me, and getting back to the MUB theory, the product represents a spin-1/2 fermion that emits three gauge bosons (in the quantum information limit where we just keep track of 1 quantum bit for the particle and don't keep track of what happens to the object that absorbed the gauge boson). The first boson emitted is one that converts a (1+x)/2 fermion to a (1+z)/2. The second converts a (1+z)/2 to a (1+y)/2, etc. As such, each of these products, i.e. "(1+x)/2 (1+z)/2", represents a gauge boson in the quantum information limit. You can think of these sort of like the six off diagonal gluons; so a red/green gluon converts a green to a red.
The gluon analogy can be taken literally. With the Pauli MUBs, there are a total of six states (1+x)/2, (1+y)/2, (1+z)/2, (1-x)/2, (1-y)/2, and (1-z)/2. These are just enough to represent red, green, blue, anti-red, anti-green, and anti-blue, respectively. Instead of representing the color states as orthogonal states with a gauge boson that interacts between them, I'm representing the color states as states that are not orthogonal; the non orthogonality automatically defines the interaction.
This use of MUBs gives a derivation of Koide's mass formulas for the hadrons and is the subject of my the paper it appears I will release first. I've finished off the mesons this morning and have just started classifying the baryons. Tommaso Dorigo has kindly undertaken to waive the submission fee at Phys Math Central.
This (e) interpretation is the one that is necessary to use to see how E8 falls out of the density matrix formalism. Calculations where we ignore the gauge bosons are not that uncommon in QM. An example is the bound state of hydrogen and Schroedinger's equation. We think of the electron as emitting photons without thinking (at first order) of how those photons are absorbed by the nucleus. This simplifies the problem to one where we watch only one of the participants in the dance. You can do the same thing with a meson, whatever is emitted by one quark has to be absorbed by the antiquark. (Of course some gluons are reabsorbed by the quark but we just sum over those and look at the quark only in terms of what it exchanges with the antiquark.)
It is in ignoring what happens to the gauge bosons that you end up with E8 from density matrix formalism. The matrices you work with have fermions down the diagonal and the gauge boson interactions off diagonal. This is simpler than it probably sounds. It's just consistency relations; there has to be just as much stuff becoming red as there is stuff that used to be red becoming something else. The E8 is only approximate because the assumption that you can ignore the gauge bosons is only approximate.
From what I can recall, the Hadamard matrices come up when you try and find complete sets of MUBs. So far, I understand that as a purely mathematical endeavor and don't have an interpretation of it. I would guess that if I was hanging around more with Marni Sheppeard she'd eventually get me to see it differently.
Part of my problem with associators is that I don't see how to give them a physical interpretation. Getting back to (c) above, the product is interpreted physically as a sequence of polarization filters. For this, clearly associativity applies and the associator is trivial.
While the octonions are not associative, my understanding of them is that they are almost associative in the sense that one only picks up a negative sign when you fool around with the parentheses on a product of basis vectors. That smells to me like something that will go away when one converts from a state vector / spinor representation to a density matrix representation. So in that sense, my feelings about octonions is better when they are used for state vector representations than density matrix as I can see a way of eventually rescuing them with a direct physical interpretation.