1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Tests MatMPIAIJSetPreallocationCSR()\n\n"; 3c4762a1bSJed Brown 4c4762a1bSJed Brown /* 5c4762a1bSJed Brown Include "petscmat.h" so that we can use matrices. Note that this file 6c4762a1bSJed Brown automatically includes: 7c4762a1bSJed Brown petscsys.h - base PETSc routines petscvec.h - vectors 8c4762a1bSJed Brown petscmat.h - matrices 9c4762a1bSJed Brown petscis.h - index sets 10c4762a1bSJed Brown petscviewer.h - viewers 11c4762a1bSJed Brown */ 12c4762a1bSJed Brown #include <petscmat.h> 13c4762a1bSJed Brown 14*9371c9d4SSatish Balay int main(int argc, char **args) { 15c4762a1bSJed Brown Mat A; 16c4762a1bSJed Brown PetscInt *ia, *ja; 17c4762a1bSJed Brown PetscMPIInt rank, size; 18c4762a1bSJed Brown 19327415f7SBarry Smith PetscFunctionBeginUser; 209566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 219566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 22be096a46SBarry Smith PetscCheck(size == 4, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run with 4 processors"); 239566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 24c4762a1bSJed Brown 259566063dSJacob Faibussowitsch PetscCall(PetscMalloc1(5, &ia)); 269566063dSJacob Faibussowitsch PetscCall(PetscMalloc1(16, &ja)); 27dd400576SPatrick Sanan if (rank == 0) { 28*9371c9d4SSatish Balay ja[0] = 1; 29*9371c9d4SSatish Balay ja[1] = 4; 30*9371c9d4SSatish Balay ja[2] = 0; 31*9371c9d4SSatish Balay ja[3] = 2; 32*9371c9d4SSatish Balay ja[4] = 5; 33*9371c9d4SSatish Balay ja[5] = 1; 34*9371c9d4SSatish Balay ja[6] = 3; 35*9371c9d4SSatish Balay ja[7] = 6; 36*9371c9d4SSatish Balay ja[8] = 2; 37*9371c9d4SSatish Balay ja[9] = 7; 38*9371c9d4SSatish Balay ia[0] = 0; 39*9371c9d4SSatish Balay ia[1] = 2; 40*9371c9d4SSatish Balay ia[2] = 5; 41*9371c9d4SSatish Balay ia[3] = 8; 42*9371c9d4SSatish Balay ia[4] = 10; 43c4762a1bSJed Brown } else if (rank == 1) { 44*9371c9d4SSatish Balay ja[0] = 0; 45*9371c9d4SSatish Balay ja[1] = 5; 46*9371c9d4SSatish Balay ja[2] = 8; 47*9371c9d4SSatish Balay ja[3] = 1; 48*9371c9d4SSatish Balay ja[4] = 4; 49*9371c9d4SSatish Balay ja[5] = 6; 50*9371c9d4SSatish Balay ja[6] = 9; 51*9371c9d4SSatish Balay ja[7] = 2; 52*9371c9d4SSatish Balay ja[8] = 5; 53*9371c9d4SSatish Balay ja[9] = 7; 54*9371c9d4SSatish Balay ja[10] = 10; 55*9371c9d4SSatish Balay ja[11] = 3; 56*9371c9d4SSatish Balay ja[12] = 6; 57*9371c9d4SSatish Balay ja[13] = 11; 58*9371c9d4SSatish Balay ia[0] = 0; 59*9371c9d4SSatish Balay ia[1] = 3; 60*9371c9d4SSatish Balay ia[2] = 7; 61*9371c9d4SSatish Balay ia[3] = 11; 62*9371c9d4SSatish Balay ia[4] = 14; 63c4762a1bSJed Brown } else if (rank == 2) { 64*9371c9d4SSatish Balay ja[0] = 4; 65*9371c9d4SSatish Balay ja[1] = 9; 66*9371c9d4SSatish Balay ja[2] = 12; 67*9371c9d4SSatish Balay ja[3] = 5; 68*9371c9d4SSatish Balay ja[4] = 8; 69*9371c9d4SSatish Balay ja[5] = 10; 70*9371c9d4SSatish Balay ja[6] = 13; 71*9371c9d4SSatish Balay ja[7] = 6; 72*9371c9d4SSatish Balay ja[8] = 9; 73*9371c9d4SSatish Balay ja[9] = 11; 74*9371c9d4SSatish Balay ja[10] = 14; 75*9371c9d4SSatish Balay ja[11] = 7; 76*9371c9d4SSatish Balay ja[12] = 10; 77*9371c9d4SSatish Balay ja[13] = 15; 78*9371c9d4SSatish Balay ia[0] = 0; 79*9371c9d4SSatish Balay ia[1] = 3; 80*9371c9d4SSatish Balay ia[2] = 7; 81*9371c9d4SSatish Balay ia[3] = 11; 82*9371c9d4SSatish Balay ia[4] = 14; 83c4762a1bSJed Brown } else { 84*9371c9d4SSatish Balay ja[0] = 8; 85*9371c9d4SSatish Balay ja[1] = 13; 86*9371c9d4SSatish Balay ja[2] = 9; 87*9371c9d4SSatish Balay ja[3] = 12; 88*9371c9d4SSatish Balay ja[4] = 14; 89*9371c9d4SSatish Balay ja[5] = 10; 90*9371c9d4SSatish Balay ja[6] = 13; 91*9371c9d4SSatish Balay ja[7] = 15; 92*9371c9d4SSatish Balay ja[8] = 11; 93*9371c9d4SSatish Balay ja[9] = 14; 94*9371c9d4SSatish Balay ia[0] = 0; 95*9371c9d4SSatish Balay ia[1] = 2; 96*9371c9d4SSatish Balay ia[2] = 5; 97*9371c9d4SSatish Balay ia[3] = 8; 98*9371c9d4SSatish Balay ia[4] = 10; 99c4762a1bSJed Brown } 100c4762a1bSJed Brown 1019566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 1029566063dSJacob Faibussowitsch PetscCall(MatSetSizes(A, 4, 4, 16, 16)); 1039566063dSJacob Faibussowitsch PetscCall(MatSetType(A, MATMPIAIJ)); 1049566063dSJacob Faibussowitsch PetscCall(MatMPIAIJSetPreallocationCSR(A, ia, ja, NULL)); 1059566063dSJacob Faibussowitsch PetscCall(PetscFree(ia)); 1069566063dSJacob Faibussowitsch PetscCall(PetscFree(ja)); 1079566063dSJacob Faibussowitsch PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 1089566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 1099566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 110b122ec5aSJacob Faibussowitsch return 0; 111c4762a1bSJed Brown } 112c4762a1bSJed Brown 113c4762a1bSJed Brown /*TEST 114c4762a1bSJed Brown 115c4762a1bSJed Brown test: 116c4762a1bSJed Brown nsize: 4 117c4762a1bSJed Brown 118c4762a1bSJed Brown TEST*/ 119