stat_util.py function truncated_multivariate_normal_fast_rvs is currently written as a rejection sampler with pre-computed cholesky factor. This rejection sampler is too slow in higher dimensions. I should implement an MCMC sampler, to speed this up.