1c4762a1bSJed Brown static char help[] = "Tests MatCreateMPISBAIJWithArrays().\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown 5d71ae5a4SJacob Faibussowitsch int main(int argc, char **args) 6d71ae5a4SJacob Faibussowitsch { 7c4762a1bSJed Brown Mat A; 8c4762a1bSJed Brown PetscInt m = 4, bs = 1, ii[5], jj[7]; 9c4762a1bSJed Brown PetscMPIInt size, rank; 10c4762a1bSJed Brown PetscScalar aa[7]; 11c4762a1bSJed Brown 12327415f7SBarry Smith PetscFunctionBeginUser; 13*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help)); 149566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 159566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 16be096a46SBarry Smith PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Only for two processes"); 17c4762a1bSJed Brown 18dd400576SPatrick Sanan if (rank == 0) { 199371c9d4SSatish Balay ii[0] = 0; 209371c9d4SSatish Balay ii[1] = 2; 219371c9d4SSatish Balay ii[2] = 5; 229371c9d4SSatish Balay ii[3] = 7; 239371c9d4SSatish Balay ii[4] = 7; 249371c9d4SSatish Balay jj[0] = 0; 259371c9d4SSatish Balay jj[1] = 1; 269371c9d4SSatish Balay jj[2] = 1; 279371c9d4SSatish Balay jj[3] = 2; 289371c9d4SSatish Balay jj[4] = 6; 299371c9d4SSatish Balay jj[5] = 3; 309371c9d4SSatish Balay jj[6] = 7; 319371c9d4SSatish Balay aa[0] = 0; 329371c9d4SSatish Balay aa[1] = 1; 339371c9d4SSatish Balay aa[2] = 2; 349371c9d4SSatish Balay aa[3] = 3; 359371c9d4SSatish Balay aa[4] = 4; 369371c9d4SSatish Balay aa[5] = 5; 379371c9d4SSatish Balay aa[6] = 6; 38c4762a1bSJed Brown /* 0 1 39c4762a1bSJed Brown 1 2 6 40c4762a1bSJed Brown 3 7 */ 41c4762a1bSJed Brown } else { 429371c9d4SSatish Balay ii[0] = 0; 439371c9d4SSatish Balay ii[1] = 2; 449371c9d4SSatish Balay ii[2] = 4; 459371c9d4SSatish Balay ii[3] = 6; 469371c9d4SSatish Balay ii[4] = 7; 479371c9d4SSatish Balay jj[0] = 4; 489371c9d4SSatish Balay jj[1] = 5; 499371c9d4SSatish Balay jj[2] = 5; 509371c9d4SSatish Balay jj[3] = 7; 519371c9d4SSatish Balay jj[4] = 6; 529371c9d4SSatish Balay jj[5] = 7; 539371c9d4SSatish Balay jj[6] = 7; 549371c9d4SSatish Balay aa[0] = 8; 559371c9d4SSatish Balay aa[1] = 9; 569371c9d4SSatish Balay aa[2] = 10; 579371c9d4SSatish Balay aa[3] = 11; 589371c9d4SSatish Balay aa[4] = 12; 599371c9d4SSatish Balay aa[5] = 13; 609371c9d4SSatish Balay aa[6] = 14; 61c4762a1bSJed Brown /* 4 5 62c4762a1bSJed Brown 5 7 63c4762a1bSJed Brown 6 7 64c4762a1bSJed Brown 7 */ 65c4762a1bSJed Brown } 669566063dSJacob Faibussowitsch PetscCall(MatCreateMPISBAIJWithArrays(PETSC_COMM_WORLD, bs, m, m, PETSC_DECIDE, PETSC_DECIDE, ii, jj, aa, &A)); 679566063dSJacob Faibussowitsch PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 689566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 699566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 70b122ec5aSJacob Faibussowitsch return 0; 71c4762a1bSJed Brown } 72c4762a1bSJed Brown 73c4762a1bSJed Brown /*TEST 74c4762a1bSJed Brown 75c4762a1bSJed Brown test: 76c4762a1bSJed Brown nsize: 2 77c4762a1bSJed Brown 78c4762a1bSJed Brown TEST*/ 79