Cholesky factorization problem

[JUNWANG at uwm.edu]

Hi, Hong:

    Really appreciate your help!
    The MatSolve() works fine. It provides the right solution. How about if i
wanna see the elements in the factorized up-triangle matix R, where A =R'R.
Since the MatView() cannot show the logical matrix elements, which function
should i call?
    Thank you !

Jun

Quoting Hong Zhang <hzhang at mcs.anl.gov>:

>
> Jun,
>
> Your code looks fine. The factored matrix
> "result" is an internal matrix to be used by MatSolve,
> thus, its actural storage doesn't match the logical
> matrix elements. MatView() only display off-diagonal blocks
> of the factored matrix (these values may not match the logical
> factored matrix values!).
>
> We'll add a note to MatView() when user display
> an internal factored matrix.
>
> If you get wrong MatSolve() with your code,
> let us know.
>
> Thanks for reporting if,
>
> Hong
>
>
> On Fri, 13 Apr 2007 junwang at uwm.edu wrote:
>
> >
> >
> > Hello:
> >    I have written a test sample program for using the Cholesky
> Factorization
> > function, here is the code i wrote, it compiled successfully. But the
> result is
> > not correct, is anyone who can help me out? Thank you!
> > ***************************************************************************
> >
> > #include "petsc.h"
> > #include "petscmat.h"
> > #include <iostream.h>
> > #include <stdlib.h>
> >
> > int main(int argc,char **args)
> > {
> >   PetscErrorCode ierr;
> >   Mat            mat;
> >   int ind1[2],ind2[2];
> >   PetscScalar temp[2*2];
> >   PetscInt *nnz=new PetscInt[3];
> >
> >   PetscInitialize(&argc,&args,0,0);
> >   nnz[0]=2;nnz[1]=1;nnz[1]=1;
> >
> >   ierr = MatCreateSeqSBAIJ(PETSC_COMM_SELF,2,6,6,0,nnz,&mat);CHKERRQ(ierr);
> >   ind1[0]=0;ind1[1]=1;
> >   temp[0]=3;temp[1]=2;temp[2]=0;temp[3]=3;
> >   ierr = MatSetValues(mat,2,ind1,2,ind1,temp,INSERT_VALUES);CHKERRQ(ierr);
> >   ind2[0]=4;ind2[1]=5;
> >   temp[0]=1;temp[1]=1;temp[2]=2;temp[3]=1;
> >   ierr = MatSetValues(mat,2,ind1,2,ind2,temp,INSERT_VALUES);CHKERRQ(ierr);
> >   ind1[0]=2;ind1[1]=3;
> >   temp[0]=4;temp[1]=1;temp[2]=1;temp[3]=5;
> >   ierr = MatSetValues(mat,2,ind1,2,ind1,temp,INSERT_VALUES);CHKERRQ(ierr);
> >   ind1[0]=4;ind1[1]=5;
> >   temp[0]=5;temp[1]=1;temp[2]=1;temp[3]=6;
> >   ierr = MatSetValues(mat,2,ind1,2,ind1,temp,INSERT_VALUES);CHKERRQ(ierr);
> >
> >    ierr = MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
> >    ierr = MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
> >
> >    MatView(mat,PETSC_VIEWER_STDOUT_SELF);
> >    // begin cholesky factorization
> >    IS perm,colp;
> >
> >    MatGetOrdering(mat,MATORDERING_NATURAL,&perm,&colp);
> >    ierr = ISDestroy(colp);CHKERRQ(ierr);
> >
> >    MatFactorInfo info;
> >    info.fill=1.0;
> >    Mat result;
> >
> >    ierr = MatCholeskyFactorSymbolic(mat,perm,&info,&result); CHKERRQ(ierr);
> >    ierr = MatCholeskyFactorNumeric(mat,&info,&result);CHKERRQ(ierr);
> >     MatView(result, PETSC_VIEWER_STDOUT_SELF);
> >
> >    ierr = ISDestroy(perm);CHKERRQ(ierr);
> >    ierr= MatDestroy(mat);CHKERRQ(ierr);
> >    PetscFinalize();
> >    return 0;
> > }
> > ***************************************************************************
> >
> > the result given below:
> > mat:
> > row0: (0,3),(1,2),(4,1),(5,1)
> > row1: (0,2),(1,3),(4,2),(5,1)
> > row2: (2,4),(3,1)
> > row3: (2,1),(3,5)
> > row4: (4,5),(5,1)
> > row5: (4,1),(5,6)
> > result:
> > row0: (4,0.2),(5,-0.2)
> > row1: (4,-0.8) (5, -0.2)
> > row2:
> > row3:
> > row4:
> > row5:
> > ***************************************************************************
> >
> >    The matrix is in SBAIJ format, block size is 2, which cholesky
> factorization
> > function should i call? really need that for the project! Thank you!
> >
> > Jun
> >
> >
>
>