Actual source code: ex230.c
petsc-3.14.6 2021-03-30
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;
17: /*
18: Create the Jacobian matrix
19: */
20: M = 10;
21: N = 8;
22: MatCreate(PETSC_COMM_WORLD,&A);
23: MatSetType(A,MATAIJ);
24: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);
25: MatSetBlockSize(A,1);
26: MatSetFromOptions(A);
28: /*
29: Get the sizes of the jacobian.
30: */
31: MatGetLocalSize(A,&m,&n);
32: MatGetBlockSize(A,&bs);
34: /*
35: Create a preallocator matrix with sizes (local and global) matching the jacobian A.
36: */
37: MatCreate(PetscObjectComm((PetscObject)A),&preallocator);
38: MatSetType(preallocator,MATPREALLOCATOR);
39: MatSetSizes(preallocator,m,n,M,N);
40: MatSetBlockSize(preallocator,bs);
41: MatSetUp(preallocator);
43: /*
44: Insert non-zero pattern (e.g. perform a sweep over the grid).
45: You can use MatSetValues(), MatSetValuesBlocked() or MatSetValue().
46: */
47: {
48: PetscInt ii,jj;
49: PetscScalar vv = 0.0;
51: ii = 3; jj = 3;
52: MatSetValues(preallocator,1,&ii,1,&jj,&vv,INSERT_VALUES);
54: ii = 7; jj = 4;
55: MatSetValues(preallocator,1,&ii,1,&jj,&vv,INSERT_VALUES);
57: ii = 9; jj = 7;
58: MatSetValue(preallocator,ii,jj,vv,INSERT_VALUES);
59: }
60: MatAssemblyBegin(preallocator,MAT_FINAL_ASSEMBLY);
61: MatAssemblyEnd(preallocator,MAT_FINAL_ASSEMBLY);
63: /*
64: Push the non-zero pattern defined within preallocator into A.
65: Internally this will call MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE).
66: The arg fill = PETSC_TRUE will insert zeros in the matrix A automatically.
67: */
68: MatPreallocatorPreallocate(preallocator,PETSC_TRUE,A);
70: /*
71: We no longer require the preallocator object so destroy it.
72: */
73: MatDestroy(&preallocator);
75: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
77: /*
78: Insert non-zero values into A.
79: */
80: {
81: PetscInt ii,jj;
82: PetscScalar vv;
84: ii = 3; jj = 3; vv = 0.3;
85: MatSetValue(A,ii,jj,vv,INSERT_VALUES);
87: ii = 7; jj = 4; vv = 3.3;
88: MatSetValues(A,1,&ii,1,&jj,&vv,INSERT_VALUES);
90: ii = 9; jj = 7; vv = 4.3;
91: MatSetValues(A,1,&ii,1,&jj,&vv,INSERT_VALUES);
92: }
93: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
94: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
96: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
98: MatDestroy(&A);
99: return(0);
100: }
102: PetscErrorCode ex2_square_bsvariable(void)
103: {
104: Mat A,preallocator;
105: PetscInt M,N,m,n,bs = 1;
109: PetscOptionsGetInt(NULL,NULL,"-block_size",&bs,NULL);
111: /*
112: Create the Jacobian matrix.
113: */
114: M = 10 * bs;
115: N = 10 * bs;
116: MatCreate(PETSC_COMM_WORLD,&A);
117: MatSetType(A,MATAIJ);
118: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);
119: MatSetBlockSize(A,bs);
120: MatSetFromOptions(A);
122: /*
123: Get the sizes of the jacobian.
124: */
125: MatGetLocalSize(A,&m,&n);
126: MatGetBlockSize(A,&bs);
128: /*
129: Create a preallocator matrix with dimensions matching the jacobian A.
130: */
131: MatCreate(PetscObjectComm((PetscObject)A),&preallocator);
132: MatSetType(preallocator,MATPREALLOCATOR);
133: MatSetSizes(preallocator,m,n,M,N);
134: MatSetBlockSize(preallocator,bs);
135: MatSetUp(preallocator);
137: /*
138: Insert non-zero pattern (e.g. perform a sweep over the grid).
139: You can use MatSetValues(), MatSetValuesBlocked() or MatSetValue().
140: */
141: {
142: PetscInt ii,jj;
143: PetscScalar *vv;
145: PetscCalloc1(bs*bs,&vv);
147: ii = 0; jj = 9;
148: MatSetValue(preallocator,ii,jj,vv[0],INSERT_VALUES);
150: ii = 0; jj = 0;
151: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
153: ii = 2; jj = 4;
154: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
156: ii = 4; jj = 2;
157: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
159: ii = 4; jj = 4;
160: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
162: ii = 9; jj = 9;
163: MatSetValuesBlocked(preallocator,1,&ii,1,&jj,vv,INSERT_VALUES);
165: PetscFree(vv);
166: }
167: MatAssemblyBegin(preallocator,MAT_FINAL_ASSEMBLY);
168: MatAssemblyEnd(preallocator,MAT_FINAL_ASSEMBLY);
170: /*
171: Push non-zero pattern defined from preallocator into A.
172: Internally this will call MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE).
173: The arg fill = PETSC_TRUE will insert zeros in the matrix A automatically.
174: */
175: MatPreallocatorPreallocate(preallocator,PETSC_TRUE,A);
177: /*
178: We no longer require the preallocator object so destroy it.
179: */
180: MatDestroy(&preallocator);
182: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
184: {
185: PetscInt ii,jj;
186: PetscScalar *vv;
188: PetscCalloc1(bs*bs,&vv);
189: for (ii=0; ii<bs*bs; ii++) {
190: vv[ii] = (PetscReal)(ii+1);
191: }
193: ii = 0; jj = 9;
194: MatSetValue(A,ii,jj,33.3,INSERT_VALUES);
196: ii = 0; jj = 0;
197: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
199: ii = 2; jj = 4;
200: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
202: ii = 4; jj = 2;
203: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
205: ii = 4; jj = 4;
206: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
208: ii = 9; jj = 9;
209: MatSetValuesBlocked(A,1,&ii,1,&jj,vv,INSERT_VALUES);
211: PetscFree(vv);
212: }
213: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
214: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
216: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
218: MatDestroy(&A);
219: return(0);
220: }
222: int main(int argc, char **args)
223: {
225: PetscInt testid = 0;
226: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
227: PetscOptionsGetInt(NULL,NULL,"-test_id",&testid,NULL);
228: switch (testid) {
229: case 0:
230: ex1_nonsquare_bs1();
231: break;
232: case 1:
233: ex2_square_bsvariable();
234: break;
235: default:
236: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Invalid value for -test_id. Must be {0,1}");
237: break;
238: }
239: PetscFinalize();
240: return ierr;
241: }
243: /*TEST
245: test:
246: suffix: t0_a_aij
247: nsize: 1
248: args: -test_id 0 -mat_type aij
250: test:
251: suffix: t0_b_aij
252: nsize: 6
253: args: -test_id 0 -mat_type aij
256: test:
257: suffix: t1_a_aij
258: nsize: 1
259: args: -test_id 1 -mat_type aij
261: test:
262: suffix: t2_a_baij
263: nsize: 1
264: args: -test_id 1 -mat_type baij
266: test:
267: suffix: t3_a_sbaij
268: nsize: 1
269: args: -test_id 1 -mat_type sbaij
273: test:
274: suffix: t4_a_aij_bs3
275: nsize: 1
276: args: -test_id 1 -mat_type aij -block_size 3
278: test:
279: suffix: t5_a_baij_bs3
280: nsize: 1
281: args: -test_id 1 -mat_type baij -block_size 3
283: test:
284: suffix: t6_a_sbaij_bs3
285: nsize: 1
286: args: -test_id 1 -mat_type sbaij -block_size 3
290: test:
291: suffix: t4_b_aij_bs3
292: nsize: 6
293: args: -test_id 1 -mat_type aij -block_size 3
295: test:
296: suffix: t5_b_baij_bs3
297: nsize: 6
298: args: -test_id 1 -mat_type baij -block_size 3
299: filter: grep -v Mat_
301: test:
302: suffix: t6_b_sbaij_bs3
303: nsize: 6
304: args: -test_id 1 -mat_type sbaij -block_size 3
305: filter: grep -v Mat_
307: TEST*/