Actual source code: ex90.c


  2: static char help[] = "Tests MatPtAP() \n";

  4: #include <petscmat.h>

  6: /*
  7:  * This is an extremely simple example to test MatPtAP. It is very useful when developing and debugging the code.
  8:  *
  9:  * A =

 11:    1   2   0   4
 12:    0   1   2   0
 13:    2   0   4   0
 14:    0   1   2   1
 15:  *
 16:  *
 17:  *
 18:  * P =

 20:    1.00000   0.00000
 21:    0.30000   0.50000
 22:    0.00000   0.80000
 23:    0.90000   0.00000
 24:  *
 25:  *
 26:  *AP =
 27:  *
 28:  * 5.20000   1.00000
 29:    0.30000   2.10000
 30:    2.00000   3.20000
 31:    1.20000   2.10000
 32:  *
 33:  * PT =

 35:    1.00000   0.30000   0.00000   0.90000
 36:    0.00000   0.50000   0.80000   0.00000

 38:  *
 39:  * C =

 41:    6.3700   3.5200
 42:    1.7500   3.6100
 43:  *
 44:  * */

 46: int main(int argc, char **argv)
 47: {
 48:   Mat         A, P, PtAP;
 49:   PetscInt    i1[] = {0, 3, 5}, i2[] = {0, 2, 5};
 50:   PetscInt    j1[] = {0, 1, 3, 1, 2}, j2[] = {0, 2, 1, 2, 3};
 51:   PetscScalar a1[] = {1, 2, 4, 1, 2}, a2[] = {2, 4, 1, 2, 1};
 52:   PetscInt    pi1[] = {0, 1, 3}, pi2[] = {0, 1, 2};
 53:   PetscInt    pj1[] = {0, 0, 1}, pj2[] = {1, 0};
 54:   PetscScalar pa1[] = {1, 0.3, 0.5}, pa2[] = {0.8, 0.9};
 55:   MPI_Comm    comm;
 56:   PetscMPIInt rank, size;

 58:   PetscFunctionBeginUser;
 59:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 60:   comm = PETSC_COMM_WORLD;
 61:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
 62:   PetscCallMPI(MPI_Comm_size(comm, &size));
 63:   PetscCheck(size == 2, comm, PETSC_ERR_WRONG_MPI_SIZE, "You have to use two processor cores to run this example ");
 64:   PetscCall(MatCreateMPIAIJWithArrays(comm, 2, 2, PETSC_DETERMINE, PETSC_DETERMINE, rank ? i2 : i1, rank ? j2 : j1, rank ? a2 : a1, &A));
 65:   PetscCall(MatCreateMPIAIJWithArrays(comm, 2, 1, PETSC_DETERMINE, PETSC_DETERMINE, rank ? pi2 : pi1, rank ? pj2 : pj1, rank ? pa2 : pa1, &P));
 66:   PetscCall(MatPtAP(A, P, MAT_INITIAL_MATRIX, 1.1, &PtAP));
 67:   PetscCall(MatView(A, NULL));
 68:   PetscCall(MatView(P, NULL));
 69:   PetscCall(MatView(PtAP, NULL));
 70:   PetscCall(MatPtAP(A, P, MAT_REUSE_MATRIX, 1.1, &PtAP));
 71:   PetscCall(MatView(A, NULL));
 72:   PetscCall(MatView(P, NULL));
 73:   PetscCall(MatView(PtAP, NULL));
 74:   PetscCall(MatDestroy(&A));
 75:   PetscCall(MatDestroy(&P));
 76:   PetscCall(MatDestroy(&PtAP));
 77:   PetscCall(PetscFinalize());
 78:   return 0;
 79: }

 81: /*TEST
 82:    test:
 83:      nsize: 2
 84:      args:   -matptap_via allatonce
 85:      output_file: output/ex90_1.out

 87:    test:
 88:      nsize: 2
 89:      suffix: merged
 90:      args:   -matptap_via allatonce_merged
 91:      output_file: output/ex90_1.out

 93: TEST*/