Actual source code: ex230.c
1: static char help[] = "Example of using MatPreallocator\n\n";
3: /*T
4: Concepts: Mat^matrix preallocation
5: Processors: n
6: T*/
8: #include <petscmat.h>
10: PetscErrorCode ex1_nonsquare_bs1(void)
11: {
12: Mat A,preallocator;
13: PetscInt M,N,m,n,bs;
15: /*
16: Create the Jacobian matrix
17: */
18: M = 10;
19: N = 8;
20: MatCreate(PETSC_COMM_WORLD,&A);
21: MatSetType(A,MATAIJ);
22: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);
23: MatSetBlockSize(A,1);
24: MatSetFromOptions(A);
26: /*
27: Get the sizes of the jacobian.
28: */
29: MatGetLocalSize(A,&m,&n);
30: MatGetBlockSize(A,&bs);
32: /*
33: Create a preallocator matrix with sizes (local and global) matching the jacobian A.
34: */
35: MatCreate(PetscObjectComm((PetscObject)A),&preallocator);
36: MatSetType(preallocator,MATPREALLOCATOR);
37: MatSetSizes(preallocator,m,n,M,N);
38: MatSetBlockSize(preallocator,bs);
39: MatSetUp(preallocator);
41: /*
42: Insert non-zero pattern (e.g. perform a sweep over the grid).
43: You can use MatSetValues(), MatSetValuesBlocked() or MatSetValue().
44: */
45: {
46: PetscInt ii,jj;
47: PetscScalar vv = 0.0;
49: ii = 3; jj = 3;
50: MatSetValues(preallocator,1,&ii,1,&jj,&vv,INSERT_VALUES);
52: ii = 7; jj = 4;
53: MatSetValues(preallocator,1,&ii,1,&jj,&vv,INSERT_VALUES);
55: ii = 9; jj = 7;
56: MatSetValue(preallocator,ii,jj,vv,INSERT_VALUES);
57: }
58: MatAssemblyBegin(preallocator,MAT_FINAL_ASSEMBLY);
59: MatAssemblyEnd(preallocator,MAT_FINAL_ASSEMBLY);
61: /*
62: Push the non-zero pattern defined within preallocator into A.
63: Internally this will call MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE).
64: The arg fill = PETSC_TRUE will insert zeros in the matrix A automatically.
65: */
66: MatPreallocatorPreallocate(preallocator,PETSC_TRUE,A);
68: /*
69: We no longer require the preallocator object so destroy it.
70: */
71: MatDestroy(&preallocator);
73: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
75: /*
76: Insert non-zero values into A.
77: */
78: {
79: PetscInt ii,jj;
80: PetscScalar vv;
82: ii = 3; jj = 3; vv = 0.3;
83: MatSetValue(A,ii,jj,vv,INSERT_VALUES);
85: ii = 7; jj = 4; vv = 3.3;
86: MatSetValues(A,1,&ii,1,&jj,&vv,INSERT_VALUES);
88: ii = 9; jj = 7; vv = 4.3;
89: MatSetValues(A,1,&ii,1,&jj,&vv,INSERT_VALUES);
90: }
91: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
92: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
94: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
96: MatDestroy(&A);
97: return 0;
98: }
100: PetscErrorCode ex2_square_bsvariable(void)
101: {
102: Mat A,preallocator;
103: PetscInt M,N,m,n,bs = 1;
105: PetscOptionsGetInt(NULL,NULL,"-block_size",&bs,NULL);
107: /*
108: Create the Jacobian matrix.
109: */
110: M = 10 * bs;
111: N = 10 * bs;
112: MatCreate(PETSC_COMM_WORLD,&A);
113: MatSetType(A,MATAIJ);
114: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);
115: MatSetBlockSize(A,bs);
116: MatSetFromOptions(A);
118: /*
119: Get the sizes of the jacobian.
120: */
121: MatGetLocalSize(A,&m,&n);
122: MatGetBlockSize(A,&bs);
124: /*
125: Create a preallocator matrix with dimensions matching the jacobian A.
126: */
127: MatCreate(PetscObjectComm((PetscObject)A),&preallocator);
128: MatSetType(preallocator,MATPREALLOCATOR);
129: MatSetSizes(preallocator,m,n,M,N);
130: MatSetBlockSize(preallocator,bs);
131: MatSetUp(preallocator);
133: /*
134: Insert non-zero pattern (e.g. perform a sweep over the grid).
135: You can use MatSetValues(), MatSetValuesBlocked() or MatSetValue().
136: */
137: {
138: PetscInt ii,jj;
139: PetscScalar *vv;
141: PetscCalloc1(bs*bs,&vv);
143: ii = 0; jj = 9;
144: MatSetValue(preallocator,ii,jj,vv[0],INSERT_VALUES);
146: ii = 0; jj = 0;
147: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
149: ii = 2; jj = 4;
150: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
152: ii = 4; jj = 2;
153: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
155: ii = 4; jj = 4;
156: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
158: ii = 9; jj = 9;
159: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
161: PetscFree(vv);
162: }
163: MatAssemblyBegin(preallocator,MAT_FINAL_ASSEMBLY);
164: MatAssemblyEnd(preallocator,MAT_FINAL_ASSEMBLY);
166: /*
167: Push non-zero pattern defined from preallocator into A.
168: Internally this will call MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE).
169: The arg fill = PETSC_TRUE will insert zeros in the matrix A automatically.
170: */
171: MatPreallocatorPreallocate(preallocator,PETSC_TRUE,A);
173: /*
174: We no longer require the preallocator object so destroy it.
175: */
176: MatDestroy(&preallocator);
178: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
180: {
181: PetscInt ii,jj;
182: PetscScalar *vv;
184: PetscCalloc1(bs*bs,&vv);
185: for (ii=0; ii<bs*bs; ii++) {
186: vv[ii] = (PetscReal)(ii+1);
187: }
189: ii = 0; jj = 9;
190: MatSetValue(A,ii,jj,33.3,INSERT_VALUES);
192: ii = 0; jj = 0;
193: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
195: ii = 2; jj = 4;
196: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
198: ii = 4; jj = 2;
199: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
201: ii = 4; jj = 4;
202: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
204: ii = 9; jj = 9;
205: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
207: PetscFree(vv);
208: }
209: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
210: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
212: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
214: MatDestroy(&A);
215: return 0;
216: }
218: int main(int argc, char **args)
219: {
220: PetscInt testid = 0;
221: PetscInitialize(&argc,&args,(char*)0,help);
222: PetscOptionsGetInt(NULL,NULL,"-test_id",&testid,NULL);
223: switch (testid) {
224: case 0:
225: ex1_nonsquare_bs1();
226: break;
227: case 1:
228: ex2_square_bsvariable();
229: break;
230: default:
231: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Invalid value for -test_id. Must be {0,1}");
232: }
233: PetscFinalize();
234: return 0;
235: }
237: /*TEST
239: test:
240: suffix: t0_a_aij
241: nsize: 1
242: args: -test_id 0 -mat_type aij
244: test:
245: suffix: t0_b_aij
246: nsize: 6
247: args: -test_id 0 -mat_type aij
249: test:
250: suffix: t1_a_aij
251: nsize: 1
252: args: -test_id 1 -mat_type aij
254: test:
255: suffix: t2_a_baij
256: nsize: 1
257: args: -test_id 1 -mat_type baij
259: test:
260: suffix: t3_a_sbaij
261: nsize: 1
262: args: -test_id 1 -mat_type sbaij
264: test:
265: suffix: t4_a_aij_bs3
266: nsize: 1
267: args: -test_id 1 -mat_type aij -block_size 3
269: test:
270: suffix: t5_a_baij_bs3
271: nsize: 1
272: args: -test_id 1 -mat_type baij -block_size 3
274: test:
275: suffix: t6_a_sbaij_bs3
276: nsize: 1
277: args: -test_id 1 -mat_type sbaij -block_size 3
279: test:
280: suffix: t4_b_aij_bs3
281: nsize: 6
282: args: -test_id 1 -mat_type aij -block_size 3
284: test:
285: suffix: t5_b_baij_bs3
286: nsize: 6
287: args: -test_id 1 -mat_type baij -block_size 3
288: filter: grep -v Mat_
290: test:
291: suffix: t6_b_sbaij_bs3
292: nsize: 6
293: args: -test_id 1 -mat_type sbaij -block_size 3
294: filter: grep -v Mat_
296: TEST*/