Actual source code: ex176.c

  1: static char help[] = "Tests MatCreateMPIAIJWithArrays() and MatUpdateMPIAIJWithArray()\n";

  3: #include <petscmat.h>

  5: /*
  6:  * This is an extremely simple example to test MatUpdateMPIAIJWithArrays()
  7:  *
  8:  * A =

 10:    1    2   0   3  0  0
 11:    0    4   5   0  0  6
 12:    7    0   8   0  9  0
 13:    0   10  11  12  0  13
 14:    0   14  15   0  0  16
 15:   17    0   0   0  0  18
 16:  *
 17:  * */

 19: int main(int argc, char **argv)
 20: {
 21:   Mat      A, B;
 22:   PetscInt i[3][3] = {
 23:     {0, 3, 6},
 24:     {0, 3, 7},
 25:     {0, 3, 5}
 26:   };
 27:   PetscInt j[3][7] = {
 28:     {0, 1, 3, 1, 2, 5,  -1},
 29:     {0, 2, 4, 1, 2, 3,  5 },
 30:     {1, 2, 5, 0, 5, -1, -1}
 31:   };
 32:   PetscScalar a[3][7] = {
 33:     {1,  2,  3,  4,  5,  6,  -1},
 34:     {7,  8,  9,  10, 11, 12, 13},
 35:     {14, 15, 16, 17, 18, -1, -1}
 36:   };
 37:   PetscScalar anew[3][7] = {
 38:     {10,  20,  30,  40,  50,  60,  -1 },
 39:     {70,  80,  90,  100, 110, 120, 130},
 40:     {140, 150, 160, 170, 180, -1,  -1 }
 41:   };
 42:   MPI_Comm    comm;
 43:   PetscMPIInt rank;
 44:   PetscBool   equal = PETSC_FALSE;

 46:   PetscFunctionBeginUser;
 47:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 48:   comm = PETSC_COMM_WORLD;
 49:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
 50:   PetscCall(MatCreateMPIAIJWithArrays(comm, 2, PETSC_DETERMINE, PETSC_DETERMINE, 6, i[rank], j[rank], a[rank], &B));

 52:   PetscCall(MatCreateMPIAIJWithArrays(comm, 2, PETSC_DETERMINE, PETSC_DETERMINE, 6, i[rank], j[rank], a[rank], &A));
 53:   PetscCall(MatSetFromOptions(A)); /* might change A's type */

 55:   PetscCall(MatEqual(A, B, &equal));
 56:   PetscCheck(equal, comm, PETSC_ERR_PLIB, "wrong results");

 58:   PetscCall(MatUpdateMPIAIJWithArray(A, anew[rank]));
 59:   PetscCall(MatUpdateMPIAIJWithArray(B, anew[rank]));
 60:   PetscCall(MatEqual(A, B, &equal));
 61:   PetscCheck(equal, comm, PETSC_ERR_PLIB, "wrong results");

 63:   PetscCall(MatUpdateMPIAIJWithArray(A, a[rank]));
 64:   PetscCall(MatUpdateMPIAIJWithArray(B, a[rank]));
 65:   PetscCall(MatEqual(A, B, &equal));
 66:   PetscCheck(equal, comm, PETSC_ERR_PLIB, "wrong results");

 68:   PetscCall(MatDestroy(&A));
 69:   PetscCall(MatDestroy(&B));
 70:   PetscCall(PetscFinalize());
 71:   return 0;
 72: }

 74: /*TEST
 75:    testset:
 76:      nsize: {{1 3}}
 77:      output_file: output/empty.out

 79:      test:
 80:        suffix: aij

 82:      test:
 83:        requires: cuda
 84:        suffix: cuda
 85:        # since the matrices are created with MatCreateMPIxxx(), users are allowed to pass 'mpiaijcusparse' even with one rank
 86:        args: -mat_type {{aijcusparse mpiaijcusparse}}

 88:      test:
 89:        requires: kokkos_kernels
 90:        suffix: kok
 91:        args: -mat_type aijkokkos
 92: TEST*/