Actual source code: ex101.c

petsc-3.3-p6 2013-02-11
  1: static char help[] = "Testing PtAP for SeqMAIJ matrix, P, with SeqAIJ matrix, A.\n\n";

  3: #include <petscmat.h>

  7: int main(int argc,char **argv) {
  8:   Mat            pA,P,aijP;
  9:   PetscScalar    pa[]={1.,-1.,0.,0.,1.,-1.,0.,0.,1.};
 10:   PetscInt       pij[]={0,1,2};
 11:   PetscInt       aij[3][3]={{0,1,2},{3,4,5},{6,7,8}};
 12:   Mat            A,mC,C;
 13:   PetscScalar    one=1.;

 16:   PetscInitialize(&argc,&argv,(char *)0,help);

 18:   /* Create MAIJ matrix, P */
 19:   MatCreate(PETSC_COMM_SELF,&pA);
 20:   MatSetSizes(pA,3,3,3,3);
 21:   MatSetType(pA,MATSEQAIJ);
 22:   MatSetUp(pA);
 23:   MatSetOption(pA,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE);
 24:   MatSetValues(pA,3,pij,3,pij,pa,ADD_VALUES);
 25:   MatAssemblyBegin(pA,MAT_FINAL_ASSEMBLY);
 26:   MatAssemblyEnd(pA,MAT_FINAL_ASSEMBLY);
 27:   MatCreateMAIJ(pA,3,&P);
 28:   MatDestroy(&pA);

 30:   /* Create AIJ equivalent matrix, aijP, for comparison testing */
 31:   MatConvert(P,MATSEQAIJ,MAT_INITIAL_MATRIX,&aijP);

 33:   /* Create AIJ matrix, A */
 34:   MatCreate(PETSC_COMM_SELF,&A);
 35:   MatSetSizes(A,9,9,9,9);
 36:   MatSetType(A,MATSEQAIJ);
 37:   MatSetUp(A);
 38:   MatSetOption(A,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE);
 39:   MatSetValues(A,3,aij[0],3,aij[0],pa,ADD_VALUES);
 40:   MatSetValues(A,3,aij[1],3,aij[1],pa,ADD_VALUES);
 41:   MatSetValues(A,3,aij[2],3,aij[2],pa,ADD_VALUES);
 42:   {int i;
 43:    for (i=0;i<9;i++) {
 44:     MatSetValue(A,i,i,one,ADD_VALUES);
 45:    }
 46:   }
 47:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 48:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 49: 
 50:   /* Perform SeqAIJ_SeqMAIJ PtAP */
 51:   MatPtAP(A,P,MAT_INITIAL_MATRIX,1.,&mC);
 52:   MatView(mC,PETSC_VIEWER_STDOUT_SELF);

 54:   /* Perform SeqAIJ_SeqAIJ PtAP for comparison testing */
 55:   MatPtAP(A,aijP,MAT_INITIAL_MATRIX,1.,&C);
 56:   MatView(C,PETSC_VIEWER_STDOUT_SELF);

 58:   /* Perform diff of two matrices */
 59:   MatAXPY(C,-1.0,mC,DIFFERENT_NONZERO_PATTERN);
 60:   /* Note: We should be able to use SAME_NONZERO_PATTERN on the line above, */
 61:   /*       but don't because this flag doesn't assist testing. */
 62:   MatView(C,PETSC_VIEWER_STDOUT_SELF);

 64:   /* Cleanup */
 65:   MatDestroy(&P);
 66:   MatDestroy(&aijP);
 67:   MatDestroy(&A);
 68:   MatDestroy(&C);
 69:   MatDestroy(&mC);
 70:   PetscFinalize();
 71:   return(0);
 72: }