Right now, we don't support extraction of the sparse Cholesky factor, but I think that it will be supported soon, see
https://github.com/JuliaLang/julia/pull/10402 I don't fully understand the structure of your problem, but I think that you can only save flops by using the Cholesky factor when B=D. In that case solving with the whole factorization and D should be okay from an efficiency point of view. Alternative, you'll probably have to wait for the pull request to be merged. 2015-04-27 15:47 GMT-04:00 matt <[email protected]>: > There are a small number of different Bs and different Ds, and I have to > compute the quadratic form for a number of different combinations. I can > probably do things the way you suggested, but I would lose the nice > symmetry that I have in my current R code. (I am trying to switch to Julia > for speed.) > > > On Monday, April 27, 2015 at 2:30:00 PM UTC-5, Andreas Noack wrote: >> >> If B and D are different then why is it not okay to calculate x = C\D and >> then B'x afterwards? >> >> 2015-04-27 15:24 GMT-04:00 matt <[email protected]>: >> >>> I would like to compute multiple quadratic forms B'*C^(-1)*D, where B, >>> C, and D are sparse, and C is always the same symmetric positive matrix. >>> For fast computation, I would like to do this by computing and re-using the >>> Cholesky factor: C=LL'. However, chol(C) does not work (because C is >>> sparse), and cholfact(C)\D returns C^(-1)*D instead of C^(-1/2)*D. >>> >>> Minimal working example: >>> C=4*speye(3); # sparse matrix >>> B=sparse(ones(1,3)) # sparse unit vector >>> chol(C) # returns "ERROR: Use cholfact()instead of chol() for sparse >>> matrices..." >>> L=cholfact(C); # works, but returns Cholesky decomposition instead of >>> cholesky factor >>> L\B' # would like this to be vector with elements .5 instead of .25 >>> >>> Thank you! >>> >> >>
