1 static char help[] = "Test MatProductReplaceMats() \n\
2 Modified from the code contributed by Pierre Jolivet \n\n";
3
4 #include <petscmat.h>
5
main(int argc,char ** args)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;
11
12 PetscFunctionBeginUser;
13 PetscCall(PetscInitialize(&argc, &args, NULL, help));
14
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));
24
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));
28
29 /* compute B = A^T * A or B = A^H * A */
30 PetscCall(MatProductCreate(A, A, NULL, &B));
31
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));
51
52 PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &Conjugate));
53 if (conjugate) PetscCall(MatConjugate(Conjugate));
54
55 /* replace input A by Conjugate */
56 PetscCall(MatProductReplaceMats(Conjugate, NULL, NULL, B));
57
58 PetscCall(MatProductNumeric(B));
59 PetscCall(MatViewFromOptions(B, NULL, "-product_view"));
60
61 PetscCall(MatDestroy(&Conjugate));
62 if (!convert) {
63 Bdense = B;
64 B = NULL;
65 }
66 }
67
68 /* Compare Bdense and B */
69 PetscCall(MatMultEqual(Bdense, B, 10, &equal));
70 PetscCheck(equal, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Bdense != B");
71
72 PetscCall(MatDestroy(&Bdense));
73 PetscCall(MatDestroy(&B));
74 PetscCall(MatDestroy(&A));
75 PetscCall(PetscFinalize());
76 return 0;
77 }
78
79 /*TEST
80
81 test:
82 suffix: 1
83 args: -conjugate false -atb
84 output_file: output/empty.out
85
86 test:
87 suffix: 2
88 args: -conjugate true -atb
89 output_file: output/empty.out
90
91 test:
92 suffix: 3
93 args: -conjugate false
94 output_file: output/empty.out
95
96 test:
97 suffix: 4
98 args: -ptap
99 output_file: output/empty.out
100
101 test:
102 suffix: 5
103 args: -abt
104 output_file: output/empty.out
105
106 test:
107 suffix: 6
108 nsize: 2
109 args: -conjugate false -atb
110 output_file: output/empty.out
111
112 test:
113 suffix: 7
114 nsize: 2
115 args: -conjugate true -atb
116 output_file: output/empty.out
117
118 test:
119 suffix: 8
120 nsize: 2
121 args: -conjugate false
122 output_file: output/empty.out
123
124 test:
125 suffix: 9
126 nsize: 2
127 args: -ptap
128 output_file: output/empty.out
129
130 test:
131 suffix: 10
132 nsize: 2
133 args: -abt
134 output_file: output/empty.out
135
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
141
142 TEST*/
143