[petsc-users] Determinant of the matrix

[Barry Smith]

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

> Thank you very much, but I need to access this value in my program and not just print it, is there a way to do this?

   You need to look at the code that prints out the determinant and call the same routine from your code directly.

>  Also, I need to use Cholesky factorization, is there a way to avoid MUMPS so I could access the diagonal?

   If the problem is small enough to run sequentially then you can just the MPIBAIJ format Cholesky and access the diagonal directly.

   Barry

> 
> 
> On Mon, Feb 25, 2013 at 7:13 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
> 
>   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
> 
>