[petsc-users] Matrix inverse and LU decomposition

[Bartłomiej W]

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/e331c338/attachment.htm>