The derivatives for log_det and chol(Q)^{-1} x can be computed analytically and reduced to sparse matrix calculations. This is why I’m a bit confused about the previous statement of needing dynamic var arrays – if you’re just exposing something like a sparse quadratic form then why can’t the calculation and it’s derivative all be done in doubles and then wrapped in a vari?
To match the mat spec, where the chol is diff’d separately.
But if we only need to get a multi normal working than that’s maybe easier. (Maybe = let’s wait for the maths!)
Because of the difficulty of full on sparse matrix support, in my opinion the fastest path to useful contributions will be functions that wrap sparse operations in which the user has to only pass the initial sparse structure. This was the motivation, for example, behind our recommendations for GP covariance functions, https://github.com/stan-dev/stan/wiki/Adding-a-Gaussian-Process-Covariance-Function.
A sparse multi_normal would be an excellent starting point.