Per @lucianopaz suggestion, I have renamed my new version as fast_sample_posterior_predictive(), and restored the old version.
This is not just nice for testing, but also is critical because we cannot currently vectorize the sampling for multivariate distributions.