Actual source code: ex122.c
petsc-3.14.6 2021-03-30
1: static char help[] = "Test MatMatMult() for AIJ and Dense matrices.\n\n";
3: #include <petscmat.h>
5: int main(int argc,char **argv)
6: {
7: Mat A,B,C;
8: PetscInt M=10,N=5;
10: PetscRandom r;
11: PetscBool equal=PETSC_FALSE;
12: PetscReal fill = 1.0;
13: PetscInt nza,am,an;
15: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
16: PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);
17: PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);
18: PetscOptionsGetReal(NULL,NULL,"-fill",&fill,NULL);
20: PetscRandomCreate(PETSC_COMM_WORLD,&r);
21: PetscRandomSetFromOptions(r);
23: /* create a aij matrix A */
24: MatCreate(PETSC_COMM_WORLD,&A);
25: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,M);
26: MatSetType(A,MATAIJ);
27: nza = (PetscInt)(.3*M); /* num of nozeros in each row of A */
28: MatSeqAIJSetPreallocation(A,nza,NULL);
29: MatMPIAIJSetPreallocation(A,nza,NULL,nza,NULL);
30: MatSetRandom(A,r);
32: /* create a dense matrix B */
33: MatGetLocalSize(A,&am,&an);
34: MatCreate(PETSC_COMM_WORLD,&B);
35: MatSetSizes(B,PETSC_DECIDE,am,N,PETSC_DECIDE);
36: MatSetType(B,MATDENSE);
37: MatSeqDenseSetPreallocation(B,NULL);
38: MatMPIDenseSetPreallocation(B,NULL);
39: MatSetRandom(B,r);
40: PetscRandomDestroy(&r);
41: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
42: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
45: /* Test MatMatMult() */
46: MatMatMult(B,A,MAT_INITIAL_MATRIX,fill,&C);
47: MatMatMult(B,A,MAT_REUSE_MATRIX,fill,&C);
48: MatMatMultEqual(B,A,C,10,&equal);
49: if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"C != B*A");
51: MatDestroy(&C);
52: MatDestroy(&B);
53: MatDestroy(&A);
54: PetscFinalize();
55: return ierr;
56: }
59: /*TEST
61: test:
62: output_file: output/ex122.out
64: TEST*/