I gave in and looked up the conditional-multivariate-normal formula, and implemented a gibbs sampler based on it

it's wrong.
of course it's wrong.
now i get to debug it.
comments and wtfs
one bug squashed: make sure to init Y to a *floating point* array
two bug squashed: it was making a giant list of pointers to the same array(), so the entire result had the value of the last sample taken
bugthree: chol() the appropriate place.

at least for j=0, this gives a distribution matching multinormal(),
which basically means I've gotten all my syntax in the right order,
but not that the gibbs part is working, since in the j=0 case the
gibbs sampler collapses to just a really slow standard multinormal sampler.

the new guards are since chol() is less tolerant of empty vectors than dot().

I still don't understand how the *derivation* of this came about
 and therefore i can't debug it any further :(
LOL SHOT MYSELF IN THE FOOT
I think multigibbs is working. This commits adds basic exercises of it.

TODO: test on bigger dimensions, and with and without covariance between the parts
and also derive what the covariance matrix should be for the hiearchical normal...
TODO:
TODO: derive the conditioning formulas; try to see if what the assignment says has any sensible relation at all
and of course, code cleanup and porting to Rcpp needs to happen
Increase the test case dosage.

Because eyeballing an 18x18 covariance matrix is stupid hard,
I use max(abs(C - E[C])) as a summary statistic.
Sometimes the gibbs sampler ends up closer to E[C] than the multinormal!!
Oooh, it apppppears to converge  even--even!--if the partition changes during sampling.
better testing
Ah, but changing partitions is a needless complication.

I suspect it slows down convergence, though i actually have no evidence for this
and I would love to know if it *doesn't*
Stubbing out the multivariate hierarchical sampler.

I have a baseline function which samples direct from the definition,
and several others which should do the same with linalg.

TODO: import my paper derivations to digital ones.
test the gibbs hierarchical normal sampler. it checks out.