1 static char help[] = "Tests resetting preallocation after filling the full sparsity pattern";
2
3 #include <petscmat.h>
4
Assemble(Mat mat)5 PetscErrorCode Assemble(Mat mat)
6 {
7 PetscInt idx[4], i;
8 PetscScalar vals[16];
9 int rank;
10
11 PetscFunctionBegin;
12 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
13 for (i = 0; i < 16; ++i) vals[i] = 1;
14 if (rank == 0) {
15 // element 0
16 idx[0] = 0;
17 idx[1] = 1;
18 idx[2] = 2;
19 idx[3] = 3;
20 PetscCall(MatSetValues(mat, 4, idx, 4, idx, vals, ADD_VALUES));
21 // element 1
22 idx[0] = 3;
23 idx[1] = 2;
24 idx[2] = 4;
25 idx[3] = 5;
26 PetscCall(MatSetValues(mat, 4, idx, 4, idx, vals, ADD_VALUES));
27 } else {
28 // element 2
29 idx[0] = 6;
30 idx[1] = 0;
31 idx[2] = 3;
32 idx[3] = 7;
33 PetscCall(MatSetValues(mat, 4, idx, 4, idx, vals, ADD_VALUES));
34 // element 3
35 idx[0] = 7;
36 idx[1] = 3;
37 idx[2] = 5;
38 idx[3] = 8;
39 PetscCall(MatSetValues(mat, 4, idx, 4, idx, vals, ADD_VALUES));
40 }
41 PetscCall(MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY));
42 PetscCall(MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY));
43 PetscFunctionReturn(PETSC_SUCCESS);
44 }
45
main(int argc,char ** argv)46 int main(int argc, char **argv)
47 {
48 Mat mat;
49 int rank;
50
51 PetscFunctionBeginUser;
52 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
53 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
54 PetscCall(MatCreate(PETSC_COMM_WORLD, &mat));
55 if (rank == 0) PetscCall(MatSetSizes(mat, 6, 6, PETSC_DETERMINE, PETSC_DETERMINE));
56 else PetscCall(MatSetSizes(mat, 3, 3, PETSC_DETERMINE, PETSC_DETERMINE));
57 PetscCall(MatSetFromOptions(mat));
58 if (rank == 0) {
59 PetscInt ndz[6], noz[6];
60 ndz[0] = 4;
61 noz[0] = 2;
62 ndz[1] = 4;
63 noz[1] = 0;
64 ndz[2] = 6;
65 noz[2] = 0;
66 ndz[3] = 6;
67 noz[3] = 3;
68 ndz[4] = 4;
69 noz[4] = 0;
70 ndz[5] = 4;
71 noz[5] = 2;
72 PetscCall(MatMPIAIJSetPreallocation(mat, 0, ndz, 0, noz));
73 } else {
74 PetscInt ndz[3], noz[3];
75 ndz[0] = 2;
76 noz[0] = 2;
77 ndz[1] = 3;
78 noz[1] = 3;
79 ndz[2] = 2;
80 noz[2] = 2;
81 PetscCall(MatMPIAIJSetPreallocation(mat, 0, ndz, 0, noz));
82 }
83 PetscCall(MatSetUp(mat));
84 PetscCall(Assemble(mat));
85 PetscCall(MatView(mat, NULL));
86 PetscCall(MatResetPreallocation(mat));
87 PetscCall(Assemble(mat));
88 PetscCall(MatView(mat, NULL));
89 PetscCall(MatDestroy(&mat));
90 PetscCall(PetscFinalize());
91 return 0;
92 }
93
94 /*TEST
95
96 test:
97 suffix: 1
98 nsize: 2
99
100 TEST*/
101