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