[petsc-users] multiply a mpibaij matrix by its block diagonal inverse

[Xiangdong]

The reason for the operation invdiag(A)*A is to have a decoupled
matrix/physics for preconditioning. For example, after the transformation,
the diagonal block is identity matrix ( e.g. [1,0,0;0,1,0;0,0,1]  for
bs=3). One can extract a submatrix (e.g. corresponding to only first
unknown) and apply special preconditioners for the extracted/decoupled
matrix. The motivation is that after the transformation, one can get a
better decoupled matrix to preserve the properties of the unknowns.

Thanks.

Xiangdong

On Tue, Feb 13, 2018 at 6:27 PM, Smith, Barry F. <bsmith at mcs.anl.gov> wrote:

>
>  In general you probably don't want to do this. Most good preconditioners
> (like AMG) rely on the matrix having the "natural" scaling that arises from
> the discretization and doing a scaling like you describe destroys that
> natural scaling. You can use PCPBJACOBI to use point block Jacobi
> preconditioner on the matrix without needing to do the scaling up front.
> The ILU preconditioners for BAIJ matrices work directly with the block
> structure so again pre-scaling the matrix buys you nothing. PETSc doesn't
> have any particularly efficient routines for computing what you desire, the
> only way to get something truly efficient is to write the code directly
> using the BAIJ data structure, doable but probably not worth it.
>
>   Barry
>
>
> > On Feb 13, 2018, at 5:21 PM, Xiangdong <epscodes at gmail.com> wrote:
> >
> > Hello everyone,
> >
> > I have a block sparse matrices A created from the DMDA3d. Before passing
> the matrix to ksp solver, I want to apply a transformation to this matrix:
> namely A:= invdiag(A)*A. Here invdiag(A) is the inverse of the block
> diagonal of A. What is the best way to get the transformed matrix?
> >
> > At this moment, I created a new mat IDA=inv(diag(A)) by looping through
> each row and call MatMatMult to get B=invdiag(A)*A, then destroy the
> temporary matrix B. However, I prefer the in-place transformation if
> possible, namely, without the additional matrix B for memory saving purpose.
> >
> > Do you have any suggestion on compute invdiag(A)*A for mpibaij matrix?
> >
> > Thanks for your help.
> >
> > Best,
> > Xiangdong
> >
> >
> >
> >
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20180214/0f87b96d/attachment.html>