[petsc-users] Determinant of the matrix

Jelena Slivka slivkaje at gmail.com
Mon Feb 25 18:47:43 CST 2013


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? Also, I need to use Cholesky
factorization, is there a way to avoid MUMPS so I could access the diagonal?


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
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130225/5c3d87dc/attachment.html>


More information about the petsc-users mailing list