How to implement a block version of MatDiagonalScale?
Stephane Aubert
stephane.aubert at fluorem.com
Tue May 6 12:20:07 CDT 2008
Hello,
I'm "playing" with badly conditioned block matrices (MPIBAIJ type)
arising from turbulent compressible fluid dynamics (RANS eqs, block size
= 5 (laminar) or 7 (turbulent)).
I tested that using bi-normalization approach to build l and r vectors,
then calling MatDiagonalScale(A,l,r), improve the accuracy and the
convergence of GMRES+ILU(0), basically by reducing the dependence to the
mesh cells size and to the variables magnitudes (in particular for k-w
turbulence model).
Now, I would like to go further and use L and R block diagonal matrices,
instead of vectors, for example to get block identity along the diagonal
of L.A.R. The sparcity of A is preserved, from a block point of view.
My first idea is to use a LU factorization of Aii, with L=Li^(-1),
R=Ui^(-1), Aii=Li.Ui.
My question is: How to compute L and R using PETSC available functions
in a clever way, instead of calling some LAPACK functions after getting
individual diagonal blocks in my own piece of code?
My actual guess is:
1. Create a new matrix D with MPIBDIAG type, from the block-diagonal
of A using MatGetSubMatrices(), to do something like MatGetDiagonal().
2. Create a PC of type PCLU using D as operators, with the sequence
PCCreate/PCSetType/PCSetOperators/PCSetUp.
3. Get Li and Ui. This is where I'm glued: I can't figure out how to
use PCGetFactoredMatrix() to get Li and Ui... Is MatLUFactor() a
better candidate?
4. Do I need to compute Li^(-1) and Ui^(-1), or is there some
"backsubstitution" functions available?
With many thanks for your suggestions,
Stef.
--
___________________________________________________________
Dr. Stephane AUBERT, CEO & CTO
FLUOREM s.a.s
Centre Scientifique Auguste MOIROUX
64 chemin des MOUILLES
F-69130 ECULLY, FRANCE
International: fax: +33 4.78.33.99.39 tel: +33 4.78.33.99.35
France: fax: 04.78.33.99.39 tel: 04.78.33.99.35
email: stephane.aubert at fluorem.com
web: www.fluorem.com
More information about the petsc-users
mailing list