Actual source code: ex258.c
1: static char help[] = "Test MatProductReplaceMats() \n\
2: Modified from the code contributed by Pierre Jolivet \n\n";
4: #include <petscmat.h>
6: int main(int argc, char **args)
7: {
8: PetscInt n = 2, convert;
9: Mat A, B, Bdense, Conjugate;
10: PetscBool conjugate = PETSC_FALSE, equal, flg;
12: PetscFunctionBeginUser;
13: PetscCall(PetscInitialize(&argc, &args, NULL, help));
15: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
16: PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n));
17: PetscCall(MatSetType(A, MATDENSE));
18: PetscCall(MatSetFromOptions(A));
19: PetscCall(MatSeqDenseSetPreallocation(A, NULL));
20: PetscCall(MatMPIDenseSetPreallocation(A, NULL));
21: PetscCall(MatSetRandom(A, NULL));
22: PetscCall(MatViewFromOptions(A, NULL, "-A_view"));
23: PetscCall(PetscOptionsGetBool(NULL, NULL, "-conjugate", &conjugate, NULL));
25: for (convert = 0; convert < 2; convert++) {
26: /* convert dense matrix A to aij format */
27: if (convert) PetscCall(MatConvert(A, MATAIJ, MAT_INPLACE_MATRIX, &A));
29: /* compute B = A^T * A or B = A^H * A */
30: PetscCall(MatProductCreate(A, A, NULL, &B));
32: flg = PETSC_FALSE;
33: PetscCall(PetscOptionsGetBool(NULL, NULL, "-atb", &flg, NULL));
34: if (flg) {
35: PetscCall(MatProductSetType(B, MATPRODUCT_AtB));
36: } else {
37: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ptap", &flg, NULL));
38: if (flg) {
39: PetscCall(MatProductSetType(B, MATPRODUCT_PtAP));
40: } else {
41: PetscCall(PetscOptionsGetBool(NULL, NULL, "-abt", &flg, NULL));
42: if (flg) {
43: PetscCall(MatProductSetType(B, MATPRODUCT_ABt));
44: } else {
45: PetscCall(MatProductSetType(B, MATPRODUCT_AB));
46: }
47: }
48: }
49: PetscCall(MatProductSetFromOptions(B));
50: PetscCall(MatProductSymbolic(B));
52: PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &Conjugate));
53: if (conjugate) PetscCall(MatConjugate(Conjugate));
55: /* replace input A by Conjugate */
56: PetscCall(MatProductReplaceMats(Conjugate, NULL, NULL, B));
58: PetscCall(MatProductNumeric(B));
59: PetscCall(MatViewFromOptions(B, NULL, "-product_view"));
61: PetscCall(MatDestroy(&Conjugate));
62: if (!convert) {
63: Bdense = B;
64: B = NULL;
65: }
66: }
68: /* Compare Bdense and B */
69: PetscCall(MatMultEqual(Bdense, B, 10, &equal));
70: PetscCheck(equal, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Bdense != B");
72: PetscCall(MatDestroy(&Bdense));
73: PetscCall(MatDestroy(&B));
74: PetscCall(MatDestroy(&A));
75: PetscCall(PetscFinalize());
76: return 0;
77: }
79: /*TEST
81: test:
82: suffix: 1
83: args: -conjugate false -atb
84: output_file: output/ex258_1.out
86: test:
87: suffix: 2
88: args: -conjugate true -atb
89: output_file: output/ex258_1.out
91: test:
92: suffix: 3
93: args: -conjugate false
94: output_file: output/ex258_1.out
96: test:
97: suffix: 4
98: args: -ptap
99: output_file: output/ex258_1.out
101: test:
102: suffix: 5
103: args: -abt
104: output_file: output/ex258_1.out
106: test:
107: suffix: 6
108: nsize: 2
109: args: -conjugate false -atb
110: output_file: output/ex258_1.out
112: test:
113: suffix: 7
114: nsize: 2
115: args: -conjugate true -atb
116: output_file: output/ex258_1.out
118: test:
119: suffix: 8
120: nsize: 2
121: args: -conjugate false
122: output_file: output/ex258_1.out
124: test:
125: suffix: 9
126: nsize: 2
127: args: -ptap
128: output_file: output/ex258_1.out
130: test:
131: suffix: 10
132: nsize: 2
133: args: -abt
134: output_file: output/ex258_1.out
136: test:
137: suffix: 11
138: nsize: 2
139: args: -conjugate true -atb -mat_product_algorithm backend
140: TODO: bug: MatProductReplaceMats() does not change the product for this test
142: TEST*/