1 static char help[] = "Testing MatCreateMPIMatConcatenateSeqMat().\n\n"; 2 3 #include <petscmat.h> 4 int main(int argc, char **argv) { 5 Mat seqmat, mpimat; 6 PetscMPIInt rank; 7 PetscScalar value[3], *vals; 8 PetscInt i, col[3], n = 5, bs = 1; 9 10 PetscFunctionBeginUser; 11 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 12 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 13 PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL)); 14 15 /* Create seqaij matrices of size (n+rank) by n */ 16 PetscCall(MatCreate(PETSC_COMM_SELF, &seqmat)); 17 PetscCall(MatSetSizes(seqmat, (n + rank) * bs, PETSC_DECIDE, PETSC_DECIDE, n * bs)); 18 PetscCall(MatSetFromOptions(seqmat)); 19 PetscCall(MatSeqAIJSetPreallocation(seqmat, 3 * bs, NULL)); 20 PetscCall(MatSeqBAIJSetPreallocation(seqmat, bs, 3, NULL)); 21 22 if (bs == 1) { 23 value[0] = -1.0; 24 value[1] = 2.0; 25 value[2] = -1.0; 26 for (i = 1; i < n - 1; i++) { 27 col[0] = i - 1; 28 col[1] = i; 29 col[2] = i + 1; 30 PetscCall(MatSetValues(seqmat, 1, &i, 3, col, value, INSERT_VALUES)); 31 } 32 i = n - 1; 33 col[0] = n - 2; 34 col[1] = n - 1; 35 PetscCall(MatSetValues(seqmat, 1, &i, 2, col, value, INSERT_VALUES)); 36 37 i = 0; 38 col[0] = 0; 39 col[1] = 1; 40 value[0] = 2.0; 41 value[1] = -1.0; 42 PetscCall(MatSetValues(seqmat, 1, &i, 2, col, value, INSERT_VALUES)); 43 } else { 44 PetscInt *rows, *cols, j; 45 PetscCall(PetscMalloc3(bs * bs, &vals, bs, &rows, bs, &cols)); 46 /* diagonal blocks */ 47 for (i = 0; i < bs * bs; i++) vals[i] = 2.0; 48 for (i = 0; i < n * bs; i += bs) { 49 for (j = 0; j < bs; j++) { 50 rows[j] = i + j; 51 cols[j] = i + j; 52 } 53 PetscCall(MatSetValues(seqmat, bs, rows, bs, cols, vals, INSERT_VALUES)); 54 } 55 /* off-diagonal blocks */ 56 for (i = 0; i < bs * bs; i++) vals[i] = -1.0; 57 for (i = 0; i < (n - 1) * bs; i += bs) { 58 for (j = 0; j < bs; j++) { 59 rows[j] = i + j; 60 cols[j] = i + bs + j; 61 } 62 PetscCall(MatSetValues(seqmat, bs, rows, bs, cols, vals, INSERT_VALUES)); 63 } 64 65 PetscCall(PetscFree3(vals, rows, cols)); 66 } 67 PetscCall(MatAssemblyBegin(seqmat, MAT_FINAL_ASSEMBLY)); 68 PetscCall(MatAssemblyEnd(seqmat, MAT_FINAL_ASSEMBLY)); 69 if (rank == 0) { 70 PetscCall(PetscPrintf(PETSC_COMM_SELF, "[%d] seqmat:\n", rank)); 71 PetscCall(MatView(seqmat, PETSC_VIEWER_STDOUT_SELF)); 72 } 73 74 PetscCall(MatCreateMPIMatConcatenateSeqMat(PETSC_COMM_WORLD, seqmat, PETSC_DECIDE, MAT_INITIAL_MATRIX, &mpimat)); 75 PetscCall(MatCreateMPIMatConcatenateSeqMat(PETSC_COMM_WORLD, seqmat, PETSC_DECIDE, MAT_REUSE_MATRIX, &mpimat)); 76 PetscCall(MatView(mpimat, PETSC_VIEWER_STDOUT_WORLD)); 77 78 PetscCall(MatDestroy(&seqmat)); 79 PetscCall(MatDestroy(&mpimat)); 80 PetscCall(PetscFinalize()); 81 return 0; 82 } 83 84 /*TEST 85 86 test: 87 nsize: 3 88 89 test: 90 suffix: 2 91 nsize: 3 92 args: -mat_type baij 93 94 test: 95 suffix: 3 96 nsize: 3 97 args: -mat_type baij -bs 2 98 99 TEST*/ 100