1 static char help[] = "Tests MatSeqAIJKron.\n\n";
2
3 #include <petscmat.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 Mat A, B, C, K, Ad, Bd;
8 const PetscScalar *Bv;
9 PetscInt n = 10, m = 20, p = 7, q = 17;
10 PetscBool flg;
11
12 PetscFunctionBeginUser;
13 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
14 PetscCall(MatCreateDense(PETSC_COMM_SELF, m, n, m, n, NULL, &Ad));
15 PetscCall(MatCreateDense(PETSC_COMM_SELF, p, q, p, q, NULL, &Bd));
16 PetscCall(MatSetRandom(Ad, NULL));
17 PetscCall(MatSetRandom(Bd, NULL));
18 PetscCall(MatFilter(Ad, 0.2, PETSC_FALSE, PETSC_FALSE));
19 PetscCall(MatFilter(Bd, 0.2, PETSC_FALSE, PETSC_FALSE));
20 PetscCall(MatConvert(Ad, MATAIJ, MAT_INITIAL_MATRIX, &A));
21 PetscCall(MatConvert(Bd, MATAIJ, MAT_INITIAL_MATRIX, &B));
22 PetscCall(MatSeqAIJKron(A, B, MAT_INITIAL_MATRIX, &C));
23 PetscCall(MatViewFromOptions(A, NULL, "-A_view"));
24 PetscCall(MatViewFromOptions(B, NULL, "-B_view"));
25 PetscCall(MatViewFromOptions(C, NULL, "-C_view"));
26 PetscCall(MatDenseGetArrayRead(Bd, &Bv));
27 PetscCall(MatCreateKAIJ(A, p, q, NULL, Bv, &K));
28 PetscCall(MatDenseRestoreArrayRead(Bd, &Bv));
29 PetscCall(MatMultEqual(C, K, 10, &flg));
30 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "K*x != C*x");
31 PetscCall(MatScale(A, 1.3));
32 PetscCall(MatScale(B, 0.3));
33 PetscCall(MatScale(Bd, 0.3));
34 PetscCall(MatSeqAIJKron(A, B, MAT_REUSE_MATRIX, &C));
35 PetscCall(MatDenseGetArrayRead(Bd, &Bv));
36 PetscCall(MatKAIJSetT(K, p, q, Bv));
37 PetscCall(MatDenseRestoreArrayRead(Bd, &Bv));
38 PetscCall(MatMultEqual(C, K, 10, &flg));
39 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "K*x != C*x");
40 PetscCall(MatDestroy(&K));
41 PetscCall(MatDestroy(&C));
42 PetscCall(MatDestroy(&B));
43 PetscCall(MatDestroy(&A));
44 PetscCall(MatDestroy(&Bd));
45 PetscCall(MatDestroy(&Ad));
46 PetscCall(PetscFinalize());
47 return 0;
48 }
49
50 /*TEST
51
52 test:
53 suffix: 1
54 nsize: 1
55 output_file: output/empty.out
56
57 TEST*/
58