1 static char help[] = "Test MatCreateRedundantMatrix for rectangular matrix.\n\ 2 Contributed by Jose E. Roman, July 2017\n\n"; 3 4 #include <petscmat.h> 5 int main(int argc, char **args) { 6 Mat A, B; 7 PetscInt m = 3, n = 4, i, nsubcomm; 8 PetscMPIInt size, rank; 9 10 PetscFunctionBeginUser; 11 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 12 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 13 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 14 15 nsubcomm = size; 16 PetscCall(PetscOptionsGetInt(NULL, NULL, "-nsubcomm", &nsubcomm, NULL)); 17 18 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 19 PetscCall(MatSetSizes(A, m, n, PETSC_DETERMINE, PETSC_DETERMINE)); 20 PetscCall(MatSetType(A, MATAIJ)); 21 PetscCall(MatSetFromOptions(A)); 22 PetscCall(MatSetUp(A)); 23 24 if (rank == 0) { 25 for (i = 0; i < size * PetscMin(m, n); i++) { PetscCall(MatSetValue(A, i, i, 1.0, INSERT_VALUES)); } 26 } 27 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 28 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 29 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 30 31 PetscCall(MatCreateRedundantMatrix(A, nsubcomm, MPI_COMM_NULL, MAT_INITIAL_MATRIX, &B)); 32 if (nsubcomm == size) { /* B is a sequential matrix */ 33 if (rank == 0) { PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF)); } 34 } else { 35 MPI_Comm comm; 36 PetscCall(PetscObjectGetComm((PetscObject)B, &comm)); 37 PetscCall(MatView(B, PETSC_VIEWER_STDOUT_(comm))); 38 } 39 40 PetscCall(MatDestroy(&A)); 41 PetscCall(MatDestroy(&B)); 42 PetscCall(PetscFinalize()); 43 return 0; 44 } 45 46 /*TEST 47 48 test: 49 50 test: 51 suffix: 2 52 nsize: 3 53 54 test: 55 suffix: baij 56 args: -mat_type baij 57 58 test: 59 suffix: baij_2 60 nsize: 3 61 args: -mat_type baij 62 63 test: 64 suffix: dense 65 args: -mat_type dense 66 67 test: 68 suffix: dense_2 69 nsize: 3 70 args: -mat_type dense 71 72 TEST*/ 73