Actual source code: ex98.c

  1: static char help[] = "Tests MatMPIAIJSetPreallocationCSR()\n\n";

  3: /*
  4:   Include "petscmat.h" so that we can use matrices.  Note that this file
  5:   automatically includes:
  6:      petscsys.h       - base PETSc routines   petscvec.h - vectors
  7:      petscmat.h - matrices
  8:      petscis.h     - index sets
  9:      petscviewer.h - viewers
 10: */
 11: #include <petscmat.h>

 13: int main(int argc, char **args)
 14: {
 15:   Mat         A;
 16:   PetscInt   *ia, *ja;
 17:   PetscMPIInt rank, size;

 19:   PetscFunctionBeginUser;
 20:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 21:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 22:   PetscCheck(size == 4, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run with 4 processors");
 23:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));

 25:   PetscCall(PetscMalloc1(5, &ia));
 26:   PetscCall(PetscMalloc1(16, &ja));
 27:   if (rank == 0) {
 28:     ja[0] = 1;
 29:     ja[1] = 4;
 30:     ja[2] = 0;
 31:     ja[3] = 2;
 32:     ja[4] = 5;
 33:     ja[5] = 1;
 34:     ja[6] = 3;
 35:     ja[7] = 6;
 36:     ja[8] = 2;
 37:     ja[9] = 7;
 38:     ia[0] = 0;
 39:     ia[1] = 2;
 40:     ia[2] = 5;
 41:     ia[3] = 8;
 42:     ia[4] = 10;
 43:   } else if (rank == 1) {
 44:     ja[0]  = 0;
 45:     ja[1]  = 5;
 46:     ja[2]  = 8;
 47:     ja[3]  = 1;
 48:     ja[4]  = 4;
 49:     ja[5]  = 6;
 50:     ja[6]  = 9;
 51:     ja[7]  = 2;
 52:     ja[8]  = 5;
 53:     ja[9]  = 7;
 54:     ja[10] = 10;
 55:     ja[11] = 3;
 56:     ja[12] = 6;
 57:     ja[13] = 11;
 58:     ia[0]  = 0;
 59:     ia[1]  = 3;
 60:     ia[2]  = 7;
 61:     ia[3]  = 11;
 62:     ia[4]  = 14;
 63:   } else if (rank == 2) {
 64:     ja[0]  = 4;
 65:     ja[1]  = 9;
 66:     ja[2]  = 12;
 67:     ja[3]  = 5;
 68:     ja[4]  = 8;
 69:     ja[5]  = 10;
 70:     ja[6]  = 13;
 71:     ja[7]  = 6;
 72:     ja[8]  = 9;
 73:     ja[9]  = 11;
 74:     ja[10] = 14;
 75:     ja[11] = 7;
 76:     ja[12] = 10;
 77:     ja[13] = 15;
 78:     ia[0]  = 0;
 79:     ia[1]  = 3;
 80:     ia[2]  = 7;
 81:     ia[3]  = 11;
 82:     ia[4]  = 14;
 83:   } else {
 84:     ja[0] = 8;
 85:     ja[1] = 13;
 86:     ja[2] = 9;
 87:     ja[3] = 12;
 88:     ja[4] = 14;
 89:     ja[5] = 10;
 90:     ja[6] = 13;
 91:     ja[7] = 15;
 92:     ja[8] = 11;
 93:     ja[9] = 14;
 94:     ia[0] = 0;
 95:     ia[1] = 2;
 96:     ia[2] = 5;
 97:     ia[3] = 8;
 98:     ia[4] = 10;
 99:   }

101:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
102:   PetscCall(MatSetSizes(A, 4, 4, 16, 16));
103:   PetscCall(MatSetType(A, MATMPIAIJ));
104:   PetscCall(MatMPIAIJSetPreallocationCSR(A, ia, ja, NULL));
105:   PetscCall(PetscFree(ia));
106:   PetscCall(PetscFree(ja));
107:   PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD));
108:   PetscCall(MatDestroy(&A));
109:   PetscCall(PetscFinalize());
110:   return 0;
111: }

113: /*TEST

115:    test:
116:       nsize: 4

118: TEST*/