Actual source code: ex247.c

  1: static char help[] = "Tests MATCENTERING matrix type.\n\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **argv)
  6: {
  7:   PetscInt    n;
  8:   Mat         C;
  9:   Vec         x, y;
 10:   PetscReal   norm;
 11:   PetscMPIInt size;

 13:   PetscFunctionBeginUser;
 14:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 15:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));

 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));

 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));

 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));

 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));

 38:   /* Clean up. */
 39:   PetscCall(VecDestroy(&x));
 40:   PetscCall(VecDestroy(&y));
 41:   PetscCall(MatDestroy(&C));
 42:   PetscCall(PetscFinalize());
 43:   return 0;
 44: }

 46: /*TEST

 48:     test:
 49:       suffix: 1
 50:       nsize: 1
 51:       output_file: output/ex247.out

 53:     test:
 54:       suffix: 2
 55:       nsize: 2
 56:       output_file: output/ex247.out

 58: TEST*/