1c1fff1c6SRichard Tran Mills 2c1fff1c6SRichard Tran Mills static char help[] = "Tests MATCENTERING matrix type.\n\n"; 3c1fff1c6SRichard Tran Mills 4c1fff1c6SRichard Tran Mills #include <petscmat.h> 5c1fff1c6SRichard Tran Mills 6*d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv) 7*d71ae5a4SJacob Faibussowitsch { 8c1fff1c6SRichard Tran Mills PetscInt n; 9c1fff1c6SRichard Tran Mills Mat C; 10c1fff1c6SRichard Tran Mills Vec x, y; 11c1fff1c6SRichard Tran Mills PetscReal norm; 12c1fff1c6SRichard Tran Mills PetscMPIInt size; 13c1fff1c6SRichard Tran Mills 14327415f7SBarry Smith PetscFunctionBeginUser; 159566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 169566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 17c1fff1c6SRichard Tran Mills 18c1fff1c6SRichard Tran Mills /* Create a parallel vector with 10*size total entries, and fill it with 1s. */ 19c1fff1c6SRichard Tran Mills n = 10 * size; 209566063dSJacob Faibussowitsch PetscCall(VecCreate(PETSC_COMM_WORLD, &x)); 219566063dSJacob Faibussowitsch PetscCall(VecSetSizes(x, PETSC_DECIDE, n)); 229566063dSJacob Faibussowitsch PetscCall(VecSetFromOptions(x)); 239566063dSJacob Faibussowitsch PetscCall(VecSet(x, 1.0)); 24c1fff1c6SRichard Tran Mills 25c1fff1c6SRichard Tran Mills /* Create a corresponding n x n centering matrix and use it to create a mean-centered y = C * x. */ 269566063dSJacob Faibussowitsch PetscCall(VecDuplicate(x, &y)); 279566063dSJacob Faibussowitsch PetscCall(MatCreateCentering(PETSC_COMM_WORLD, PETSC_DECIDE, n, &C)); 289566063dSJacob Faibussowitsch PetscCall(MatMult(C, x, y)); 29c1fff1c6SRichard Tran Mills 30c1fff1c6SRichard Tran Mills /* Verify that the centered vector y has norm 0. */ 319566063dSJacob Faibussowitsch PetscCall(VecNorm(y, NORM_2, &norm)); 329566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Vector norm after MatMult() with centering matrix applied to vector of ones is %f.\n", (double)norm)); 33c1fff1c6SRichard Tran Mills 34c1fff1c6SRichard Tran Mills /* Now repeat, but using MatMultTranspose(). */ 359566063dSJacob Faibussowitsch PetscCall(MatMultTranspose(C, x, y)); 369566063dSJacob Faibussowitsch PetscCall(VecNorm(y, NORM_2, &norm)); 379566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Vector norm after MatMultTranspose() with centering matrix applied to vector of ones is %f.\n", (double)norm)); 38c1fff1c6SRichard Tran Mills 39c1fff1c6SRichard Tran Mills /* Clean up. */ 409566063dSJacob Faibussowitsch PetscCall(VecDestroy(&x)); 419566063dSJacob Faibussowitsch PetscCall(VecDestroy(&y)); 429566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 439566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 44b122ec5aSJacob Faibussowitsch return 0; 45c1fff1c6SRichard Tran Mills } 46c1fff1c6SRichard Tran Mills 47c1fff1c6SRichard Tran Mills /*TEST 48c1fff1c6SRichard Tran Mills 49c1fff1c6SRichard Tran Mills test: 50c1fff1c6SRichard Tran Mills suffix: 1 51c1fff1c6SRichard Tran Mills nsize: 1 52c1fff1c6SRichard Tran Mills output_file: output/ex247.out 53c1fff1c6SRichard Tran Mills 54c1fff1c6SRichard Tran Mills test: 55c1fff1c6SRichard Tran Mills suffix: 2 56c1fff1c6SRichard Tran Mills nsize: 2 57c1fff1c6SRichard Tran Mills output_file: output/ex247.out 58c1fff1c6SRichard Tran Mills 59c1fff1c6SRichard Tran Mills TEST*/ 60