Actual source code: ex34.c
petsc-3.14.6 2021-03-30
1: static char help[] = "Test MatMatMult() and MatTransposeMatMult() for MPIAIJ and MPIDENSE matrices. \n\
2: Sequential part of mpidense matrix allows changes made by MatDenseSetLDA(). \n\n";
4: #include <petsc.h>
6: int main(int argc, char ** argv)
7: {
8: Mat A, B, C, C1;
9: PetscMPIInt size;
11: PetscInt i,ia[2] = { 0, 2 }, ja[2] = { 0, 1 }, lda = 4;
12: PetscScalar a[2] = { 1.0, 1.0 }, *data;
13: PetscBool flg;
15: PetscInitialize(&argc, &argv, (char*)0,help);if (ierr) return ierr;
16: MPI_Comm_size(PETSC_COMM_WORLD, &size);
17: if (size != 2) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Must use 2 processors");
19: MatCreate(PETSC_COMM_WORLD, &A);
20: MatSetType(A, MATMPIAIJ);
21: MatSetSizes(A, 1, 1, 2, 2);
22: MatMPIAIJSetPreallocationCSR(A, ia, ja, a);
24: PetscCalloc1(4 * lda,&data);
25: for (i = 0; i < 4; ++i) data[lda * i] = i * 1.0;
27: MatCreateDense(PETSC_COMM_WORLD, 1, PETSC_DECIDE, 2, 4, data, &B);
28: MatSetOptionsPrefix(B,"b_");
29: MatSetFromOptions(B);
30: MatDenseSetLDA(B, lda);
32: /* Test MatMatMult() */
33: MatMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C);
34: MatMatMult(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C);
36: MatMatMultEqual(A,B,C,10,&flg);
37: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in MatMatMult() for C");
39: /* Test user-provided mpidense matrix product */
40: MatDuplicate(C,MAT_COPY_VALUES,&C1);
41: MatMatMult(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C1);
42: MatMatMultEqual(A,B,C1,10,&flg);
43: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in MatMatMult() for C1");
45: MatDestroy(&C1);
46: MatDestroy(&C);
48: /* Test MatTransposeMatMult() */
49: MatTransposeMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C);
50: MatTransposeMatMult(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C);
52: MatTransposeMatMultEqual(A,B,C,10,&flg);
53: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in MatTransposeMatMult()");
54: MatDestroy(&C);
56: MatDestroy(&B);
57: MatDestroy(&A);
58: PetscFree(data);
60: PetscFinalize();
61: return ierr;
62: }
64: /*TEST
66: test:
67: suffix: 1
68: nsize: 2
69: output_file: output/ex34.out
71: test:
72: suffix: 1_cuda
73: requires: cuda
74: nsize: 2
75: args: -b_mat_type mpidensecuda
76: output_file: output/ex34.out
77: TEST*/