Problems with LDL factorization

8 次查看（过去 30 天）
Kitic
Kitic 2014-6-7
评论： Kitic 2014-6-7
Hello,
I have a very large scale, overdetermined linear system (10^7 variables), arising from the discretization of a time-dependent PDE. I need to compute least squares estimate (A'*Ax = A'b) several times with a different right hand side b. To do this on modest scale problem (10^4 to 10^5) I factorize the LDL decomosition of A'A and then the solution is fast enough. I do not use Cholesky to avoid potential rounding errors. However, when the dimensions increase, LDL does not prouduce accurate decomposition (even with threshold set to 0.5). The system is increasingly ill-conditioned, so this may be the source of the problem. Do you have any suggestions?
Best, Srdan

请先登录，再回答此问题。

回答（2 个）

John D'Errico
John D'Errico 2014-6-7
Yes, BUT if the problem is ill-conditioned, then computing A'*A is crazy, since that will make it significantly more so. If it cannot succeed because of the ill-conditioning, what does it matter if it is faster?
I would suggest trying an iterative scheme. Perhaps a version of PCCGLS (preconditioned conjugate gradient least squares) There are codes for it in MATLAB. This algorithm does not force you to form A'*A, although you will probably need a preconditioner.
  1 个评论
Kitic
Kitic 2014-6-7
Sure, but it did make sense for the smaller scale problem (which is not as badly conditioned). I have tried warm-started iterative solvers (only the stuff provided in Matlab), and this is one of the approaches I plan to use if factorization becomes impossible. The issue is that the initial point (taken as the estimate of the previous iteration) need not be close enough for the new problem - mainly due to conditioning. However I haven't tried PCCGLS (haven't heard of it), and except Jacobi preconditioner, everything else seems rather expensive.
Thanks all for the advices.

请先登录，再进行评论。

Bill Greene
Bill Greene 2014-6-7
Have you tried solving the over-determined system directly:
x = A\b;
Computing using the normal equations, A'*A, is almost certainly making your problem more ill-conditioned.
Bill
  3 个评论
显示 1更早的评论
Bill Greene
Bill Greene 2014-6-7
I am very surprised to hear that. I would expect that mldivide is doing a sparse qr factorization for this problem and I would expect that to be faster than an ldl composition of A'*A. You can see what mldivide is doing by putting this statement before the call:
spparms('spumoni',1);
The next experiment I would suggest is to call qr directly to solve your problem. The reason I first suggested mldivide is that it takes only one line and I thought the performance would be essentially the same as calling sparse qr and then doing a back solution with R. The qr doc page has an example of this:
http://www.mathworks.com/help/matlab/ref/qr.html
Bill
Kitic
Kitic 2014-6-7
I didn't mean that solving A\b once is slower than LDL. As I said, it is a problem that needs to be solved once per (external) iteration, thus even though factorizing LDL costs more than A\b, the overall cost is less since in the subsequent iterations the problem is solved by forward-backward approach. Sparse QR is the only thing I haven't considered - so I will give it a try, thanks a lot.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Linear Algebra 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by