フィルターのクリア

Quick method in solving symmetric linear equation

11 ビュー (過去 30 日間)
Ziwen Gu
Ziwen Gu 2023 年 10 月 22 日
コメント済み: Ziwen Gu 2023 年 10 月 23 日
In simulation of PDE, I have to solve a linear equation: Ax = b for many times( with different b, iteratively ) .
The matrix A is symmetric with size of about 1e6*1e6 , and doesn't change in time. So I use cholesky decomposition
R = chol(A);
and solve two triangle equation
y = linsolve(R',b,opts.LT = true);
x = linsolve(R ,y,opts.UT = true);
at each iteration.
More detailed, matrix A is like:
A = M + S*inv(M)*S + S ;
where M is a diagonal matrix( lump mass matrix in finite element method ) , S is a symmetric stiff matrix.
My question is : is there a better way to solve this equation?( quick and stable ) All numerical method are acceptable, for example preconditioner, iterative method...
Thanks for your help!

採用された回答

Bruno Luong
Bruno Luong 2023 年 10 月 22 日
編集済み: Bruno Luong 2023 年 10 月 22 日
If your S is sparse you should look at iterative methods. For symmtric case it is recommended to use pcg. To take further advantage you might want to find a cheap preconditioner. The problem is you have to iterate for new rhs.
The cholesky decomposition will fill the entire half matrix, so it will not be cheap. I'm even surprised you can even do it with 1e6 matrix.
  1 件のコメント
Ziwen Gu
Ziwen Gu 2023 年 10 月 22 日
Yes it's very stupid, I find Matlab uses cholesky decomposition when I solve A\b, so I adapt this method...

サインインしてコメントする。

その他の回答 (1 件)

Christine Tobler
Christine Tobler 2023 年 10 月 23 日
You can definitely use iterative methods such as pcg. The success of this usually depends on how well-conditioned your matrix is, or how good of a preconditioner you can come up with.
In terms of direct methods, consider using permuted Cholesky (the three-output syntax). The permutation allows reduction of fill-in in a sparse direct factorization in many cases - meaning you need less storage and less computation for the factor R, as it contains fewer nonzero elements.
You can also try the decomposition object (which uses permuted Cholesky, and additionally stores the factorization of A in the same internal format used in backslash, instead of converting it to a generic sparse matrix):
dA = decomposition(A);
x = dA \ b;
  1 件のコメント
Ziwen Gu
Ziwen Gu 2023 年 10 月 23 日
Thanks a lot for your reply! I further study the structure of A, its condition number is about 10(very well-conditioned), and I use incomplete cholesky as preconditioner, it works well.
I will try another two methods in spare time. Thanks again for your help!

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeMathematics and Optimization についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by