Actual source code: ex49.c
2: static char help[] = "Tests MatTranspose(), MatNorm(), and MatAXPY().\n\n";
4: #include <petscmat.h>
6: int main(int argc, char **argv)
7: {
8: Mat mat, tmat = 0;
9: PetscInt m = 4, n, i, j;
10: PetscMPIInt size, rank;
11: PetscInt rstart, rend, rect = 0;
12: PetscBool flg;
13: PetscScalar v;
14: PetscReal normf, normi, norm1;
15: MatInfo info;
17: PetscFunctionBeginUser;
18: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
19: PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
20: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
21: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
22: n = m;
23: PetscCall(PetscOptionsHasName(NULL, NULL, "-rect1", &flg));
24: if (flg) {
25: n += 2;
26: rect = 1;
27: }
28: PetscCall(PetscOptionsHasName(NULL, NULL, "-rect2", &flg));
29: if (flg) {
30: n -= 2;
31: rect = 1;
32: }
34: /* Create and assemble matrix */
35: PetscCall(MatCreate(PETSC_COMM_WORLD, &mat));
36: PetscCall(MatSetSizes(mat, PETSC_DECIDE, PETSC_DECIDE, m, n));
37: PetscCall(MatSetFromOptions(mat));
38: PetscCall(MatSetUp(mat));
39: PetscCall(MatGetOwnershipRange(mat, &rstart, &rend));
40: for (i = rstart; i < rend; i++) {
41: for (j = 0; j < n; j++) {
42: v = 10 * i + j;
43: PetscCall(MatSetValues(mat, 1, &i, 1, &j, &v, INSERT_VALUES));
44: }
45: }
46: PetscCall(MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY));
47: PetscCall(MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY));
49: /* Print info about original matrix */
50: PetscCall(MatGetInfo(mat, MAT_GLOBAL_SUM, &info));
51: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "original matrix nonzeros = %" PetscInt_FMT ", allocated nonzeros = %" PetscInt_FMT "\n", (PetscInt)info.nz_used, (PetscInt)info.nz_allocated));
52: PetscCall(MatNorm(mat, NORM_FROBENIUS, &normf));
53: PetscCall(MatNorm(mat, NORM_1, &norm1));
54: PetscCall(MatNorm(mat, NORM_INFINITY, &normi));
55: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "original: Frobenious norm = %g, one norm = %g, infinity norm = %g\n", (double)normf, (double)norm1, (double)normi));
56: PetscCall(MatView(mat, PETSC_VIEWER_STDOUT_WORLD));
58: /* Form matrix transpose */
59: PetscCall(PetscOptionsHasName(NULL, NULL, "-in_place", &flg));
60: if (flg) {
61: PetscCall(MatTranspose(mat, MAT_INPLACE_MATRIX, &mat)); /* in-place transpose */
62: tmat = mat;
63: mat = 0;
64: } else { /* out-of-place transpose */
65: PetscCall(MatTranspose(mat, MAT_INITIAL_MATRIX, &tmat));
66: }
68: /* Print info about transpose matrix */
69: PetscCall(MatGetInfo(tmat, MAT_GLOBAL_SUM, &info));
70: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "transpose matrix nonzeros = %" PetscInt_FMT ", allocated nonzeros = %" PetscInt_FMT "\n", (PetscInt)info.nz_used, (PetscInt)info.nz_allocated));
71: PetscCall(MatNorm(tmat, NORM_FROBENIUS, &normf));
72: PetscCall(MatNorm(tmat, NORM_1, &norm1));
73: PetscCall(MatNorm(tmat, NORM_INFINITY, &normi));
74: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "transpose: Frobenious norm = %g, one norm = %g, infinity norm = %g\n", (double)normf, (double)norm1, (double)normi));
75: PetscCall(MatView(tmat, PETSC_VIEWER_STDOUT_WORLD));
77: /* Test MatAXPY */
78: if (mat && !rect) {
79: PetscScalar alpha = 1.0;
80: PetscCall(PetscOptionsGetScalar(NULL, NULL, "-alpha", &alpha, NULL));
81: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "matrix addition: B = B + alpha * A\n"));
82: PetscCall(MatAXPY(tmat, alpha, mat, DIFFERENT_NONZERO_PATTERN));
83: PetscCall(MatView(tmat, PETSC_VIEWER_STDOUT_WORLD));
84: }
86: /* Free data structures */
87: PetscCall(MatDestroy(&tmat));
88: if (mat) PetscCall(MatDestroy(&mat));
90: PetscCall(PetscFinalize());
91: return 0;
92: }
94: /*TEST
96: test:
98: testset:
99: args: -rect1
100: test:
101: suffix: r1
102: output_file: output/ex49_r1.out
103: test:
104: suffix: r1_inplace
105: args: -in_place
106: output_file: output/ex49_r1.out
107: test:
108: suffix: r1_par
109: nsize: 2
110: output_file: output/ex49_r1_par.out
111: test:
112: suffix: r1_par_inplace
113: args: -in_place
114: nsize: 2
115: output_file: output/ex49_r1_par.out
117: TEST*/