Just got burned by a silly Matlab gotcha. Suppose you have a diagonal matrix D and you want to compute the inverse square root matrix:
Disqrt = diag(1./sqrt(diag(D))
But this will be dense!
Instead use
Disqrt = diag(sqrt(diag(D).^-1)
Or maybe
Disqrt = diag(diag(D).^-0.5)
Not sure if there's an accuracy difference (hopefully not).