[petsc-users] Determinant of the matrix

[Barry Smith]

  Jelena,

     MUMPS has its own (complicated) data structure for storing the factored matrix, you cannot easily look inside it. 

      However, fortunately MUMPS has a way of accessing the determinate of the matrix directly. Look at src/mat/impls/aij/mpi/mumps.c and search for the word determinant, you will find the command line argument to have it computed and where it will be printed out if you call MatView() on the factored matrix. You mush pass the PetscViewerFormat PETSC_VIEWER_ASCII_INFO to PetscViewerSetFormat() before viewing the matrix.

   Good luck,

   Barry

On Feb 25, 2013, at 6:02 PM, Jelena Slivka <slivkaje at gmail.com> wrote:

> Hello,
> I am trying to find the determinant of the matrix. I have found this thread:
> 
> http://lists.mcs.anl.gov/pipermail/petsc-users/2010-July/006716.html
> 
> But, could you please explain how to access the data structure directly?
> So far I have tried the following:
> 
> #include <petsc-private/matimpl.h>
> #include <../src/mat/impls/aij/seq/aij.h>
> ...
> // (Solving the system A*x = rhs
> MatSetOption(A, MAT_SPD, PETSC_TRUE);
> Mat F;
> MatGetFactor(A, MATSOLVERMUMPS, MAT_FACTOR_CHOLESKY, &F);
> MatFactorInfo factinfo;
> MatFactorInfoInitialize(&factinfo);
> IS perm, iscol;
> MatGetOrdering(A, MATORDERINGNATURAL, &perm, &iscol);
> MatCholeskyFactorSymbolic(F, A, perm, &factinfo);
> MatCholeskyFactorNumeric(F, A, &factinfo);        
> ISDestroy(&iscol);
> ISDestroy(&perm);
>         
> MatSolve(F, rhs, x);    
> 
> // Now I would like to access the diagonal part of the matrix F, can I do something like:
> 
> Mat_SeqAIJ *aa = (Mat_SeqAIJ*)F->data;
> 
> Do you maybe have an example I could look at?
> 
> Grateful in advance,
> Jelena