1 static char help[] = "Tests MatComputeOperator() and MatComputeOperatorTranspose()\n\n";
2
3 #include <petscmat.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 Mat A, Ae, Aet;
8 char filename[PETSC_MAX_PATH_LEN];
9 char expltype[128], *etype = NULL;
10 PetscInt bs = 1;
11 PetscBool flg, check = PETSC_TRUE;
12
13 PetscFunctionBeginUser;
14 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
15
16 PetscCall(PetscOptionsGetString(NULL, NULL, "-expl_type", expltype, sizeof(expltype), &flg));
17 if (flg) PetscCall(PetscStrallocpy(expltype, &etype));
18 PetscCall(PetscOptionsGetString(NULL, NULL, "-f", filename, sizeof(filename), &flg));
19 PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL));
20 if (!flg) {
21 PetscInt M = 13, N = 6;
22 PetscBool sbaij;
23
24 PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL));
25 PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL));
26 PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, N, NULL, &A));
27 PetscCall(MatSetBlockSize(A, bs));
28 PetscCall(MatSetRandom(A, NULL));
29 PetscCall(PetscStrcmp(etype, MATSBAIJ, &sbaij));
30 if (sbaij) {
31 Mat At;
32 PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &At));
33 PetscCall(MatAXPY(A, 1.0, At, DIFFERENT_NONZERO_PATTERN));
34 PetscCall(MatDestroy(&At));
35 }
36 } else {
37 PetscViewer viewer;
38
39 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, filename, FILE_MODE_READ, &viewer));
40 PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
41 PetscCall(MatSetBlockSize(A, bs));
42 PetscCall(MatSetFromOptions(A));
43 PetscCall(MatLoad(A, viewer));
44 PetscCall(PetscViewerDestroy(&viewer));
45 }
46 PetscCall(PetscObjectSetName((PetscObject)A, "Matrix"));
47 PetscCall(MatViewFromOptions(A, NULL, "-view_expl"));
48
49 PetscCall(MatComputeOperator(A, etype, &Ae));
50 PetscCall(PetscObjectSetName((PetscObject)Ae, "Explicit matrix"));
51 PetscCall(MatViewFromOptions(Ae, NULL, "-view_expl"));
52
53 PetscCall(PetscOptionsGetBool(NULL, NULL, "-check", &check, NULL));
54 if (check) {
55 Mat A2;
56 PetscReal err, tol = PETSC_SMALL;
57
58 PetscCall(PetscOptionsGetReal(NULL, NULL, "-tol", &tol, NULL));
59 PetscCall(MatConvert(A, etype, MAT_INITIAL_MATRIX, &A2));
60 PetscCall(MatAXPY(A2, -1.0, Ae, DIFFERENT_NONZERO_PATTERN));
61 PetscCall(MatNorm(A2, NORM_FROBENIUS, &err));
62 if (err > tol) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error %g > %g (type %s)\n", (double)err, (double)tol, etype));
63 PetscCall(MatDestroy(&A2));
64 }
65
66 PetscCall(MatComputeOperatorTranspose(A, etype, &Aet));
67 PetscCall(PetscObjectSetName((PetscObject)Aet, "Explicit matrix transpose"));
68 PetscCall(MatViewFromOptions(Aet, NULL, "-view_expl"));
69
70 PetscCall(PetscFree(etype));
71 PetscCall(MatDestroy(&Ae));
72 PetscCall(MatDestroy(&Aet));
73 PetscCall(MatDestroy(&A));
74 PetscCall(PetscFinalize());
75 return 0;
76 }
77
78 /*TEST
79
80 test:
81 output_file: output/empty.out
82
83 testset:
84 suffix: matexpl_rect
85 output_file: output/empty.out
86 nsize: {{1 3}}
87 args: -expl_type {{dense aij baij}}
88
89 testset:
90 suffix: matexpl_square
91 output_file: output/empty.out
92 nsize: {{1 3}}
93 args: -bs {{1 2 3}} -M 36 -N 36 -expl_type {{dense aij baij sbaij}}
94
95 TEST*/
96