1c4762a1bSJed Brown static char help[] = "Tests MatPtAP() \n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown /* 6c4762a1bSJed Brown * This is an extremely simple example to test MatPtAP. It is very useful when developing and debugging the code. 7c4762a1bSJed Brown * 8c4762a1bSJed Brown * A = 9c4762a1bSJed Brown 10c4762a1bSJed Brown 1 2 0 4 11c4762a1bSJed Brown 0 1 2 0 12c4762a1bSJed Brown 2 0 4 0 13c4762a1bSJed Brown 0 1 2 1 14c4762a1bSJed Brown * 15c4762a1bSJed Brown * 16c4762a1bSJed Brown * 17c4762a1bSJed Brown * P = 18c4762a1bSJed Brown 19c4762a1bSJed Brown 1.00000 0.00000 20c4762a1bSJed Brown 0.30000 0.50000 21c4762a1bSJed Brown 0.00000 0.80000 22c4762a1bSJed Brown 0.90000 0.00000 23c4762a1bSJed Brown * 24c4762a1bSJed Brown * 25c4762a1bSJed Brown *AP = 26c4762a1bSJed Brown * 27c4762a1bSJed Brown * 5.20000 1.00000 28c4762a1bSJed Brown 0.30000 2.10000 29c4762a1bSJed Brown 2.00000 3.20000 30c4762a1bSJed Brown 1.20000 2.10000 31c4762a1bSJed Brown * 32c4762a1bSJed Brown * PT = 33c4762a1bSJed Brown 34c4762a1bSJed Brown 1.00000 0.30000 0.00000 0.90000 35c4762a1bSJed Brown 0.00000 0.50000 0.80000 0.00000 36c4762a1bSJed Brown 37c4762a1bSJed Brown * 38c4762a1bSJed Brown * C = 39c4762a1bSJed Brown 40c4762a1bSJed Brown 6.3700 3.5200 41c4762a1bSJed Brown 1.7500 3.6100 42c4762a1bSJed Brown * 43c4762a1bSJed Brown * */ 44c4762a1bSJed Brown 45*d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv) 46*d71ae5a4SJacob Faibussowitsch { 47c4762a1bSJed Brown Mat A, P, PtAP; 48c4762a1bSJed Brown PetscInt i1[] = {0, 3, 5}, i2[] = {0, 2, 5}; 49c4762a1bSJed Brown PetscInt j1[] = {0, 1, 3, 1, 2}, j2[] = {0, 2, 1, 2, 3}; 50c4762a1bSJed Brown PetscScalar a1[] = {1, 2, 4, 1, 2}, a2[] = {2, 4, 1, 2, 1}; 51c4762a1bSJed Brown PetscInt pi1[] = {0, 1, 3}, pi2[] = {0, 1, 2}; 52c4762a1bSJed Brown PetscInt pj1[] = {0, 0, 1}, pj2[] = {1, 0}; 53c4762a1bSJed Brown PetscScalar pa1[] = {1, 0.3, 0.5}, pa2[] = {0.8, 0.9}; 54c4762a1bSJed Brown MPI_Comm comm; 55c4762a1bSJed Brown PetscMPIInt rank, size; 56c4762a1bSJed Brown 57327415f7SBarry Smith PetscFunctionBeginUser; 589566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 59c4762a1bSJed Brown comm = PETSC_COMM_WORLD; 609566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(comm, &rank)); 619566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(comm, &size)); 62be096a46SBarry Smith PetscCheck(size == 2, comm, PETSC_ERR_WRONG_MPI_SIZE, "You have to use two processor cores to run this example "); 639566063dSJacob Faibussowitsch PetscCall(MatCreateMPIAIJWithArrays(comm, 2, 2, PETSC_DETERMINE, PETSC_DETERMINE, rank ? i2 : i1, rank ? j2 : j1, rank ? a2 : a1, &A)); 649566063dSJacob Faibussowitsch PetscCall(MatCreateMPIAIJWithArrays(comm, 2, 1, PETSC_DETERMINE, PETSC_DETERMINE, rank ? pi2 : pi1, rank ? pj2 : pj1, rank ? pa2 : pa1, &P)); 659566063dSJacob Faibussowitsch PetscCall(MatPtAP(A, P, MAT_INITIAL_MATRIX, 1.1, &PtAP)); 669566063dSJacob Faibussowitsch PetscCall(MatView(A, NULL)); 679566063dSJacob Faibussowitsch PetscCall(MatView(P, NULL)); 689566063dSJacob Faibussowitsch PetscCall(MatView(PtAP, NULL)); 699566063dSJacob Faibussowitsch PetscCall(MatPtAP(A, P, MAT_REUSE_MATRIX, 1.1, &PtAP)); 709566063dSJacob Faibussowitsch PetscCall(MatView(A, NULL)); 719566063dSJacob Faibussowitsch PetscCall(MatView(P, NULL)); 729566063dSJacob Faibussowitsch PetscCall(MatView(PtAP, NULL)); 739566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 749566063dSJacob Faibussowitsch PetscCall(MatDestroy(&P)); 759566063dSJacob Faibussowitsch PetscCall(MatDestroy(&PtAP)); 769566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 77b122ec5aSJacob Faibussowitsch return 0; 78c4762a1bSJed Brown } 79c4762a1bSJed Brown 80c4762a1bSJed Brown /*TEST 81c4762a1bSJed Brown test: 82c4762a1bSJed Brown nsize: 2 83c4762a1bSJed Brown args: -matptap_via allatonce 84c4762a1bSJed Brown output_file: output/ex90_1.out 85c4762a1bSJed Brown 86c4762a1bSJed Brown test: 87c4762a1bSJed Brown nsize: 2 88c4762a1bSJed Brown suffix: merged 89c4762a1bSJed Brown args: -matptap_via allatonce_merged 90c4762a1bSJed Brown output_file: output/ex90_1.out 91c4762a1bSJed Brown 92c4762a1bSJed Brown TEST*/ 93