xref: /petsc/src/mat/tests/ex23.c (revision 0be0d8bde7cbed9986facbcd3689d3188b7894db)
1 static char help[] = "Tests the use of interface functions for MATIS matrices and conversion routines.\n";
2 
3 #include <petscmat.h>
4 
5 PetscErrorCode TestMatZeroRows(Mat, Mat, PetscBool, IS, PetscScalar);
6 PetscErrorCode CheckMat(Mat, Mat, PetscBool, const char *);
7 
8 int main(int argc, char **args)
9 {
10   Mat                    A, B, A2, B2, T;
11   Mat                    Aee, Aeo, Aoe, Aoo;
12   Mat                   *mats, *Asub, *Bsub;
13   Vec                    x, y;
14   MatInfo                info;
15   ISLocalToGlobalMapping cmap, rmap;
16   IS                     is, is2, reven, rodd, ceven, codd;
17   IS                    *rows, *cols;
18   IS                     irow[2], icol[2];
19   PetscLayout            rlayout, clayout;
20   const PetscInt        *rrange, *crange;
21   MatType                lmtype;
22   PetscScalar            diag = 2.;
23   PetscInt               n, m, i, lm, ln;
24   PetscInt               rst, ren, cst, cen, nr, nc;
25   PetscMPIInt            rank, size, lrank, rrank;
26   PetscBool              testT, squaretest, isaij;
27   PetscBool              permute = PETSC_FALSE, negmap = PETSC_FALSE, repmap = PETSC_FALSE, allow_repeated = PETSC_TRUE;
28   PetscBool              diffmap = PETSC_TRUE, symmetric = PETSC_FALSE, issymmetric, test_matlab = PETSC_FALSE;
29 
30   PetscFunctionBeginUser;
31   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
32   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
33   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
34   m = n = 2 * size;
35   PetscCall(PetscOptionsGetBool(NULL, NULL, "-symmetric", &symmetric, NULL));
36   PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
37   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
38   PetscCall(PetscOptionsGetBool(NULL, NULL, "-negmap", &negmap, NULL));
39   PetscCall(PetscOptionsGetBool(NULL, NULL, "-repmap", &repmap, NULL));
40   PetscCall(PetscOptionsGetBool(NULL, NULL, "-permmap", &permute, NULL));
41   PetscCall(PetscOptionsGetBool(NULL, NULL, "-diffmap", &diffmap, NULL));
42   PetscCall(PetscOptionsGetBool(NULL, NULL, "-allow_repeated", &allow_repeated, NULL));
43   PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_matlab", &test_matlab, NULL));
44   PetscCheck(size == 1 || m >= 4, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Number of rows should be larger or equal 4 for parallel runs");
45   PetscCheck(size != 1 || m >= 2, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Number of rows should be larger or equal 2 for uniprocessor runs");
46   PetscCheck(n >= 2, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Number of cols should be larger or equal 2");
47 
48   /* create a MATIS matrix */
49   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
50   PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m, n));
51   PetscCall(MatSetType(A, MATIS));
52   PetscCall(MatSetFromOptions(A));
53   if (!negmap && !repmap) {
54     /* This is not the proper setting for MATIS for finite elements, it is just used to test the routines
55        Here we use a one-to-one correspondence between local row/column spaces and global row/column spaces
56        Equivalent to passing NULL for the mapping */
57     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, 0, 1, &is));
58   } else if (negmap && !repmap) { /* non repeated but with negative indices */
59     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n + 2, -2, 1, &is));
60   } else if (!negmap && repmap) { /* non negative but repeated indices */
61     IS isl[2];
62 
63     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, 0, 1, &isl[0]));
64     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, n - 1, -1, &isl[1]));
65     PetscCall(ISConcatenate(PETSC_COMM_WORLD, 2, isl, &is));
66     PetscCall(ISDestroy(&isl[0]));
67     PetscCall(ISDestroy(&isl[1]));
68   } else { /* negative and repeated indices */
69     IS isl[2];
70 
71     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n + 1, -1, 1, &isl[0]));
72     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n + 1, n - 1, -1, &isl[1]));
73     PetscCall(ISConcatenate(PETSC_COMM_WORLD, 2, isl, &is));
74     PetscCall(ISDestroy(&isl[0]));
75     PetscCall(ISDestroy(&isl[1]));
76   }
77   PetscCall(ISLocalToGlobalMappingCreateIS(is, &cmap));
78   PetscCall(ISDestroy(&is));
79 
80   if (m != n || diffmap) {
81     PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, permute ? m - 1 : 0, permute ? -1 : 1, &is));
82     PetscCall(ISLocalToGlobalMappingCreateIS(is, &rmap));
83     PetscCall(ISDestroy(&is));
84   } else {
85     PetscCall(PetscObjectReference((PetscObject)cmap));
86     rmap = cmap;
87   }
88 
89   PetscCall(MatISSetAllowRepeated(A, allow_repeated));
90   PetscCall(MatSetLocalToGlobalMapping(A, rmap, cmap));
91   PetscCall(MatISStoreL2L(A, PETSC_FALSE));
92   PetscCall(MatISSetPreallocation(A, 3, NULL, 3, NULL));
93   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, (PetscBool) !(repmap || negmap))); /* I do not want to precompute the pattern */
94   PetscCall(ISLocalToGlobalMappingGetSize(rmap, &lm));
95   PetscCall(ISLocalToGlobalMappingGetSize(cmap, &ln));
96   for (i = 0; i < lm; i++) {
97     PetscScalar v[3];
98     PetscInt    cols[3];
99 
100     cols[0] = (i - 1 + n) % n;
101     cols[1] = i % n;
102     cols[2] = (i + 1) % n;
103     v[0]    = -1. * (symmetric ? PetscMin(i + 1, cols[0] + 1) : i + 1);
104     v[1]    = 2. * (symmetric ? PetscMin(i + 1, cols[1] + 1) : i + 1);
105     v[2]    = -1. * (symmetric ? PetscMin(i + 1, cols[2] + 1) : i + 1);
106     PetscCall(ISGlobalToLocalMappingApply(cmap, IS_GTOLM_MASK, 3, cols, NULL, cols));
107     PetscCall(MatSetValuesLocal(A, 1, &i, 3, cols, v, ADD_VALUES));
108   }
109   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
110   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
111 
112   /* activate tests for square matrices with same maps only */
113   PetscCall(MatHasCongruentLayouts(A, &squaretest));
114   if (squaretest && rmap != cmap) {
115     PetscInt nr, nc;
116 
117     PetscCall(ISLocalToGlobalMappingGetSize(rmap, &nr));
118     PetscCall(ISLocalToGlobalMappingGetSize(cmap, &nc));
119     if (nr != nc) squaretest = PETSC_FALSE;
120     else {
121       const PetscInt *idxs1, *idxs2;
122 
123       PetscCall(ISLocalToGlobalMappingGetIndices(rmap, &idxs1));
124       PetscCall(ISLocalToGlobalMappingGetIndices(cmap, &idxs2));
125       PetscCall(PetscArraycmp(idxs1, idxs2, nr, &squaretest));
126       PetscCall(ISLocalToGlobalMappingRestoreIndices(rmap, &idxs1));
127       PetscCall(ISLocalToGlobalMappingRestoreIndices(cmap, &idxs2));
128     }
129     PetscCall(MPIU_Allreduce(MPI_IN_PLACE, &squaretest, 1, MPIU_BOOL, MPI_LAND, PetscObjectComm((PetscObject)A)));
130   }
131   if (negmap && repmap) squaretest = PETSC_FALSE;
132 
133   /* test MatISGetLocalMat */
134   PetscCall(MatISGetLocalMat(A, &B));
135   PetscCall(MatGetType(B, &lmtype));
136 
137   /* test MatGetInfo */
138   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatGetInfo\n"));
139   PetscCall(MatGetInfo(A, MAT_LOCAL, &info));
140   PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD));
141   PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "Process  %2d: %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", PetscGlobalRank, (PetscInt)info.nz_used, (PetscInt)info.nz_allocated,
142                                                (PetscInt)info.nz_unneeded, (PetscInt)info.assemblies, (PetscInt)info.mallocs));
143   PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
144   PetscCall(MatGetInfo(A, MAT_GLOBAL_MAX, &info));
145   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "GlobalMax  : %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", (PetscInt)info.nz_used, (PetscInt)info.nz_allocated, (PetscInt)info.nz_unneeded,
146                                    (PetscInt)info.assemblies, (PetscInt)info.mallocs));
147   PetscCall(MatGetInfo(A, MAT_GLOBAL_SUM, &info));
148   PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "GlobalSum  : %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", (PetscInt)info.nz_used, (PetscInt)info.nz_allocated, (PetscInt)info.nz_unneeded,
149                                    (PetscInt)info.assemblies, (PetscInt)info.mallocs));
150 
151   /* test MatIsSymmetric */
152   PetscCall(MatIsSymmetric(A, 0.0, &issymmetric));
153   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatIsSymmetric: %d\n", issymmetric));
154 
155   /* Create a MPIAIJ matrix, same as A */
156   PetscCall(MatCreate(PETSC_COMM_WORLD, &B));
157   PetscCall(MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, m, n));
158   PetscCall(MatSetType(B, MATAIJ));
159   PetscCall(MatSetFromOptions(B));
160   PetscCall(MatSetLocalToGlobalMapping(B, rmap, cmap));
161   PetscCall(MatMPIAIJSetPreallocation(B, 3, NULL, 3, NULL));
162   PetscCall(MatMPIBAIJSetPreallocation(B, 1, 3, NULL, 3, NULL));
163 #if defined(PETSC_HAVE_HYPRE)
164   PetscCall(MatHYPRESetPreallocation(B, 3, NULL, 3, NULL));
165 #endif
166   PetscCall(MatISSetPreallocation(B, 3, NULL, 3, NULL));
167   PetscCall(MatSetOption(B, MAT_NEW_NONZERO_ALLOCATION_ERR, (PetscBool) !(repmap || negmap))); /* I do not want to precompute the pattern */
168   for (i = 0; i < lm; i++) {
169     PetscScalar v[3];
170     PetscInt    cols[3];
171 
172     cols[0] = (i - 1 + n) % n;
173     cols[1] = i % n;
174     cols[2] = (i + 1) % n;
175     v[0]    = -1. * (symmetric ? PetscMin(i + 1, cols[0] + 1) : i + 1);
176     v[1]    = 2. * (symmetric ? PetscMin(i + 1, cols[1] + 1) : i + 1);
177     v[2]    = -1. * (symmetric ? PetscMin(i + 1, cols[2] + 1) : i + 1);
178     PetscCall(ISGlobalToLocalMappingApply(cmap, IS_GTOLM_MASK, 3, cols, NULL, cols));
179     PetscCall(MatSetValuesLocal(B, 1, &i, 3, cols, v, ADD_VALUES));
180   }
181   PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
182   PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
183 
184   /* test MatView */
185   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatView\n"));
186   PetscCall(MatView(A, NULL));
187   PetscCall(MatView(B, NULL));
188 
189   /* test MATLAB ASCII view */
190   if (test_matlab) { /* output is different when using real or complex numbers */
191     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatView ASCII MATLAB\n"));
192     PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_MATLAB));
193     PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD));
194     PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
195   }
196 
197   /* test CheckMat */
198   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test CheckMat\n"));
199   PetscCall(CheckMat(A, B, PETSC_FALSE, "CheckMat"));
200 
201   /* test binary MatView/MatLoad */
202   {
203     PetscMPIInt color = rank % 2;
204     MPI_Comm    comm;
205     char        name[PETSC_MAX_PATH_LEN];
206     PetscViewer wview, cview, sview, view;
207     Mat         A2;
208 
209     PetscCallMPI(MPI_Comm_split(PETSC_COMM_WORLD, color, rank, &comm));
210 
211     PetscCall(PetscSNPrintf(name, PETSC_STATIC_ARRAY_LENGTH(name), "world_is"));
212     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_WRITE, &wview));
213     PetscCall(PetscSNPrintf(name, PETSC_STATIC_ARRAY_LENGTH(name), "seq_is_%d", rank));
214     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, name, FILE_MODE_WRITE, &sview));
215     PetscCall(PetscSNPrintf(name, PETSC_STATIC_ARRAY_LENGTH(name), "color_is_%d", color));
216     PetscCall(PetscViewerBinaryOpen(comm, name, FILE_MODE_WRITE, &cview));
217     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatView on binary world\n"));
218     PetscCall(MatView(A, wview));
219     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatView on binary self\n"));
220     PetscCall(MatView(A, sview));
221     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatView on binary subcomm\n"));
222     PetscCall(MatView(A, cview));
223     PetscCall(PetscViewerDestroy(&wview));
224     PetscCall(PetscViewerDestroy(&cview));
225     PetscCall(PetscViewerDestroy(&sview));
226 
227     /* Load a world matrix */
228     PetscCall(MatCreate(PETSC_COMM_WORLD, &A2));
229     PetscCall(MatSetType(A2, MATIS));
230     PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_INFO_DETAIL));
231 
232     /* Read back the same matrix and check */
233     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatLoad from world\n"));
234     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "world_is", FILE_MODE_READ, &view));
235     PetscCall(MatLoad(A2, view));
236     PetscCall(CheckMat(A, A2, PETSC_TRUE, "Load"));
237     PetscCall(MatView(A2, PETSC_VIEWER_STDOUT_WORLD));
238     PetscCall(PetscViewerDestroy(&view));
239 
240     /* Read the matrix from rank 0 only */
241     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatLoad from self\n"));
242     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "seq_is_0", FILE_MODE_READ, &view));
243     PetscCall(MatLoad(A2, view));
244     PetscCall(MatView(A2, PETSC_VIEWER_STDOUT_WORLD));
245     PetscCall(PetscViewerDestroy(&view));
246 
247     /* Read the matrix from subcomm */
248     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatLoad from subcomm\n"));
249     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, "color_is_0", FILE_MODE_READ, &view));
250     PetscCall(MatLoad(A2, view));
251     PetscCall(MatView(A2, PETSC_VIEWER_STDOUT_WORLD));
252     PetscCall(PetscViewerDestroy(&view));
253 
254     PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
255     PetscCall(MatDestroy(&A2));
256 
257     /* now load the original matrix from color 0 only processes */
258     if (!color) {
259       PetscCall(PetscPrintf(comm, "Test subcomm MatLoad from world\n"));
260       PetscCall(MatCreate(comm, &A2));
261       PetscCall(MatSetType(A2, MATIS));
262       PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(comm), PETSC_VIEWER_ASCII_INFO_DETAIL));
263       PetscCall(PetscViewerBinaryOpen(comm, "world_is", FILE_MODE_READ, &view));
264       PetscCall(MatLoad(A2, view));
265       PetscCall(MatView(A2, PETSC_VIEWER_STDOUT_(comm)));
266       PetscCall(PetscViewerDestroy(&view));
267       PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(comm)));
268       PetscCall(MatDestroy(&A2));
269     }
270 
271     PetscCallMPI(MPI_Comm_free(&comm));
272   }
273 
274   /* test MatDuplicate and MatAXPY */
275   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatDuplicate and MatAXPY\n"));
276   PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
277   PetscCall(CheckMat(A, A2, PETSC_FALSE, "MatDuplicate and MatAXPY"));
278 
279   /* test MatConvert */
280   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatConvert_IS_XAIJ\n"));
281   PetscCall(MatConvert(A2, MATAIJ, MAT_INITIAL_MATRIX, &B2));
282   PetscCall(CheckMat(B, B2, PETSC_TRUE, "MatConvert_IS_XAIJ MAT_INITIAL_MATRIX"));
283   PetscCall(MatConvert(A2, MATAIJ, MAT_REUSE_MATRIX, &B2));
284   PetscCall(CheckMat(B, B2, PETSC_TRUE, "MatConvert_IS_XAIJ MAT_REUSE_MATRIX"));
285   PetscCall(MatConvert(A2, MATAIJ, MAT_INPLACE_MATRIX, &A2));
286   PetscCall(CheckMat(B, A2, PETSC_TRUE, "MatConvert_IS_XAIJ MAT_INPLACE_MATRIX"));
287   PetscCall(MatDestroy(&A2));
288   PetscCall(MatDestroy(&B2));
289   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatConvert_XAIJ_IS\n"));
290   PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
291   PetscCall(MatConvert(B2, MATIS, MAT_INITIAL_MATRIX, &A2));
292   PetscCall(CheckMat(A, A2, PETSC_TRUE, "MatConvert_XAIJ_IS MAT_INITIAL_MATRIX"));
293   PetscCall(MatConvert(B2, MATIS, MAT_REUSE_MATRIX, &A2));
294   PetscCall(CheckMat(A, A2, PETSC_TRUE, "MatConvert_XAIJ_IS MAT_REUSE_MATRIX"));
295   PetscCall(MatConvert(B2, MATIS, MAT_INPLACE_MATRIX, &B2));
296   PetscCall(CheckMat(A, B2, PETSC_TRUE, "MatConvert_XAIJ_IS MAT_INPLACE_MATRIX"));
297   PetscCall(MatDestroy(&A2));
298   PetscCall(MatDestroy(&B2));
299   PetscCall(PetscStrcmp(lmtype, MATSEQAIJ, &isaij));
300   if (size == 1 && isaij) { /* tests special code paths in MatConvert_IS_XAIJ */
301     PetscInt               ri, ci, rr[3] = {0, 1, 0}, cr[4] = {1, 2, 0, 1}, rk[3] = {0, 2, 1}, ck[4] = {1, 0, 3, 2};
302     ISLocalToGlobalMapping tcmap, trmap;
303 
304     for (ri = 0; ri < 2; ri++) {
305       PetscInt *r;
306 
307       r = (PetscInt *)(ri == 0 ? rr : rk);
308       for (ci = 0; ci < 2; ci++) {
309         PetscInt *c, rb, cb;
310 
311         c = (PetscInt *)(ci == 0 ? cr : ck);
312         for (rb = 1; rb < 4; rb++) {
313           PetscCall(ISCreateBlock(PETSC_COMM_SELF, rb, 3, r, PETSC_COPY_VALUES, &is));
314           PetscCall(ISLocalToGlobalMappingCreateIS(is, &trmap));
315           PetscCall(ISDestroy(&is));
316           for (cb = 1; cb < 4; cb++) {
317             Mat  T, lT, T2;
318             char testname[256];
319 
320             PetscCall(PetscSNPrintf(testname, sizeof(testname), "MatConvert_IS_XAIJ special case (%" PetscInt_FMT " %" PetscInt_FMT ", bs %" PetscInt_FMT " %" PetscInt_FMT ")", ri, ci, rb, cb));
321             PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test %s\n", testname));
322 
323             PetscCall(ISCreateBlock(PETSC_COMM_SELF, cb, 4, c, PETSC_COPY_VALUES, &is));
324             PetscCall(ISLocalToGlobalMappingCreateIS(is, &tcmap));
325             PetscCall(ISDestroy(&is));
326 
327             PetscCall(MatCreate(PETSC_COMM_SELF, &T));
328             PetscCall(MatSetSizes(T, PETSC_DECIDE, PETSC_DECIDE, rb * 3, cb * 4));
329             PetscCall(MatSetType(T, MATIS));
330             PetscCall(MatSetLocalToGlobalMapping(T, trmap, tcmap));
331             PetscCall(ISLocalToGlobalMappingDestroy(&tcmap));
332             PetscCall(MatISGetLocalMat(T, &lT));
333             PetscCall(MatSetType(lT, MATSEQAIJ));
334             PetscCall(MatSeqAIJSetPreallocation(lT, cb * 4, NULL));
335             PetscCall(MatSetRandom(lT, NULL));
336             PetscCall(MatConvert(lT, lmtype, MAT_INPLACE_MATRIX, &lT));
337             PetscCall(MatISRestoreLocalMat(T, &lT));
338             PetscCall(MatAssemblyBegin(T, MAT_FINAL_ASSEMBLY));
339             PetscCall(MatAssemblyEnd(T, MAT_FINAL_ASSEMBLY));
340 
341             PetscCall(MatConvert(T, MATAIJ, MAT_INITIAL_MATRIX, &T2));
342             PetscCall(CheckMat(T, T2, PETSC_TRUE, "MAT_INITIAL_MATRIX"));
343             PetscCall(MatConvert(T, MATAIJ, MAT_REUSE_MATRIX, &T2));
344             PetscCall(CheckMat(T, T2, PETSC_TRUE, "MAT_REUSE_MATRIX"));
345             PetscCall(MatDestroy(&T2));
346             PetscCall(MatDuplicate(T, MAT_COPY_VALUES, &T2));
347             PetscCall(MatConvert(T2, MATAIJ, MAT_INPLACE_MATRIX, &T2));
348             PetscCall(CheckMat(T, T2, PETSC_TRUE, "MAT_INPLACE_MATRIX"));
349             PetscCall(MatDestroy(&T));
350             PetscCall(MatDestroy(&T2));
351           }
352           PetscCall(ISLocalToGlobalMappingDestroy(&trmap));
353         }
354       }
355     }
356   }
357 
358   /* test MatDiagonalScale */
359   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatDiagonalScale\n"));
360   PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
361   PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
362   PetscCall(MatCreateVecs(A, &x, &y));
363   PetscCall(VecSetRandom(x, NULL));
364   if (issymmetric) {
365     PetscCall(VecCopy(x, y));
366   } else {
367     PetscCall(VecSetRandom(y, NULL));
368     PetscCall(VecScale(y, 8.));
369   }
370   PetscCall(MatDiagonalScale(A2, y, x));
371   PetscCall(MatDiagonalScale(B2, y, x));
372   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatDiagonalScale"));
373   PetscCall(MatDestroy(&A2));
374   PetscCall(MatDestroy(&B2));
375   PetscCall(VecDestroy(&x));
376   PetscCall(VecDestroy(&y));
377 
378   /* test MatPtAP (A IS and B AIJ) */
379   if (isaij && m == n) {
380     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatPtAP\n"));
381     /* There's a bug in MatCreateSubMatrices_MPIAIJ I cannot figure out */
382     if (!allow_repeated || !repmap || size == 1) {
383       PetscCall(MatISStoreL2L(A, PETSC_TRUE));
384       PetscCall(MatPtAP(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &A2));
385       PetscCall(MatPtAP(B, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &B2));
386       PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatPtAP MAT_INITIAL_MATRIX"));
387       PetscCall(MatPtAP(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &A2));
388       PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatPtAP MAT_REUSE_MATRIX"));
389       PetscCall(MatDestroy(&A2));
390       PetscCall(MatDestroy(&B2));
391     }
392   }
393 
394   /* test MatGetLocalSubMatrix */
395   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatGetLocalSubMatrix\n"));
396   PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &A2));
397   PetscCall(ISCreateStride(PETSC_COMM_SELF, lm / 2 + lm % 2, 0, 2, &reven));
398   PetscCall(ISComplement(reven, 0, lm, &rodd));
399   PetscCall(ISCreateStride(PETSC_COMM_SELF, ln / 2 + ln % 2, 0, 2, &ceven));
400   PetscCall(ISComplement(ceven, 0, ln, &codd));
401   PetscCall(MatGetLocalSubMatrix(A2, reven, ceven, &Aee));
402   PetscCall(MatGetLocalSubMatrix(A2, reven, codd, &Aeo));
403   PetscCall(MatGetLocalSubMatrix(A2, rodd, ceven, &Aoe));
404   PetscCall(MatGetLocalSubMatrix(A2, rodd, codd, &Aoo));
405   for (i = 0; i < lm; i++) {
406     PetscInt    j, je, jo, colse[3], colso[3];
407     PetscScalar ve[3], vo[3];
408     PetscScalar v[3];
409     PetscInt    cols[3];
410     PetscInt    row = i / 2;
411 
412     cols[0] = (i - 1 + n) % n;
413     cols[1] = i % n;
414     cols[2] = (i + 1) % n;
415     v[0]    = -1. * (symmetric ? PetscMin(i + 1, cols[0] + 1) : i + 1);
416     v[1]    = 2. * (symmetric ? PetscMin(i + 1, cols[1] + 1) : i + 1);
417     v[2]    = -1. * (symmetric ? PetscMin(i + 1, cols[2] + 1) : i + 1);
418     PetscCall(ISGlobalToLocalMappingApply(cmap, IS_GTOLM_MASK, 3, cols, NULL, cols));
419     for (j = 0, je = 0, jo = 0; j < 3; j++) {
420       if (cols[j] % 2) {
421         vo[jo]      = v[j];
422         colso[jo++] = cols[j] / 2;
423       } else {
424         ve[je]      = v[j];
425         colse[je++] = cols[j] / 2;
426       }
427     }
428     if (i % 2) {
429       PetscCall(MatSetValuesLocal(Aoe, 1, &row, je, colse, ve, ADD_VALUES));
430       PetscCall(MatSetValuesLocal(Aoo, 1, &row, jo, colso, vo, ADD_VALUES));
431     } else {
432       PetscCall(MatSetValuesLocal(Aee, 1, &row, je, colse, ve, ADD_VALUES));
433       PetscCall(MatSetValuesLocal(Aeo, 1, &row, jo, colso, vo, ADD_VALUES));
434     }
435   }
436   PetscCall(MatRestoreLocalSubMatrix(A2, reven, ceven, &Aee));
437   PetscCall(MatRestoreLocalSubMatrix(A2, reven, codd, &Aeo));
438   PetscCall(MatRestoreLocalSubMatrix(A2, rodd, ceven, &Aoe));
439   PetscCall(MatRestoreLocalSubMatrix(A2, rodd, codd, &Aoo));
440   PetscCall(ISDestroy(&reven));
441   PetscCall(ISDestroy(&ceven));
442   PetscCall(ISDestroy(&rodd));
443   PetscCall(ISDestroy(&codd));
444   PetscCall(MatAssemblyBegin(A2, MAT_FINAL_ASSEMBLY));
445   PetscCall(MatAssemblyEnd(A2, MAT_FINAL_ASSEMBLY));
446   PetscCall(MatAXPY(A2, -1., A, SAME_NONZERO_PATTERN));
447   PetscCall(CheckMat(A2, NULL, PETSC_FALSE, "MatGetLocalSubMatrix"));
448   PetscCall(MatDestroy(&A2));
449 
450   /* test MatConvert_Nest_IS */
451   testT = PETSC_FALSE;
452   PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_trans", &testT, NULL));
453 
454   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatConvert_Nest_IS\n"));
455   nr = 2;
456   nc = 2;
457   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nr", &nr, NULL));
458   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nc", &nc, NULL));
459   if (testT) {
460     PetscCall(MatGetOwnershipRange(A, &cst, &cen));
461     PetscCall(MatGetOwnershipRangeColumn(A, &rst, &ren));
462   } else {
463     PetscCall(MatGetOwnershipRange(A, &rst, &ren));
464     PetscCall(MatGetOwnershipRangeColumn(A, &cst, &cen));
465   }
466   PetscCall(PetscMalloc3(nr, &rows, nc, &cols, 2 * nr * nc, &mats));
467   for (i = 0; i < nr * nc; i++) {
468     if (testT) {
469       PetscCall(MatCreateTranspose(A, &mats[i]));
470       PetscCall(MatTranspose(B, MAT_INITIAL_MATRIX, &mats[i + nr * nc]));
471     } else {
472       PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &mats[i]));
473       PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &mats[i + nr * nc]));
474     }
475   }
476   for (i = 0; i < nr; i++) PetscCall(ISCreateStride(PETSC_COMM_WORLD, ren - rst, i + rst, nr, &rows[i]));
477   for (i = 0; i < nc; i++) PetscCall(ISCreateStride(PETSC_COMM_WORLD, cen - cst, i + cst, nc, &cols[i]));
478   PetscCall(MatCreateNest(PETSC_COMM_WORLD, nr, rows, nc, cols, mats, &A2));
479   PetscCall(MatCreateNest(PETSC_COMM_WORLD, nr, rows, nc, cols, mats + nr * nc, &B2));
480   for (i = 0; i < nr; i++) PetscCall(ISDestroy(&rows[i]));
481   for (i = 0; i < nc; i++) PetscCall(ISDestroy(&cols[i]));
482   for (i = 0; i < 2 * nr * nc; i++) PetscCall(MatDestroy(&mats[i]));
483   PetscCall(PetscFree3(rows, cols, mats));
484   PetscCall(MatConvert(B2, MATAIJ, MAT_INITIAL_MATRIX, &T));
485   PetscCall(MatDestroy(&B2));
486   PetscCall(MatConvert(A2, MATIS, MAT_INITIAL_MATRIX, &B2));
487   PetscCall(CheckMat(B2, T, PETSC_TRUE, "MatConvert_Nest_IS MAT_INITIAL_MATRIX"));
488   PetscCall(MatConvert(A2, MATIS, MAT_REUSE_MATRIX, &B2));
489   PetscCall(CheckMat(B2, T, PETSC_TRUE, "MatConvert_Nest_IS MAT_REUSE_MATRIX"));
490   PetscCall(MatDestroy(&B2));
491   PetscCall(MatConvert(A2, MATIS, MAT_INPLACE_MATRIX, &A2));
492   PetscCall(CheckMat(A2, T, PETSC_TRUE, "MatConvert_Nest_IS MAT_INPLACE_MATRIX"));
493   PetscCall(MatDestroy(&T));
494   PetscCall(MatDestroy(&A2));
495 
496   /* test MatCreateSubMatrix */
497   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatCreateSubMatrix\n"));
498   if (rank == 0) {
499     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 1, 1, 1, &is));
500     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 2, 0, 1, &is2));
501   } else if (rank == 1) {
502     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 1, 0, 1, &is));
503     if (n > 3) {
504       PetscCall(ISCreateStride(PETSC_COMM_WORLD, 1, 3, 1, &is2));
505     } else {
506       PetscCall(ISCreateStride(PETSC_COMM_WORLD, 0, 0, 1, &is2));
507     }
508   } else if (rank == 2 && n > 4) {
509     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 0, 0, 1, &is));
510     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n - 4, 4, 1, &is2));
511   } else {
512     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 0, 0, 1, &is));
513     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 0, 0, 1, &is2));
514   }
515   PetscCall(MatCreateSubMatrix(A, is, is, MAT_INITIAL_MATRIX, &A2));
516   PetscCall(MatCreateSubMatrix(B, is, is, MAT_INITIAL_MATRIX, &B2));
517   PetscCall(CheckMat(A2, B2, PETSC_TRUE, "first MatCreateSubMatrix"));
518 
519   PetscCall(MatCreateSubMatrix(A, is, is, MAT_REUSE_MATRIX, &A2));
520   PetscCall(MatCreateSubMatrix(B, is, is, MAT_REUSE_MATRIX, &B2));
521   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "reuse MatCreateSubMatrix"));
522   PetscCall(MatDestroy(&A2));
523   PetscCall(MatDestroy(&B2));
524 
525   if (!issymmetric) {
526     PetscCall(MatCreateSubMatrix(A, is, is2, MAT_INITIAL_MATRIX, &A2));
527     PetscCall(MatCreateSubMatrix(B, is, is2, MAT_INITIAL_MATRIX, &B2));
528     PetscCall(MatCreateSubMatrix(A, is, is2, MAT_REUSE_MATRIX, &A2));
529     PetscCall(MatCreateSubMatrix(B, is, is2, MAT_REUSE_MATRIX, &B2));
530     PetscCall(CheckMat(A2, B2, PETSC_FALSE, "second MatCreateSubMatrix"));
531   }
532 
533   PetscCall(MatDestroy(&A2));
534   PetscCall(MatDestroy(&B2));
535   PetscCall(ISDestroy(&is));
536   PetscCall(ISDestroy(&is2));
537 
538   /* test MatCreateSubMatrices */
539   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatCreateSubMatrices\n"));
540   PetscCall(MatGetLayouts(A, &rlayout, &clayout));
541   PetscCall(PetscLayoutGetRanges(rlayout, &rrange));
542   PetscCall(PetscLayoutGetRanges(clayout, &crange));
543   lrank = (size + rank - 1) % size;
544   rrank = (rank + 1) % size;
545   PetscCall(ISCreateStride(PETSC_COMM_SELF, (rrange[lrank + 1] - rrange[lrank]), rrange[lrank], 1, &irow[0]));
546   PetscCall(ISCreateStride(PETSC_COMM_SELF, (crange[rrank + 1] - crange[rrank]), crange[rrank], 1, &icol[0]));
547   PetscCall(ISCreateStride(PETSC_COMM_SELF, (rrange[rrank + 1] - rrange[rrank]), rrange[rrank], 1, &irow[1]));
548   PetscCall(ISCreateStride(PETSC_COMM_SELF, (crange[lrank + 1] - crange[lrank]), crange[lrank], 1, &icol[1]));
549   PetscCall(MatCreateSubMatrices(A, 2, irow, icol, MAT_INITIAL_MATRIX, &Asub));
550   PetscCall(MatCreateSubMatrices(B, 2, irow, icol, MAT_INITIAL_MATRIX, &Bsub));
551   PetscCall(CheckMat(Asub[0], Bsub[0], PETSC_FALSE, "MatCreateSubMatrices[0]"));
552   PetscCall(CheckMat(Asub[1], Bsub[1], PETSC_FALSE, "MatCreateSubMatrices[1]"));
553   PetscCall(MatCreateSubMatrices(A, 2, irow, icol, MAT_REUSE_MATRIX, &Asub));
554   PetscCall(MatCreateSubMatrices(B, 2, irow, icol, MAT_REUSE_MATRIX, &Bsub));
555   PetscCall(CheckMat(Asub[0], Bsub[0], PETSC_FALSE, "MatCreateSubMatrices[0]"));
556   PetscCall(CheckMat(Asub[1], Bsub[1], PETSC_FALSE, "MatCreateSubMatrices[1]"));
557   PetscCall(MatDestroySubMatrices(2, &Asub));
558   PetscCall(MatDestroySubMatrices(2, &Bsub));
559   PetscCall(ISDestroy(&irow[0]));
560   PetscCall(ISDestroy(&irow[1]));
561   PetscCall(ISDestroy(&icol[0]));
562   PetscCall(ISDestroy(&icol[1]));
563 
564   /* Create an IS required by MatZeroRows(): just rank zero provides the rows to be eliminated */
565   if (size > 1) {
566     if (rank == 0) {
567       PetscInt st, len;
568 
569       st  = (m + 1) / 2;
570       len = PetscMin(m / 2, PetscMax(m - (m + 1) / 2 - 1, 0));
571       PetscCall(ISCreateStride(PETSC_COMM_WORLD, len, st, 1, &is));
572     } else {
573       PetscCall(ISCreateStride(PETSC_COMM_WORLD, 0, 0, 1, &is));
574     }
575   } else {
576     PetscCall(ISCreateStride(PETSC_COMM_WORLD, 1, 0, 1, &is));
577   }
578 
579   if (squaretest) { /* tests for square matrices only, with same maps for rows and columns */
580     PetscInt *idx0, *idx1, n0, n1;
581     IS        Ais[2], Bis[2];
582 
583     /* test MatDiagonalSet */
584     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatDiagonalSet\n"));
585     PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
586     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
587     PetscCall(MatCreateVecs(A, NULL, &x));
588     PetscCall(VecSetRandom(x, NULL));
589     PetscCall(MatDiagonalSet(A2, x, allow_repeated ? ADD_VALUES : INSERT_VALUES));
590     PetscCall(MatDiagonalSet(B2, x, allow_repeated ? ADD_VALUES : INSERT_VALUES));
591     PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatDiagonalSet"));
592     PetscCall(VecDestroy(&x));
593     PetscCall(MatDestroy(&A2));
594     PetscCall(MatDestroy(&B2));
595 
596     /* test MatShift (MatShift_IS internally uses MatDiagonalSet_IS with ADD_VALUES) */
597     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatShift\n"));
598     PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
599     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
600     PetscCall(MatShift(A2, 2.0));
601     PetscCall(MatShift(B2, 2.0));
602     PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatShift"));
603     PetscCall(MatDestroy(&A2));
604     PetscCall(MatDestroy(&B2));
605 
606     /* nonzero diag value is supported for square matrices only */
607     PetscCall(TestMatZeroRows(A, B, PETSC_TRUE, is, diag));
608 
609     /* test MatIncreaseOverlap */
610     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatIncreaseOverlap\n"));
611     PetscCall(MatGetOwnershipRange(A, &rst, &ren));
612     n0 = (ren - rst) / 2;
613     n1 = (ren - rst) / 3;
614     PetscCall(PetscMalloc1(n0, &idx0));
615     PetscCall(PetscMalloc1(n1, &idx1));
616     for (i = 0; i < n0; i++) idx0[i] = ren - i - 1;
617     for (i = 0; i < n1; i++) idx1[i] = rst + i;
618     PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, n0, idx0, PETSC_OWN_POINTER, &Ais[0]));
619     PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, n1, idx1, PETSC_OWN_POINTER, &Ais[1]));
620     PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, n0, idx0, PETSC_COPY_VALUES, &Bis[0]));
621     PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, n1, idx1, PETSC_COPY_VALUES, &Bis[1]));
622     PetscCall(MatIncreaseOverlap(A, 2, Ais, 3));
623     PetscCall(MatIncreaseOverlap(B, 2, Bis, 3));
624     /* Non deterministic output! */
625     PetscCall(ISSort(Ais[0]));
626     PetscCall(ISSort(Ais[1]));
627     PetscCall(ISSort(Bis[0]));
628     PetscCall(ISSort(Bis[1]));
629     PetscCall(ISView(Ais[0], NULL));
630     PetscCall(ISView(Bis[0], NULL));
631     PetscCall(ISView(Ais[1], NULL));
632     PetscCall(ISView(Bis[1], NULL));
633     PetscCall(MatCreateSubMatrices(A, 2, Ais, Ais, MAT_INITIAL_MATRIX, &Asub));
634     PetscCall(MatCreateSubMatrices(B, 2, Bis, Ais, MAT_INITIAL_MATRIX, &Bsub));
635     PetscCall(CheckMat(Asub[0], Bsub[0], PETSC_FALSE, "MatIncreaseOverlap[0]"));
636     PetscCall(CheckMat(Asub[1], Bsub[1], PETSC_FALSE, "MatIncreaseOverlap[1]"));
637     PetscCall(MatDestroySubMatrices(2, &Asub));
638     PetscCall(MatDestroySubMatrices(2, &Bsub));
639     PetscCall(ISDestroy(&Ais[0]));
640     PetscCall(ISDestroy(&Ais[1]));
641     PetscCall(ISDestroy(&Bis[0]));
642     PetscCall(ISDestroy(&Bis[1]));
643   }
644   PetscCall(TestMatZeroRows(A, B, squaretest, is, 0.0));
645   PetscCall(ISDestroy(&is));
646 
647   /* test MatTranspose */
648   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatTranspose\n"));
649   PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &A2));
650   PetscCall(MatTranspose(B, MAT_INITIAL_MATRIX, &B2));
651   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "initial matrix MatTranspose"));
652 
653   PetscCall(MatTranspose(A, MAT_REUSE_MATRIX, &A2));
654   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "reuse matrix (not in place) MatTranspose"));
655   PetscCall(MatDestroy(&A2));
656 
657   PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
658   PetscCall(MatTranspose(A2, MAT_INPLACE_MATRIX, &A2));
659   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "reuse matrix (in place) MatTranspose"));
660   PetscCall(MatDestroy(&A2));
661 
662   PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &A2));
663   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "reuse matrix (different type) MatTranspose"));
664   PetscCall(MatDestroy(&A2));
665   PetscCall(MatDestroy(&B2));
666 
667   /* test MatISFixLocalEmpty */
668   if (isaij) {
669     PetscInt r[2];
670 
671     r[0] = 0;
672     r[1] = PetscMin(m, n) - 1;
673     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatISFixLocalEmpty\n"));
674     PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
675 
676     PetscCall(MatISFixLocalEmpty(A2, PETSC_TRUE));
677     PetscCall(MatAssemblyBegin(A2, MAT_FINAL_ASSEMBLY));
678     PetscCall(MatAssemblyEnd(A2, MAT_FINAL_ASSEMBLY));
679     PetscCall(CheckMat(A2, B, PETSC_FALSE, "MatISFixLocalEmpty (null)"));
680 
681     PetscCall(MatZeroRows(A2, 2, r, 0.0, NULL, NULL));
682     PetscCall(MatViewFromOptions(A2, NULL, "-fixempty_view"));
683     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
684     PetscCall(MatZeroRows(B2, 2, r, 0.0, NULL, NULL));
685     PetscCall(MatISFixLocalEmpty(A2, PETSC_TRUE));
686     PetscCall(MatAssemblyBegin(A2, MAT_FINAL_ASSEMBLY));
687     PetscCall(MatAssemblyEnd(A2, MAT_FINAL_ASSEMBLY));
688     PetscCall(MatViewFromOptions(A2, NULL, "-fixempty_view"));
689     PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatISFixLocalEmpty (rows)"));
690     PetscCall(MatDestroy(&A2));
691 
692     PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
693     PetscCall(MatZeroRows(A2, 2, r, 0.0, NULL, NULL));
694     PetscCall(MatTranspose(A2, MAT_INPLACE_MATRIX, &A2));
695     PetscCall(MatTranspose(B2, MAT_INPLACE_MATRIX, &B2));
696     PetscCall(MatViewFromOptions(A2, NULL, "-fixempty_view"));
697     PetscCall(MatISFixLocalEmpty(A2, PETSC_TRUE));
698     PetscCall(MatAssemblyBegin(A2, MAT_FINAL_ASSEMBLY));
699     PetscCall(MatAssemblyEnd(A2, MAT_FINAL_ASSEMBLY));
700     PetscCall(MatViewFromOptions(A2, NULL, "-fixempty_view"));
701     PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatISFixLocalEmpty (cols)"));
702 
703     PetscCall(MatDestroy(&A2));
704     PetscCall(MatDestroy(&B2));
705 
706     if (squaretest) {
707       PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
708       PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
709       PetscCall(MatZeroRowsColumns(A2, 2, r, 0.0, NULL, NULL));
710       PetscCall(MatZeroRowsColumns(B2, 2, r, 0.0, NULL, NULL));
711       PetscCall(MatViewFromOptions(A2, NULL, "-fixempty_view"));
712       PetscCall(MatISFixLocalEmpty(A2, PETSC_TRUE));
713       PetscCall(MatAssemblyBegin(A2, MAT_FINAL_ASSEMBLY));
714       PetscCall(MatAssemblyEnd(A2, MAT_FINAL_ASSEMBLY));
715       PetscCall(MatViewFromOptions(A2, NULL, "-fixempty_view"));
716       PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatISFixLocalEmpty (rows+cols)"));
717       PetscCall(MatDestroy(&A2));
718       PetscCall(MatDestroy(&B2));
719     }
720   }
721 
722   /* test MatInvertBlockDiagonal
723        special cases for block-diagonal matrices */
724   if (m == n) {
725     ISLocalToGlobalMapping map;
726     Mat                    Abd, Bbd;
727     IS                     is, bis;
728     const PetscScalar     *isbd, *aijbd;
729     PetscScalar           *vals;
730     const PetscInt        *sts, *idxs;
731     PetscInt              *idxs2, diff, perm, nl, bs, st, en, in;
732     PetscBool              ok;
733 
734     for (diff = 0; diff < 3; diff++) {
735       for (perm = 0; perm < 3; perm++) {
736         for (bs = 1; bs < 4; bs++) {
737           PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatInvertBlockDiagonal blockdiag %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", n, diff, perm, bs));
738           PetscCall(PetscMalloc1(bs * bs, &vals));
739           PetscCall(MatGetOwnershipRanges(A, &sts));
740           switch (diff) {
741           case 1: /* inverted layout by processes */
742             in = 1;
743             st = sts[size - rank - 1];
744             en = sts[size - rank];
745             nl = en - st;
746             break;
747           case 2: /* round-robin layout */
748             in = size;
749             st = rank;
750             nl = n / size;
751             if (rank < n % size) nl++;
752             break;
753           default: /* same layout */
754             in = 1;
755             st = sts[rank];
756             en = sts[rank + 1];
757             nl = en - st;
758             break;
759           }
760           PetscCall(ISCreateStride(PETSC_COMM_WORLD, nl, st, in, &is));
761           PetscCall(ISGetLocalSize(is, &nl));
762           PetscCall(ISGetIndices(is, &idxs));
763           PetscCall(PetscMalloc1(nl, &idxs2));
764           for (i = 0; i < nl; i++) {
765             switch (perm) { /* invert some of the indices */
766             case 2:
767               idxs2[i] = rank % 2 ? idxs[i] : idxs[nl - i - 1];
768               break;
769             case 1:
770               idxs2[i] = rank % 2 ? idxs[nl - i - 1] : idxs[i];
771               break;
772             default:
773               idxs2[i] = idxs[i];
774               break;
775             }
776           }
777           PetscCall(ISRestoreIndices(is, &idxs));
778           PetscCall(ISCreateBlock(PETSC_COMM_WORLD, bs, nl, idxs2, PETSC_OWN_POINTER, &bis));
779           PetscCall(ISLocalToGlobalMappingCreateIS(bis, &map));
780           PetscCall(MatCreateIS(PETSC_COMM_WORLD, bs, PETSC_DECIDE, PETSC_DECIDE, bs * n, bs * n, map, map, &Abd));
781           PetscCall(ISLocalToGlobalMappingDestroy(&map));
782           PetscCall(MatISSetPreallocation(Abd, bs, NULL, 0, NULL));
783           for (i = 0; i < nl; i++) {
784             PetscInt b1, b2;
785 
786             for (b1 = 0; b1 < bs; b1++) {
787               for (b2 = 0; b2 < bs; b2++) vals[b1 * bs + b2] = i * bs * bs + b1 * bs + b2 + 1 + (b1 == b2 ? 1.0 : 0);
788             }
789             PetscCall(MatSetValuesBlockedLocal(Abd, 1, &i, 1, &i, vals, INSERT_VALUES));
790           }
791           PetscCall(MatAssemblyBegin(Abd, MAT_FINAL_ASSEMBLY));
792           PetscCall(MatAssemblyEnd(Abd, MAT_FINAL_ASSEMBLY));
793           PetscCall(MatConvert(Abd, MATAIJ, MAT_INITIAL_MATRIX, &Bbd));
794           PetscCall(MatInvertBlockDiagonal(Abd, &isbd));
795           PetscCall(MatInvertBlockDiagonal(Bbd, &aijbd));
796           PetscCall(MatGetLocalSize(Bbd, &nl, NULL));
797           ok = PETSC_TRUE;
798           for (i = 0; i < nl / bs; i++) {
799             PetscInt b1, b2;
800 
801             for (b1 = 0; b1 < bs; b1++) {
802               for (b2 = 0; b2 < bs; b2++) {
803                 if (PetscAbsScalar(isbd[i * bs * bs + b1 * bs + b2] - aijbd[i * bs * bs + b1 * bs + b2]) > PETSC_SMALL) ok = PETSC_FALSE;
804                 if (!ok) {
805                   PetscCall(PetscPrintf(PETSC_COMM_SELF, "[%d] ERROR block %" PetscInt_FMT ", entry %" PetscInt_FMT " %" PetscInt_FMT ": %g %g\n", rank, i, b1, b2, (double)PetscAbsScalar(isbd[i * bs * bs + b1 * bs + b2]), (double)PetscAbsScalar(aijbd[i * bs * bs + b1 * bs + b2])));
806                   break;
807                 }
808               }
809               if (!ok) break;
810             }
811             if (!ok) break;
812           }
813           PetscCall(MatDestroy(&Abd));
814           PetscCall(MatDestroy(&Bbd));
815           PetscCall(PetscFree(vals));
816           PetscCall(ISDestroy(&is));
817           PetscCall(ISDestroy(&bis));
818         }
819       }
820     }
821   }
822 
823   /* test MatGetDiagonalBlock */
824   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatGetDiagonalBlock\n"));
825   PetscCall(MatGetDiagonalBlock(A, &A2));
826   PetscCall(MatGetDiagonalBlock(B, &B2));
827   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatGetDiagonalBlock"));
828   PetscCall(MatScale(A, 2.0));
829   PetscCall(MatScale(B, 2.0));
830   PetscCall(MatGetDiagonalBlock(A, &A2));
831   PetscCall(MatGetDiagonalBlock(B, &B2));
832   PetscCall(CheckMat(A2, B2, PETSC_FALSE, "MatGetDiagonalBlock"));
833 
834   /* test MatISSetAllowRepeated on a MATIS */
835   PetscCall(MatISSetAllowRepeated(A, allow_repeated));
836   if (allow_repeated) { /* original MATIS maybe with repeated entries, test assembling of local matrices */
837     Mat lA, lA2;
838 
839     for (PetscInt i = 0; i < 1; i++) { /* TODO: make MatScatter inherit from MATSHELL and support MatProducts */
840       PetscBool usemult = PETSC_FALSE;
841 
842       PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
843       if (i) {
844         Mat tA;
845 
846         usemult = PETSC_TRUE;
847         PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatISSetAllowRepeated(false) with possibly repeated entries and local shell matrices\n"));
848         PetscCall(MatISGetLocalMat(A2, &lA2));
849         PetscCall(MatConvert(lA2, MATSHELL, MAT_INITIAL_MATRIX, &tA));
850         PetscCall(MatISRestoreLocalMat(A2, &lA2));
851         PetscCall(MatISSetLocalMat(A2, tA));
852         PetscCall(MatDestroy(&tA));
853       } else {
854         PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatISSetAllowRepeated(false) with possibly repeated entries\n"));
855       }
856       PetscCall(MatISSetAllowRepeated(A2, PETSC_FALSE));
857       PetscCall(MatISGetLocalMat(A, &lA));
858       PetscCall(MatISGetLocalMat(A2, &lA2));
859       if (!repmap) PetscCall(CheckMat(lA, lA2, usemult, "MatISSetAllowRepeated(false) with non-repeated entries"));
860       PetscCall(MatISRestoreLocalMat(A, &lA));
861       PetscCall(MatISRestoreLocalMat(A2, &lA2));
862       if (repmap) PetscCall(CheckMat(A2, B, usemult, "MatISSetAllowRepeated(false) with repeated entries"));
863       else PetscCall(CheckMat(A2, B, usemult, "MatISSetAllowRepeated(false) with non-repeated entries"));
864       PetscCall(MatDestroy(&A2));
865     }
866   } else { /* original matis with non-repeated entries, this should only recreate the local matrices */
867     Mat       lA;
868     PetscBool flg;
869 
870     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatISSetAllowRepeated(true) with non repeated entries\n"));
871     PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &A2));
872     PetscCall(MatISSetAllowRepeated(A2, PETSC_TRUE));
873     PetscCall(MatISGetLocalMat(A2, &lA));
874     PetscCall(MatAssembled(lA, &flg));
875     PetscCheck(!flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Local mat should be unassembled");
876     PetscCall(MatISRestoreLocalMat(A2, &lA));
877     PetscCall(MatDestroy(&A2));
878   }
879 
880   /* free testing matrices */
881   PetscCall(ISLocalToGlobalMappingDestroy(&cmap));
882   PetscCall(ISLocalToGlobalMappingDestroy(&rmap));
883   PetscCall(MatDestroy(&A));
884   PetscCall(MatDestroy(&B));
885   PetscCall(PetscFinalize());
886   return 0;
887 }
888 
889 PetscErrorCode CheckMat(Mat A, Mat B, PetscBool usemult, const char *func)
890 {
891   Mat       Bcheck;
892   PetscReal error;
893 
894   PetscFunctionBeginUser;
895   if (!usemult && B) {
896     PetscBool hasnorm;
897 
898     PetscCall(MatHasOperation(B, MATOP_NORM, &hasnorm));
899     if (!hasnorm) usemult = PETSC_TRUE;
900   }
901   if (!usemult) {
902     if (B) {
903       MatType Btype;
904 
905       PetscCall(MatGetType(B, &Btype));
906       PetscCall(MatConvert(A, Btype, MAT_INITIAL_MATRIX, &Bcheck));
907     } else {
908       PetscCall(MatConvert(A, MATAIJ, MAT_INITIAL_MATRIX, &Bcheck));
909     }
910     if (B) { /* if B is present, subtract it */
911       PetscCall(MatAXPY(Bcheck, -1., B, DIFFERENT_NONZERO_PATTERN));
912     }
913     PetscCall(MatNorm(Bcheck, NORM_INFINITY, &error));
914     if (error > PETSC_SQRT_MACHINE_EPSILON) {
915       ISLocalToGlobalMapping rl2g, cl2g;
916 
917       PetscCall(PetscObjectSetName((PetscObject)Bcheck, "Bcheck"));
918       PetscCall(MatView(Bcheck, NULL));
919       if (B) {
920         PetscCall(PetscObjectSetName((PetscObject)B, "B"));
921         PetscCall(MatView(B, NULL));
922         PetscCall(MatDestroy(&Bcheck));
923         PetscCall(MatConvert(A, MATAIJ, MAT_INITIAL_MATRIX, &Bcheck));
924         PetscCall(PetscObjectSetName((PetscObject)Bcheck, "Assembled A"));
925         PetscCall(MatView(Bcheck, NULL));
926       }
927       PetscCall(MatDestroy(&Bcheck));
928       PetscCall(PetscObjectSetName((PetscObject)A, "A"));
929       PetscCall(MatView(A, NULL));
930       PetscCall(MatGetLocalToGlobalMapping(A, &rl2g, &cl2g));
931       if (rl2g) PetscCall(ISLocalToGlobalMappingView(rl2g, NULL));
932       if (cl2g) PetscCall(ISLocalToGlobalMappingView(cl2g, NULL));
933       SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "ERROR ON %s: %g", func, (double)error);
934     }
935     PetscCall(MatDestroy(&Bcheck));
936   } else {
937     PetscBool ok, okt;
938 
939     PetscCall(MatMultEqual(A, B, 3, &ok));
940     PetscCall(MatMultTransposeEqual(A, B, 3, &okt));
941     PetscCheck(ok && okt, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "ERROR ON %s: mult ok ?  %d, multtranspose ok ? %d", func, ok, okt);
942   }
943   PetscFunctionReturn(PETSC_SUCCESS);
944 }
945 
946 PetscErrorCode TestMatZeroRows(Mat A, Mat Afull, PetscBool squaretest, IS is, PetscScalar diag)
947 {
948   Mat                    B, Bcheck, B2 = NULL, lB;
949   Vec                    x = NULL, b = NULL, b2 = NULL;
950   ISLocalToGlobalMapping l2gr, l2gc;
951   PetscReal              error;
952   char                   diagstr[16];
953   const PetscInt        *idxs;
954   PetscInt               rst, ren, i, n, N, d;
955   PetscMPIInt            rank;
956   PetscBool              miss, haszerorows;
957 
958   PetscFunctionBeginUser;
959   if (diag == 0.) {
960     PetscCall(PetscStrncpy(diagstr, "zero", sizeof(diagstr)));
961   } else {
962     PetscCall(PetscStrncpy(diagstr, "nonzero", sizeof(diagstr)));
963   }
964   PetscCall(ISView(is, NULL));
965   PetscCall(MatGetLocalToGlobalMapping(A, &l2gr, &l2gc));
966   /* tests MatDuplicate and MatCopy */
967   if (diag == 0.) {
968     PetscCall(MatDuplicate(A, MAT_COPY_VALUES, &B));
969   } else {
970     PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &B));
971     PetscCall(MatCopy(A, B, SAME_NONZERO_PATTERN));
972   }
973   PetscCall(MatISGetLocalMat(B, &lB));
974   PetscCall(MatHasOperation(lB, MATOP_ZERO_ROWS, &haszerorows));
975   if (squaretest && haszerorows) {
976     PetscCall(MatCreateVecs(B, &x, &b));
977     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &B2));
978     PetscCall(VecSetLocalToGlobalMapping(b, l2gr));
979     PetscCall(VecSetLocalToGlobalMapping(x, l2gc));
980     PetscCall(VecSetRandom(x, NULL));
981     PetscCall(VecSetRandom(b, NULL));
982     /* mimic b[is] = x[is] */
983     PetscCall(VecDuplicate(b, &b2));
984     PetscCall(VecSetLocalToGlobalMapping(b2, l2gr));
985     PetscCall(VecCopy(b, b2));
986     PetscCall(ISGetLocalSize(is, &n));
987     PetscCall(ISGetIndices(is, &idxs));
988     PetscCall(VecGetSize(x, &N));
989     for (i = 0; i < n; i++) {
990       if (0 <= idxs[i] && idxs[i] < N) {
991         PetscCall(VecSetValue(b2, idxs[i], diag, INSERT_VALUES));
992         PetscCall(VecSetValue(x, idxs[i], 1., INSERT_VALUES));
993       }
994     }
995     PetscCall(VecAssemblyBegin(b2));
996     PetscCall(VecAssemblyEnd(b2));
997     PetscCall(VecAssemblyBegin(x));
998     PetscCall(VecAssemblyEnd(x));
999     PetscCall(ISRestoreIndices(is, &idxs));
1000     /*  test ZeroRows on MATIS */
1001     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatZeroRows (diag %s)\n", diagstr));
1002     PetscCall(MatZeroRowsIS(B, is, diag, x, b));
1003     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatZeroRowsColumns (diag %s)\n", diagstr));
1004     PetscCall(MatZeroRowsColumnsIS(B2, is, diag, NULL, NULL));
1005   } else if (haszerorows) {
1006     /*  test ZeroRows on MATIS */
1007     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatZeroRows (diag %s)\n", diagstr));
1008     PetscCall(MatZeroRowsIS(B, is, diag, NULL, NULL));
1009     b = b2 = x = NULL;
1010   } else {
1011     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Skipping MatZeroRows (diag %s)\n", diagstr));
1012     b = b2 = x = NULL;
1013   }
1014 
1015   if (squaretest && haszerorows) {
1016     PetscCall(VecAXPY(b2, -1., b));
1017     PetscCall(VecNorm(b2, NORM_INFINITY, &error));
1018     PetscCheck(error <= PETSC_SQRT_MACHINE_EPSILON, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "ERROR IN ZEROROWS ON B %g (diag %s)", (double)error, diagstr);
1019   }
1020 
1021   /* test MatMissingDiagonal */
1022   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test MatMissingDiagonal\n"));
1023   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
1024   PetscCall(MatMissingDiagonal(B, &miss, &d));
1025   PetscCall(MatGetOwnershipRange(B, &rst, &ren));
1026   PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD));
1027   PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] [%" PetscInt_FMT ",%" PetscInt_FMT ") Missing %d, row %" PetscInt_FMT " (diag %s)\n", rank, rst, ren, (int)miss, d, diagstr));
1028   PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
1029   PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD));
1030 
1031   PetscCall(VecDestroy(&x));
1032   PetscCall(VecDestroy(&b));
1033   PetscCall(VecDestroy(&b2));
1034 
1035   /* check the result of ZeroRows with that from MPIAIJ routines
1036      assuming that MatConvert_IS_XAIJ and MatZeroRows_MPIAIJ work fine */
1037   if (haszerorows) {
1038     PetscCall(MatDuplicate(Afull, MAT_COPY_VALUES, &Bcheck));
1039     PetscCall(MatSetOption(Bcheck, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
1040     PetscCall(MatZeroRowsIS(Bcheck, is, diag, NULL, NULL));
1041     PetscCall(CheckMat(B, Bcheck, PETSC_FALSE, "Zerorows"));
1042     PetscCall(MatDestroy(&Bcheck));
1043   }
1044   PetscCall(MatDestroy(&B));
1045 
1046   if (B2) { /* test MatZeroRowsColumns */
1047     PetscCall(MatDuplicate(Afull, MAT_COPY_VALUES, &B));
1048     PetscCall(MatSetOption(B, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
1049     PetscCall(MatZeroRowsColumnsIS(B, is, diag, NULL, NULL));
1050     PetscCall(CheckMat(B2, B, PETSC_FALSE, "MatZeroRowsColumns"));
1051     PetscCall(MatDestroy(&B));
1052     PetscCall(MatDestroy(&B2));
1053   }
1054   PetscFunctionReturn(PETSC_SUCCESS);
1055 }
1056 
1057 /*TEST
1058 
1059    test:
1060       requires: !complex
1061       args: -test_matlab -test_trans
1062 
1063    test:
1064       suffix: 2
1065       nsize: 4
1066       args: -mat_is_convert_local_nest -nr 3 -nc 4
1067 
1068    test:
1069       suffix: 3
1070       nsize: 5
1071       args: -m 11 -n 10 -mat_is_convert_local_nest -nr 2 -nc 1
1072 
1073    test:
1074       suffix: 4
1075       nsize: 6
1076       args: -m 9 -n 12 -test_trans -nr 2 -nc 7
1077 
1078    test:
1079       suffix: 5
1080       nsize: 6
1081       args: -m 12 -n 12 -test_trans -nr 3 -nc 1
1082 
1083    test:
1084       suffix: 6
1085       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap
1086 
1087    test:
1088       suffix: 7
1089       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap -permmap
1090 
1091    test:
1092       suffix: 8
1093       args: -m 12 -n 17 -test_trans -nr 2 -nc 3 -permmap
1094 
1095    test:
1096       suffix: 9
1097       nsize: 5
1098       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap
1099 
1100    test:
1101       suffix: 10
1102       nsize: 5
1103       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap -permmap
1104 
1105    test:
1106       suffix: vscat_default
1107       nsize: 5
1108       args: -m 12 -n 17 -test_trans -nr 2 -nc 3 -permmap
1109       output_file: output/ex23_11.out
1110 
1111    test:
1112       suffix: 12
1113       nsize: 3
1114       args: -m 12 -n 12 -symmetric -mat_is_localmat_type sbaij -test_trans -nr 2 -nc 3
1115 
1116    testset:
1117       output_file: output/ex23_13.out
1118       nsize: 3
1119       args: -m 12 -n 17 -test_trans -nr 2 -nc 3 -diffmap -permmap
1120       filter: grep -v "type:" | grep -v "block size is 1" | grep -v "not using I-node routines"
1121       test:
1122         suffix: baij
1123         args: -mat_is_localmat_type baij
1124       test:
1125         requires: viennacl
1126         suffix: viennacl
1127         args: -mat_is_localmat_type aijviennacl
1128       test:
1129         requires: cuda
1130         suffix: cusparse
1131         args: -mat_is_localmat_type aijcusparse
1132       test:
1133         requires: kokkos_kernels
1134         suffix: kokkos
1135         args: -mat_is_localmat_type aijkokkos
1136 
1137    test:
1138       suffix: negrep
1139       nsize: {{1 3}separate output}
1140       args: -m {{5 7}separate output} -n {{5 7}separate output} -test_trans -nr 2 -nc 3 -negmap {{0 1}separate output} -repmap {{0 1}separate output} -permmap -diffmap {{0 1}separate output} -allow_repeated 0
1141 
1142    test:
1143       suffix: negrep_allowrep
1144       nsize: {{1 3}separate output}
1145       args: -m {{5 7}separate output} -n {{5 7}separate output} -test_trans -nr 2 -nc 3 -negmap {{0 1}separate output} -repmap {{0 1}separate output} -permmap -diffmap {{0 1}separate output} -allow_repeated
1146 
1147 TEST*/
1148