Actual source code: ex226.c
petsc-3.13.6 2020-09-29
1: static char help[] = "Benchmark for MatMatMult() of AIJ matrices using different 2d finite-difference stencils.\n\n";
2:
3: #include <petscmat.h>
5: /* Converts 3d grid coordinates (i,j,k) for a grid of size m \times n to global indexing. Pass k = 0 for a 2d grid. */
6: int global_index(PetscInt i,PetscInt j,PetscInt k, PetscInt m, PetscInt n) { return i + j * m + k * m * n; }
8: int main(int argc,char **argv)
9: {
10: Mat A,B,C,PtAP,PtAP_copy,PtAP_squared;
11: PetscInt i,M,N,Istart,Iend,n=7,j,J,Ii,m=8,k,o=1;
12: PetscScalar v;
14: PetscBool equal=PETSC_FALSE,mat_view=PETSC_FALSE;
15: char stencil[PETSC_MAX_PATH_LEN];
16: #if defined(PETSC_USE_LOG)
17: PetscLogStage fullMatMatMultStage;
18: #endif
20: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
21: PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
22: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
23: PetscOptionsGetInt(NULL,NULL,"-o",&o,NULL);
24: PetscOptionsHasName(NULL,NULL,"-result_view",&mat_view);
25: PetscOptionsGetString(NULL,NULL,"-stencil",stencil,PETSC_MAX_PATH_LEN,NULL);
27: /* Create a aij matrix A */
28: M = N = m*n*o;
29: MatCreate(PETSC_COMM_WORLD,&A);
30: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);
31: MatSetType(A,MATAIJ);
32: MatSetFromOptions(A);
34: /* Consistency checks */
35: if (o < 1 || m < 1 || n < 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER,"Dimensions need to be larger than zero!");
37: /************ 2D stencils ***************/
38: PetscStrcmp(stencil, "2d5point", &equal);
39: if (equal) { /* 5-point stencil, 2D */
40: MatMPIAIJSetPreallocation(A,5,NULL,5,NULL);
41: MatSeqAIJSetPreallocation(A,5,NULL);
42: MatGetOwnershipRange(A,&Istart,&Iend);
43: for (Ii=Istart; Ii<Iend; Ii++) {
44: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
45: if (i>0) {J = global_index(i-1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
46: if (i<m-1) {J = global_index(i+1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
47: if (j>0) {J = global_index(i,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
48: if (j<n-1) {J = global_index(i,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
49: v = 4.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
50: }
51: }
52: PetscStrcmp(stencil, "2d9point", &equal);
53: if (equal) { /* 9-point stencil, 2D */
54: MatMPIAIJSetPreallocation(A,9,NULL,9,NULL);
55: MatSeqAIJSetPreallocation(A,9,NULL);
56: MatGetOwnershipRange(A,&Istart,&Iend);
57: for (Ii=Istart; Ii<Iend; Ii++) {
58: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
59: if (i>0) {J = global_index(i-1,j, k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
60: if (i>0 && j>0) {J = global_index(i-1,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
61: if ( j>0) {J = global_index(i, j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
62: if (i<m-1 && j>0) {J = global_index(i+1,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
63: if (i<m-1) {J = global_index(i+1,j, k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
64: if (i<m-1 && j<n-1) {J = global_index(i+1,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
65: if (j<n-1) {J = global_index(i, j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
66: if (i>0 && j<n-1) {J = global_index(i-1,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
67: v = 8.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
68: }
69: }
70: PetscStrcmp(stencil, "2d9point2", &equal);
71: if (equal) { /* 9-point Cartesian stencil (width 2 per coordinate), 2D */
72: MatMPIAIJSetPreallocation(A,9,NULL,9,NULL);
73: MatSeqAIJSetPreallocation(A,9,NULL);
74: MatGetOwnershipRange(A,&Istart,&Iend);
75: for (Ii=Istart; Ii<Iend; Ii++) {
76: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
77: if (i>0) {J = global_index(i-1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
78: if (i>1) {J = global_index(i-2,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
79: if (i<m-1) {J = global_index(i+1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
80: if (i<m-2) {J = global_index(i+2,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
81: if (j>0) {J = global_index(i,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
82: if (j>1) {J = global_index(i,j-2,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
83: if (j<n-1) {J = global_index(i,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
84: if (j<n-2) {J = global_index(i,j+2,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
85: v = 8.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
86: }
87: }
88: PetscStrcmp(stencil, "2d13point", &equal);
89: if (equal) { /* 13-point Cartesian stencil (width 3 per coordinate), 2D */
90: MatMPIAIJSetPreallocation(A,13,NULL,13,NULL);
91: MatSeqAIJSetPreallocation(A,13,NULL);
92: MatGetOwnershipRange(A,&Istart,&Iend);
93: for (Ii=Istart; Ii<Iend; Ii++) {
94: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
95: if (i>0) {J = global_index(i-1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
96: if (i>1) {J = global_index(i-2,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
97: if (i>2) {J = global_index(i-3,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
98: if (i<m-1) {J = global_index(i+1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
99: if (i<m-2) {J = global_index(i+2,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
100: if (i<m-3) {J = global_index(i+3,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
101: if (j>0) {J = global_index(i,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
102: if (j>1) {J = global_index(i,j-2,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
103: if (j>2) {J = global_index(i,j-3,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
104: if (j<n-1) {J = global_index(i,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
105: if (j<n-2) {J = global_index(i,j+2,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
106: if (j<n-3) {J = global_index(i,j+3,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
107: v = 12.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
108: }
109: }
110: /************ 3D stencils ***************/
111: PetscStrcmp(stencil, "3d7point", &equal);
112: if (equal) { /* 7-point stencil, 3D */
113: MatMPIAIJSetPreallocation(A,7,NULL,7,NULL);
114: MatSeqAIJSetPreallocation(A,7,NULL);
115: MatGetOwnershipRange(A,&Istart,&Iend);
116: for (Ii=Istart; Ii<Iend; Ii++) {
117: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
118: if (i>0) {J = global_index(i-1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
119: if (i<m-1) {J = global_index(i+1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
120: if (j>0) {J = global_index(i,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
121: if (j<n-1) {J = global_index(i,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
122: if (k>0) {J = global_index(i,j,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
123: if (k<o-1) {J = global_index(i,j,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
124: v = 6.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
125: }
126: }
127: PetscStrcmp(stencil, "3d13point", &equal);
128: if (equal) { /* 13-point stencil, 3D */
129: MatMPIAIJSetPreallocation(A,13,NULL,13,NULL);
130: MatSeqAIJSetPreallocation(A,13,NULL);
131: MatGetOwnershipRange(A,&Istart,&Iend);
132: for (Ii=Istart; Ii<Iend; Ii++) {
133: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
134: if (i>0) {J = global_index(i-1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
135: if (i>1) {J = global_index(i-2,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
136: if (i<m-1) {J = global_index(i+1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
137: if (i<m-2) {J = global_index(i+2,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
138: if (j>0) {J = global_index(i,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
139: if (j>1) {J = global_index(i,j-2,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
140: if (j<n-1) {J = global_index(i,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
141: if (j<n-2) {J = global_index(i,j+2,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
142: if (k>0) {J = global_index(i,j,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
143: if (k>1) {J = global_index(i,j,k-2,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
144: if (k<o-1) {J = global_index(i,j,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
145: if (k<o-2) {J = global_index(i,j,k+2,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
146: v = 12.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
147: }
148: }
149: PetscStrcmp(stencil, "3d19point", &equal);
150: if (equal) { /* 19-point stencil, 3D */
151: MatMPIAIJSetPreallocation(A,19,NULL,19,NULL);
152: MatSeqAIJSetPreallocation(A,19,NULL);
153: MatGetOwnershipRange(A,&Istart,&Iend);
154: for (Ii=Istart; Ii<Iend; Ii++) {
155: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
156: /* one hop */
157: if (i>0) {J = global_index(i-1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
158: if (i<m-1) {J = global_index(i+1,j,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
159: if (j>0) {J = global_index(i,j-1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
160: if (j<n-1) {J = global_index(i,j+1,k,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
161: if (k>0) {J = global_index(i,j,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
162: if (k<o-1) {J = global_index(i,j,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
163: /* two hops */
164: if (i>0 && j>0) {J = global_index(i-1,j-1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
165: if (i>0 && k>0) {J = global_index(i-1,j, k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
166: if (i>0 && j<n-1) {J = global_index(i-1,j+1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
167: if (i>0 && k<o-1) {J = global_index(i-1,j, k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
168: if (i<m-1 && j>0) {J = global_index(i+1,j-1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
169: if (i<m-1 && k>0) {J = global_index(i+1,j, k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
170: if (i<m-1 && j<n-1) {J = global_index(i+1,j+1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
171: if (i<m-1 && k<o-1) {J = global_index(i+1,j, k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
172: if (j>0 && k>0) {J = global_index(i, j-1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
173: if (j>0 && k<o-1) {J = global_index(i, j-1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
174: if (j<n-1 && k>0) {J = global_index(i, j+1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
175: if (j<n-1 && k<o-1) {J = global_index(i, j+1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
176: v = 18.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
177: }
178: }
179: PetscStrcmp(stencil, "3d27point", &equal);
180: if (equal) { /* 27-point stencil, 3D */
181: MatMPIAIJSetPreallocation(A,27,NULL,27,NULL);
182: MatSeqAIJSetPreallocation(A,27,NULL);
183: MatGetOwnershipRange(A,&Istart,&Iend);
184: for (Ii=Istart; Ii<Iend; Ii++) {
185: v = -1.0; k = Ii / (m*n); j = (Ii - k * m * n) / m; i = (Ii - k * m * n - j * m);
186: if (k>0) {
187: if (j>0) {
188: if (i>0) {J = global_index(i-1,j-1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
189: J = global_index(i, j-1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
190: if (i<m-1) {J = global_index(i+1,j-1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
191: }
192: {
193: if (i>0) {J = global_index(i-1,j, k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
194: J = global_index(i, j, k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
195: if (i<m-1) {J = global_index(i+1,j, k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
196: }
197: if (j<n-1) {
198: if (i>0) {J = global_index(i-1,j+1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
199: J = global_index(i, j+1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
200: if (i<m-1) {J = global_index(i+1,j+1,k-1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
201: }
202: }
203: {
204: if (j>0) {
205: if (i>0) {J = global_index(i-1,j-1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
206: J = global_index(i, j-1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
207: if (i<m-1) {J = global_index(i+1,j-1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
208: }
209: {
210: if (i>0) {J = global_index(i-1,j, k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
211: J = global_index(i, j, k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
212: if (i<m-1) {J = global_index(i+1,j, k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
213: }
214: if (j<n-1) {
215: if (i>0) {J = global_index(i-1,j+1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
216: J = global_index(i, j+1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
217: if (i<m-1) {J = global_index(i+1,j+1,k ,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
218: }
219: }
220: if (k<o-1) {
221: if (j>0) {
222: if (i>0) {J = global_index(i-1,j-1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
223: J = global_index(i, j-1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
224: if (i<m-1) {J = global_index(i+1,j-1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
225: }
226: {
227: if (i>0) {J = global_index(i-1,j, k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
228: J = global_index(i, j, k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
229: if (i<m-1) {J = global_index(i+1,j, k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
230: }
231: if (j<n-1) {
232: if (i>0) {J = global_index(i-1,j+1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
233: J = global_index(i, j+1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
234: if (i<m-1) {J = global_index(i+1,j+1,k+1,m,n); MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
235: }
236: }
237: v = 26.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
238: }
239: }
240: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
241: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
243: /* Copy A into B in order to have a more representative benchmark (A*A has more cache hits than A*B) */
244: MatDuplicate(A,MAT_COPY_VALUES,&B);
246: PetscLogStageRegister("Full MatMatMult",&fullMatMatMultStage);
248: /* Test C = A*B */
249: PetscLogStagePush(fullMatMatMultStage);
250: MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C);
252: /* Test PtAP_squared = PtAP(C,C)*PtAP(C,C) */
253: MatPtAP(C,C,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&PtAP);
254: MatDuplicate(PtAP,MAT_COPY_VALUES,&PtAP_copy);
255: MatMatMult(PtAP,PtAP_copy,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&PtAP_squared);
257: MatView(C,PETSC_VIEWER_STDOUT_WORLD);
258: MatView(PtAP_squared,PETSC_VIEWER_STDOUT_WORLD);
261: MatDestroy(&PtAP_squared);
262: MatDestroy(&PtAP_copy);
263: MatDestroy(&PtAP);
264: MatDestroy(&C);
265: MatDestroy(&B);
266: MatDestroy(&A);
267: PetscFinalize();
268: return ierr;
269: }
272: /*TEST
274: test:
275: suffix: 1
276: nsize: 1
277: args: -m 8 -n 8 -stencil 2d5point -matmatmult_via sorted
279: test:
280: suffix: 2
281: nsize: 1
282: args: -m 5 -n 5 -o 5 -stencil 3d27point -matmatmult_via rowmerge
284: test:
285: suffix: 3
286: nsize: 4
287: args: -m 6 -n 6 -stencil 2d5point -matmatmult_via seqmpi
291: TEST*/