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