[petsc-users] Matrix inverse and LU decomposition

[Hong Zhang]

You need call

MatSetType(I,MATDENSE);
MatSetType(Dinv,MATDENSE);

See petsc-dev/src/mat/examples/tests/ex1.c (attached) on how to use
MatMatSolve().
Can you send me your complete code so I can check why Dinv is empty?

Hong

On Wed, Mar 30, 2011 at 8:31 AM, Bartłomiej W <bartlomiej.wach at yahoo.pl>wrote:

> Hello,
>
> I want to obtain the inverse of a matrix, even though it's not advised by
> the FAQ ;-) and I have a problem that earlier responses didn't solve for me.
> Simply speaking : the inverse matrix that is a result of  DI=Dinv as told by
> the FAQ, is empty for me. Here is my code(the significant part) :
>
>   ierr = MatCreate(PETSC_COMM_WORLD,&D);CHKERRQ(ierr);
>   ierr = MatSetSizes(D,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr);
>   ierr = MatSetFromOptions(D);CHKERRQ(ierr);
>
>   ierr = MatCreate(PETSC_COMM_WORLD,&Dinv);CHKERRQ(ierr);
>   ierr = MatSetSizes(Dinv,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr);
>   ierr = MatSetFromOptions(Dinv);CHKERRQ(ierr);
>
>   ierr = MatCreate(PETSC_COMM_WORLD,&I);CHKERRQ(ierr);
>   ierr = MatSetSizes(I,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr);
>   ierr = MatSetFromOptions(I);CHKERRQ(ierr);
>
> (...)
>
>   ierr = MatAssemblyEnd(D,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
>   ierr = MatAssemblyEnd(I,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
>   ierr = MatAssemblyEnd(Dinv,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
>
>   MatView(D,PETSC_VIEWER_STDOUT_SELF);
>
>   ierr = MatGetOrdering(D,  MATORDERING_NATURAL,  &perm,  &iperm);
> CHKERRQ(ierr);
>   ierr = MatFactorInfoInitialize(&info); CHKERRQ(ierr);
>   ierr = MatLUFactor(D, perm, iperm, &info); CHKERRQ(ierr);
>
>   MatView(D,PETSC_VIEWER_STDOUT_SELF);
>   MatView(I,PETSC_VIEWER_STDOUT_SELF);
>
>   ierr = MatMatSolve(D,I,Dinv); CHKERRQ(ierr);
>
>   MatView(Dinv,PETSC_VIEWER_STDOUT_SELF);
>
> // Results in such output:
>
> //==== MATRIX D
> row 0: (0, 2)
> row 1: (1, 1)
> row 2: (2, 1)
> row 3: (3, 1)
>
> //==== MATRIX D after LU decomposition
> row 0: (0, 0.5)
> row 1: (1, 1)
> row 2: (2, 1)
> row 3: (3, 1)
>
> //==== Identity matrix
> row 0: (0, 1)
> row 1: (1, 1)
> row 2: (2, 1)
> row 3: (3, 1)
>
> //==== inverse of D  (result of MatMatSolve(D,I,Dinv))
> row 0:
> row 1:
> row 2:
> row 3:
>
>
> The approach with :
>
>     ierr = MatGetOrdering(D,  MATORDERING_NATURAL,  &perm,  &iperm);
> CHKERRQ(ierr);
>     ierr = MatGetFactor(D, MAT_SOLVER_PETSC, MAT_FACTOR_LU,&F);
>     ierr = MatFactorInfoInitialize(&info); CHKERRQ(ierr);
>     ierr = MatLUFactorSymbolic(F,D,perm,iperm,&info);CHKERRQ(ierr);
>     ierr = MatLUFactorNumeric(F,D,&info);CHKERRQ(ierr);
>     ierr = MatLUFactor(D, perm, iperm, &info); CHKERRQ(ierr);
>
> Provides same output
>
> Thank You for Your help.
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20110330/a01eaaf3/attachment.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: ex1.c
Type: text/x-csrc
Size: 6160 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20110330/a01eaaf3/attachment.c>