1 static char help[] = "Tests MATCENTERING matrix type.\n\n";
2
3 #include <petscmat.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 PetscInt n;
8 Mat C;
9 Vec x, y;
10 PetscReal norm;
11 PetscMPIInt size;
12
13 PetscFunctionBeginUser;
14 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
15 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
16
17 /* Create a parallel vector with 10*size total entries, and fill it with 1s. */
18 n = 10 * size;
19 PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
20 PetscCall(VecSetSizes(x, PETSC_DECIDE, n));
21 PetscCall(VecSetFromOptions(x));
22 PetscCall(VecSet(x, 1.0));
23
24 /* Create a corresponding n x n centering matrix and use it to create a mean-centered y = C * x. */
25 PetscCall(VecDuplicate(x, &y));
26 PetscCall(MatCreateCentering(PETSC_COMM_WORLD, PETSC_DECIDE, n, &C));
27 PetscCall(MatMult(C, x, y));
28
29 /* Verify that the centered vector y has norm 0. */
30 PetscCall(VecNorm(y, NORM_2, &norm));
31 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Vector norm after MatMult() with centering matrix applied to vector of ones is %f.\n", (double)norm));
32
33 /* Now repeat, but using MatMultTranspose(). */
34 PetscCall(MatMultTranspose(C, x, y));
35 PetscCall(VecNorm(y, NORM_2, &norm));
36 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Vector norm after MatMultTranspose() with centering matrix applied to vector of ones is %f.\n", (double)norm));
37
38 /* Clean up. */
39 PetscCall(VecDestroy(&x));
40 PetscCall(VecDestroy(&y));
41 PetscCall(MatDestroy(&C));
42 PetscCall(PetscFinalize());
43 return 0;
44 }
45
46 /*TEST
47
48 test:
49 suffix: 1
50 nsize: 1
51 output_file: output/ex247.out
52
53 test:
54 suffix: 2
55 nsize: 2
56 output_file: output/ex247.out
57
58 TEST*/
59