Actual source code: ex137.c
1: static char help[] = "Tests MatCreateMPISBAIJWithArrays().\n\n";
3: #include <petscmat.h>
5: int main(int argc, char **args)
6: {
7: Mat A;
8: PetscInt m = 4, bs = 1, ii[5], jj[7];
9: PetscMPIInt size, rank;
10: PetscScalar aa[7];
12: PetscFunctionBeginUser;
13: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
14: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
15: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
16: PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Only for two processes");
18: if (rank == 0) {
19: ii[0] = 0;
20: ii[1] = 2;
21: ii[2] = 5;
22: ii[3] = 7;
23: ii[4] = 7;
24: jj[0] = 0;
25: jj[1] = 1;
26: jj[2] = 1;
27: jj[3] = 2;
28: jj[4] = 6;
29: jj[5] = 3;
30: jj[6] = 7;
31: aa[0] = 0;
32: aa[1] = 1;
33: aa[2] = 2;
34: aa[3] = 3;
35: aa[4] = 4;
36: aa[5] = 5;
37: aa[6] = 6;
38: /* 0 1
39: 1 2 6
40: 3 7 */
41: } else {
42: ii[0] = 0;
43: ii[1] = 2;
44: ii[2] = 4;
45: ii[3] = 6;
46: ii[4] = 7;
47: jj[0] = 4;
48: jj[1] = 5;
49: jj[2] = 5;
50: jj[3] = 7;
51: jj[4] = 6;
52: jj[5] = 7;
53: jj[6] = 7;
54: aa[0] = 8;
55: aa[1] = 9;
56: aa[2] = 10;
57: aa[3] = 11;
58: aa[4] = 12;
59: aa[5] = 13;
60: aa[6] = 14;
61: /* 4 5
62: 5 7
63: 6 7
64: 7 */
65: }
66: PetscCall(MatCreateMPISBAIJWithArrays(PETSC_COMM_WORLD, bs, m, m, PETSC_DECIDE, PETSC_DECIDE, ii, jj, aa, &A));
67: PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD));
68: PetscCall(MatDestroy(&A));
69: PetscCall(PetscFinalize());
70: return 0;
71: }
73: /*TEST
75: test:
76: nsize: 2
78: TEST*/