xref: /petsc/src/mat/tests/ex101.c (revision b122ec5aa1bd4469eb4e0673542fb7de3f411254)
1c4762a1bSJed Brown static char help[] = "Testing PtAP for SeqMAIJ matrix, P, with SeqAIJ matrix, A.\n\n";
2c4762a1bSJed Brown 
3c4762a1bSJed Brown #include <petscmat.h>
4c4762a1bSJed Brown 
5c4762a1bSJed Brown int main(int argc,char **argv)
6c4762a1bSJed Brown {
7c4762a1bSJed Brown   Mat            pA,P,aijP;
8c4762a1bSJed Brown   PetscScalar    pa[]={1.,-1.,0.,0.,1.,-1.,0.,0.,1.};
9c4762a1bSJed Brown   PetscInt       i,pij[]={0,1,2};
10c4762a1bSJed Brown   PetscInt       aij[3][3]={{0,1,2},{3,4,5},{6,7,8}};
11c4762a1bSJed Brown   Mat            A,mC,C;
12c4762a1bSJed Brown   PetscScalar    one=1.;
13c4762a1bSJed Brown   PetscMPIInt    size;
14c20d7725SJed Brown   PetscBool      flg;
15c4762a1bSJed Brown 
16*b122ec5aSJacob Faibussowitsch   CHKERRQ(PetscInitialize(&argc,&argv,(char*)0,help));
175f80ce2aSJacob Faibussowitsch   CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size));
182c71b3e2SJacob Faibussowitsch   PetscCheckFalse(size != 1,PETSC_COMM_WORLD,PETSC_ERR_SUP,"This is a uniprocessor example only!");
19c4762a1bSJed Brown 
20c4762a1bSJed Brown   /* Create MAIJ matrix, P */
215f80ce2aSJacob Faibussowitsch   CHKERRQ(MatCreate(PETSC_COMM_SELF,&pA));
225f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetSizes(pA,3,3,3,3));
235f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetType(pA,MATSEQAIJ));
245f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetUp(pA));
255f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetOption(pA,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE));
265f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetValues(pA,3,pij,3,pij,pa,ADD_VALUES));
275f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAssemblyBegin(pA,MAT_FINAL_ASSEMBLY));
285f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAssemblyEnd(pA,MAT_FINAL_ASSEMBLY));
295f80ce2aSJacob Faibussowitsch   CHKERRQ(MatCreateMAIJ(pA,3,&P));
305f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&pA));
31c4762a1bSJed Brown 
32c4762a1bSJed Brown   /* Create AIJ equivalent matrix, aijP, for comparison testing */
335f80ce2aSJacob Faibussowitsch   CHKERRQ(MatConvert(P,MATSEQAIJ,MAT_INITIAL_MATRIX,&aijP));
34c4762a1bSJed Brown 
35c20d7725SJed Brown   /* Create AIJ matrix A */
365f80ce2aSJacob Faibussowitsch   CHKERRQ(MatCreate(PETSC_COMM_SELF,&A));
375f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetSizes(A,9,9,9,9));
385f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetType(A,MATSEQAIJ));
395f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetUp(A));
405f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetOption(A,MAT_IGNORE_ZERO_ENTRIES,PETSC_TRUE));
415f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetValues(A,3,aij[0],3,aij[0],pa,ADD_VALUES));
425f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetValues(A,3,aij[1],3,aij[1],pa,ADD_VALUES));
435f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetValues(A,3,aij[2],3,aij[2],pa,ADD_VALUES));
44c4762a1bSJed Brown   for (i=0; i<9; i++) {
455f80ce2aSJacob Faibussowitsch     CHKERRQ(MatSetValue(A,i,i,one,ADD_VALUES));
46c4762a1bSJed Brown   }
475f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
485f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
49c4762a1bSJed Brown 
50c4762a1bSJed Brown   /* Perform PtAP_SeqAIJ_SeqAIJ for comparison testing */
515f80ce2aSJacob Faibussowitsch   CHKERRQ(MatPtAP(A,aijP,MAT_INITIAL_MATRIX,1.,&C));
52c20d7725SJed Brown 
53c20d7725SJed Brown   /* Perform PtAP_SeqAIJ_SeqMAIJ */
54c20d7725SJed Brown   /* Developer API */
555f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductCreate(A,P,NULL,&mC));
565f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductSetType(mC,MATPRODUCT_PtAP));
575f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductSetAlgorithm(mC,"default"));
585f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductSetFill(mC,PETSC_DEFAULT));
595f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductSetFromOptions(mC));
605f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductSymbolic(mC));
615f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductNumeric(mC));
625f80ce2aSJacob Faibussowitsch   CHKERRQ(MatProductNumeric(mC));
63c4762a1bSJed Brown 
64c4762a1bSJed Brown   /* Check mC = C */
655f80ce2aSJacob Faibussowitsch   CHKERRQ(MatEqual(C,mC,&flg));
6628b400f6SJacob Faibussowitsch   PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"MatProduct C != mC");
675f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&mC));
68c20d7725SJed Brown 
69c20d7725SJed Brown   /* User API */
705f80ce2aSJacob Faibussowitsch   CHKERRQ(MatPtAP(A,P,MAT_INITIAL_MATRIX,1.,&mC));
715f80ce2aSJacob Faibussowitsch   CHKERRQ(MatPtAP(A,P,MAT_REUSE_MATRIX,1.,&mC));
72c20d7725SJed Brown 
73c20d7725SJed Brown   /* Check mC = C */
745f80ce2aSJacob Faibussowitsch   CHKERRQ(MatEqual(C,mC,&flg));
7528b400f6SJacob Faibussowitsch   PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"MatPtAP C != mC");
765f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&mC));
77c4762a1bSJed Brown 
78c4762a1bSJed Brown   /* Cleanup */
795f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&P));
805f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&aijP));
815f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&A));
825f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&C));
83*b122ec5aSJacob Faibussowitsch   CHKERRQ(PetscFinalize());
84*b122ec5aSJacob Faibussowitsch   return 0;
85c4762a1bSJed Brown }
86c4762a1bSJed Brown 
87c4762a1bSJed Brown /*TEST
88c4762a1bSJed Brown 
89c4762a1bSJed Brown    test:
90c4762a1bSJed Brown       output_file: output/ex101.out
91c4762a1bSJed Brown 
92c4762a1bSJed Brown TEST*/
93