1c4762a1bSJed Brown static char help[] = "Tests various routines in MatMPIBAIJ format.\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 = 2, bs = 1, M, row, col, start, end, i, j, k; 9c4762a1bSJed Brown PetscMPIInt rank, size; 10c4762a1bSJed Brown PetscScalar data = 100; 11c4762a1bSJed Brown PetscBool flg; 12c4762a1bSJed Brown 13327415f7SBarry Smith PetscFunctionBeginUser; 14*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help)); 159566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 169566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 17c4762a1bSJed Brown 18c4762a1bSJed Brown /* Test MatSetValues() and MatGetValues() */ 199566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-mat_block_size", &bs, NULL)); 209566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-mat_size", &m, NULL)); 21c4762a1bSJed Brown 22c4762a1bSJed Brown M = m * bs * size; 239566063dSJacob Faibussowitsch PetscCall(MatCreateBAIJ(PETSC_COMM_WORLD, bs, PETSC_DECIDE, PETSC_DECIDE, M, M, PETSC_DECIDE, NULL, PETSC_DECIDE, NULL, &A)); 24c4762a1bSJed Brown 259566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRange(A, &start, &end)); 269566063dSJacob Faibussowitsch PetscCall(PetscOptionsHasName(NULL, NULL, "-column_oriented", &flg)); 271baa6e33SBarry Smith if (flg) PetscCall(MatSetOption(A, MAT_ROW_ORIENTED, PETSC_FALSE)); 28c4762a1bSJed Brown 29c4762a1bSJed Brown /* inproc assembly */ 30c4762a1bSJed Brown for (row = start; row < end; row++) { 3148a46eb9SPierre Jolivet for (col = start; col < end; col++, data += 1) PetscCall(MatSetValues(A, 1, &row, 1, &col, &data, INSERT_VALUES)); 32c4762a1bSJed Brown } 339566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 349566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 35c4762a1bSJed Brown 36c4762a1bSJed Brown /* offproc assembly */ 37c4762a1bSJed Brown data = 5.0; 38c4762a1bSJed Brown row = (M + start - 1) % M; 3948a46eb9SPierre Jolivet for (col = 0; col < M; col++) PetscCall(MatSetValues(A, 1, &row, 1, &col, &data, ADD_VALUES)); 409566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 419566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 42c4762a1bSJed Brown 43c4762a1bSJed Brown /* Test MatSetValuesBlocked() */ 449566063dSJacob Faibussowitsch PetscCall(PetscOptionsHasName(NULL, NULL, "-test_setvaluesblocked", &flg)); 45c4762a1bSJed Brown if (flg) { 46c4762a1bSJed Brown PetscScalar *bval; 47c4762a1bSJed Brown row /= bs; 48c4762a1bSJed Brown col = start / bs; 499566063dSJacob Faibussowitsch PetscCall(PetscMalloc1(bs * bs, &bval)); 50c4762a1bSJed Brown k = 1; 515e116b59SBarry Smith /* row-oriented - default */ 52c4762a1bSJed Brown for (i = 0; i < bs; i++) { 53c4762a1bSJed Brown for (j = 0; j < bs; j++) { 549371c9d4SSatish Balay bval[i * bs + j] = (PetscScalar)k; 559371c9d4SSatish Balay k++; 56c4762a1bSJed Brown } 57c4762a1bSJed Brown } 589566063dSJacob Faibussowitsch PetscCall(MatSetValuesBlocked(A, 1, &row, 1, &col, bval, INSERT_VALUES)); 599566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 609566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 619566063dSJacob Faibussowitsch PetscCall(PetscFree(bval)); 62c4762a1bSJed Brown } 63c4762a1bSJed Brown 649566063dSJacob Faibussowitsch PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 659566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 669566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 67b122ec5aSJacob Faibussowitsch return 0; 68c4762a1bSJed Brown } 69c4762a1bSJed Brown 70c4762a1bSJed Brown /*TEST 71c4762a1bSJed Brown 72c4762a1bSJed Brown test: 73c4762a1bSJed Brown suffix: 1 74c4762a1bSJed Brown nsize: 3 75c4762a1bSJed Brown args: -mat_block_size 2 -test_setvaluesblocked 76c4762a1bSJed Brown 77c4762a1bSJed Brown test: 78c4762a1bSJed Brown suffix: 2 79c4762a1bSJed Brown nsize: 3 80c4762a1bSJed Brown args: -mat_block_size 2 -test_setvaluesblocked -column_oriented 81c4762a1bSJed Brown 82c4762a1bSJed Brown test: 83c4762a1bSJed Brown suffix: 3 84c4762a1bSJed Brown nsize: 3 85c4762a1bSJed Brown args: -mat_block_size 1 -test_setvaluesblocked 86c4762a1bSJed Brown 87c4762a1bSJed Brown test: 88c4762a1bSJed Brown suffix: 4 89c4762a1bSJed Brown nsize: 3 90c4762a1bSJed Brown args: -mat_block_size 1 -test_setvaluesblocked -column_oriented 91c4762a1bSJed Brown 92c4762a1bSJed Brown TEST*/ 93