Actual source code: ex23.c


  2: static char help[] = "Tests the use of interface functions for MATIS matrices.\n\
  3: This example tests: MatZeroRows(), MatZeroRowsLocal(), MatView(), MatDuplicate(),\n\
  4: MatCopy(), MatCreateSubMatrix(), MatGetLocalSubMatrix(), MatAXPY(), MatShift()\n\
  5: MatDiagonalSet(), MatTranspose() and MatPtAP(). It also tests some\n\
  6: conversion routines.\n";

  8: #include <petscmat.h>

 10: PetscErrorCode TestMatZeroRows(Mat,Mat,PetscBool,IS,PetscScalar);
 11: PetscErrorCode CheckMat(Mat,Mat,PetscBool,const char*);

 13: int main(int argc,char **args)
 14: {
 15:   Mat                    A,B,A2,B2,T;
 16:   Mat                    Aee,Aeo,Aoe,Aoo;
 17:   Mat                    *mats;
 18:   Vec                    x,y;
 19:   MatInfo                info;
 20:   ISLocalToGlobalMapping cmap,rmap;
 21:   IS                     is,is2,reven,rodd,ceven,codd;
 22:   IS                     *rows,*cols;
 23:   MatType                lmtype;
 24:   PetscScalar            diag = 2.;
 25:   PetscInt               n,m,i,lm,ln;
 26:   PetscInt               rst,ren,cst,cen,nr,nc;
 27:   PetscMPIInt            rank,size;
 28:   PetscBool              testT,squaretest,isaij;
 29:   PetscBool              permute = PETSC_FALSE, negmap = PETSC_FALSE, repmap = PETSC_FALSE;
 30:   PetscBool              diffmap = PETSC_TRUE, symmetric = PETSC_FALSE, issymmetric;
 31:   PetscErrorCode         ierr;

 33:   PetscInitialize(&argc,&args,(char*)0,help);
 34:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 35:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 36:   m = n = 2*size;
 37:   PetscOptionsGetBool(NULL,NULL,"-symmetric",&symmetric,NULL);
 38:   PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
 39:   PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
 40:   PetscOptionsGetBool(NULL,NULL,"-negmap",&negmap,NULL);
 41:   PetscOptionsGetBool(NULL,NULL,"-repmap",&repmap,NULL);
 42:   PetscOptionsGetBool(NULL,NULL,"-permmap",&permute,NULL);
 43:   PetscOptionsGetBool(NULL,NULL,"-diffmap",&diffmap,NULL);

 48:   /* create a MATIS matrix */
 49:   MatCreate(PETSC_COMM_WORLD,&A);
 50:   MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);
 51:   MatSetType(A,MATIS);
 52:   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:     ISCreateStride(PETSC_COMM_WORLD,n,0,1,&is);
 58:   } else if (negmap && !repmap) { /* non repeated but with negative indices */
 59:     ISCreateStride(PETSC_COMM_WORLD,n+2,-2,1,&is);
 60:   } else if (!negmap && repmap) { /* non negative but repeated indices */
 61:     IS isl[2];

 63:     ISCreateStride(PETSC_COMM_WORLD,n,0,1,&isl[0]);
 64:     ISCreateStride(PETSC_COMM_WORLD,n,n-1,-1,&isl[1]);
 65:     ISConcatenate(PETSC_COMM_WORLD,2,isl,&is);
 66:     ISDestroy(&isl[0]);
 67:     ISDestroy(&isl[1]);
 68:   } else { /* negative and repeated indices */
 69:     IS isl[2];

 71:     ISCreateStride(PETSC_COMM_WORLD,n+1,-1,1,&isl[0]);
 72:     ISCreateStride(PETSC_COMM_WORLD,n+1,n-1,-1,&isl[1]);
 73:     ISConcatenate(PETSC_COMM_WORLD,2,isl,&is);
 74:     ISDestroy(&isl[0]);
 75:     ISDestroy(&isl[1]);
 76:   }
 77:   ISLocalToGlobalMappingCreateIS(is,&cmap);
 78:   ISDestroy(&is);

 80:   if (m != n || diffmap) {
 81:     ISCreateStride(PETSC_COMM_WORLD,m,permute ? m-1 : 0,permute ? -1 : 1,&is);
 82:     ISLocalToGlobalMappingCreateIS(is,&rmap);
 83:     ISDestroy(&is);
 84:   } else {
 85:     PetscObjectReference((PetscObject)cmap);
 86:     rmap = cmap;
 87:   }

 89:   MatSetLocalToGlobalMapping(A,rmap,cmap);
 90:   MatISStoreL2L(A,PETSC_FALSE);
 91:   MatISSetPreallocation(A,3,NULL,3,NULL);
 92:   MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,(PetscBool)!(repmap || negmap)); /* I do not want to precompute the pattern */
 93:   ISLocalToGlobalMappingGetSize(rmap,&lm);
 94:   ISLocalToGlobalMappingGetSize(cmap,&ln);
 95:   for (i=0; i<lm; i++) {
 96:     PetscScalar v[3];
 97:     PetscInt    cols[3];

 99:     cols[0] = (i-1+n)%n;
100:     cols[1] = i%n;
101:     cols[2] = (i+1)%n;
102:     v[0] = -1.*(symmetric ? PetscMin(i+1,cols[0]+1) : i+1);
103:     v[1] =  2.*(symmetric ? PetscMin(i+1,cols[1]+1) : i+1);
104:     v[2] = -1.*(symmetric ? PetscMin(i+1,cols[2]+1) : i+1);
105:     ISGlobalToLocalMappingApply(cmap,IS_GTOLM_MASK,3,cols,NULL,cols);
106:     MatSetValuesLocal(A,1,&i,3,cols,v,ADD_VALUES);
107:   }
108:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
109:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

111:   /* activate tests for square matrices with same maps only */
112:   MatHasCongruentLayouts(A,&squaretest);
113:   if (squaretest && rmap != cmap) {
114:     PetscInt nr, nc;

116:     ISLocalToGlobalMappingGetSize(rmap,&nr);
117:     ISLocalToGlobalMappingGetSize(cmap,&nc);
118:     if (nr != nc) squaretest = PETSC_FALSE;
119:     else {
120:       const PetscInt *idxs1,*idxs2;

122:       ISLocalToGlobalMappingGetIndices(rmap,&idxs1);
123:       ISLocalToGlobalMappingGetIndices(cmap,&idxs2);
124:       PetscArraycmp(idxs1,idxs2,nr,&squaretest);
125:       ISLocalToGlobalMappingRestoreIndices(rmap,&idxs1);
126:       ISLocalToGlobalMappingRestoreIndices(cmap,&idxs2);
127:     }
128:     MPIU_Allreduce(MPI_IN_PLACE,&squaretest,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)A));
129:   }

131:   /* test MatISGetLocalMat */
132:   MatISGetLocalMat(A,&B);
133:   MatGetType(B,&lmtype);

135:   /* test MatGetInfo */
136:   PetscPrintf(PETSC_COMM_WORLD,"Test MatGetInfo\n");
137:   MatGetInfo(A,MAT_LOCAL,&info);
138:   PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD);
139:   PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"Process  %2d: %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",PetscGlobalRank,(PetscInt)info.nz_used,
140:                                             (PetscInt)info.nz_allocated,(PetscInt)info.nz_unneeded,(PetscInt)info.assemblies,(PetscInt)info.mallocs);
141:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
142:   MatGetInfo(A,MAT_GLOBAL_MAX,&info);
143:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"GlobalMax  : %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",(PetscInt)info.nz_used,
144:                                 (PetscInt)info.nz_allocated,(PetscInt)info.nz_unneeded,(PetscInt)info.assemblies,(PetscInt)info.mallocs);
145:   MatGetInfo(A,MAT_GLOBAL_SUM,&info);
146:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"GlobalSum  : %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",(PetscInt)info.nz_used,
147:                                 (PetscInt)info.nz_allocated,(PetscInt)info.nz_unneeded,(PetscInt)info.assemblies,(PetscInt)info.mallocs);

149:   /* test MatIsSymmetric */
150:   MatIsSymmetric(A,0.0,&issymmetric);
151:   PetscPrintf(PETSC_COMM_WORLD,"Test MatIsSymmetric: %d\n",issymmetric);

153:   /* Create a MPIAIJ matrix, same as A */
154:   MatCreate(PETSC_COMM_WORLD,&B);
155:   MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,m,n);
156:   MatSetType(B,MATAIJ);
157:   MatSetFromOptions(B);
158:   MatSetUp(B);
159:   MatSetLocalToGlobalMapping(B,rmap,cmap);
160:   MatMPIAIJSetPreallocation(B,3,NULL,3,NULL);
161:   MatMPIBAIJSetPreallocation(B,1,3,NULL,3,NULL);
162: #if defined(PETSC_HAVE_HYPRE)
163:   MatHYPRESetPreallocation(B,3,NULL,3,NULL);
164: #endif
165:   MatISSetPreallocation(B,3,NULL,3,NULL);
166:   MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,(PetscBool)!(repmap || negmap)); /* I do not want to precompute the pattern */
167:   for (i=0; i<lm; i++) {
168:     PetscScalar v[3];
169:     PetscInt    cols[3];

171:     cols[0] = (i-1+n)%n;
172:     cols[1] = i%n;
173:     cols[2] = (i+1)%n;
174:     v[0] = -1.*(symmetric ? PetscMin(i+1,cols[0]+1) : i+1);
175:     v[1] =  2.*(symmetric ? PetscMin(i+1,cols[1]+1) : i+1);
176:     v[2] = -1.*(symmetric ? PetscMin(i+1,cols[2]+1) : i+1);
177:     ISGlobalToLocalMappingApply(cmap,IS_GTOLM_MASK,3,cols,NULL,cols);
178:     MatSetValuesLocal(B,1,&i,3,cols,v,ADD_VALUES);
179:   }
180:   MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
181:   MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);

183:   /* test MatView */
184:   PetscPrintf(PETSC_COMM_WORLD,"Test MatView\n");
185:   MatView(A,NULL);
186:   MatView(B,NULL);

188:   /* test CheckMat */
189:   PetscPrintf(PETSC_COMM_WORLD,"Test CheckMat\n");
190:   CheckMat(A,B,PETSC_FALSE,"CheckMat");

192:   /* test MatDuplicate and MatAXPY */
193:   PetscPrintf(PETSC_COMM_WORLD,"Test MatDuplicate and MatAXPY\n");
194:   MatDuplicate(A,MAT_COPY_VALUES,&A2);
195:   CheckMat(A,A2,PETSC_FALSE,"MatDuplicate and MatAXPY");

197:   /* test MatConvert */
198:   PetscPrintf(PETSC_COMM_WORLD,"Test MatConvert_IS_XAIJ\n");
199:   MatConvert(A2,MATAIJ,MAT_INITIAL_MATRIX,&B2);
200:   CheckMat(B,B2,PETSC_TRUE,"MatConvert_IS_XAIJ MAT_INITIAL_MATRIX");
201:   MatConvert(A2,MATAIJ,MAT_REUSE_MATRIX,&B2);
202:   CheckMat(B,B2,PETSC_TRUE,"MatConvert_IS_XAIJ MAT_REUSE_MATRIX");
203:   MatConvert(A2,MATAIJ,MAT_INPLACE_MATRIX,&A2);
204:   CheckMat(B,A2,PETSC_TRUE,"MatConvert_IS_XAIJ MAT_INPLACE_MATRIX");
205:   MatDestroy(&A2);
206:   MatDestroy(&B2);
207:   PetscPrintf(PETSC_COMM_WORLD,"Test MatConvert_XAIJ_IS\n");
208:   MatDuplicate(B,MAT_COPY_VALUES,&B2);
209:   MatConvert(B2,MATIS,MAT_INITIAL_MATRIX,&A2);
210:   CheckMat(A,A2,PETSC_TRUE,"MatConvert_XAIJ_IS MAT_INITIAL_MATRIX");
211:   MatConvert(B2,MATIS,MAT_REUSE_MATRIX,&A2);
212:   CheckMat(A,A2,PETSC_TRUE,"MatConvert_XAIJ_IS MAT_REUSE_MATRIX");
213:   MatConvert(B2,MATIS,MAT_INPLACE_MATRIX,&B2);
214:   CheckMat(A,B2,PETSC_TRUE,"MatConvert_XAIJ_IS MAT_INPLACE_MATRIX");
215:   MatDestroy(&A2);
216:   MatDestroy(&B2);
217:   PetscStrcmp(lmtype,MATSEQAIJ,&isaij);
218:   if (size == 1 && isaij) { /* tests special code paths in MatConvert_IS_XAIJ */
219:     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};
220:     ISLocalToGlobalMapping tcmap,trmap;

222:     for (ri = 0; ri < 2; ri++) {
223:       PetscInt *r;

225:       r = (PetscInt*)(ri == 0 ? rr : rk);
226:       for (ci = 0; ci < 2; ci++) {
227:         PetscInt *c,rb,cb;

229:         c = (PetscInt*)(ci == 0 ? cr : ck);
230:         for (rb = 1; rb < 4; rb++) {
231:           ISCreateBlock(PETSC_COMM_SELF,rb,3,r,PETSC_COPY_VALUES,&is);
232:           ISLocalToGlobalMappingCreateIS(is,&trmap);
233:           ISDestroy(&is);
234:           for (cb = 1; cb < 4; cb++) {
235:             Mat  T,lT,T2;
236:             char testname[256];

238:             PetscSNPrintf(testname,sizeof(testname),"MatConvert_IS_XAIJ special case (%" PetscInt_FMT " %" PetscInt_FMT ", bs %" PetscInt_FMT " %" PetscInt_FMT ")",ri,ci,rb,cb);
239:             PetscPrintf(PETSC_COMM_WORLD,"Test %s\n",testname);

241:             ISCreateBlock(PETSC_COMM_SELF,cb,4,c,PETSC_COPY_VALUES,&is);
242:             ISLocalToGlobalMappingCreateIS(is,&tcmap);
243:             ISDestroy(&is);

245:             MatCreate(PETSC_COMM_SELF,&T);
246:             MatSetSizes(T,PETSC_DECIDE,PETSC_DECIDE,rb*3,cb*4);
247:             MatSetType(T,MATIS);
248:             MatSetLocalToGlobalMapping(T,trmap,tcmap);
249:             ISLocalToGlobalMappingDestroy(&tcmap);
250:             MatISGetLocalMat(T,&lT);
251:             MatSetType(lT,MATSEQAIJ);
252:             MatSeqAIJSetPreallocation(lT,cb*4,NULL);
253:             MatSetRandom(lT,NULL);
254:             MatConvert(lT,lmtype,MAT_INPLACE_MATRIX,&lT);
255:             MatISRestoreLocalMat(T,&lT);
256:             MatAssemblyBegin(T,MAT_FINAL_ASSEMBLY);
257:             MatAssemblyEnd(T,MAT_FINAL_ASSEMBLY);

259:             MatConvert(T,MATAIJ,MAT_INITIAL_MATRIX,&T2);
260:             CheckMat(T,T2,PETSC_TRUE,"MAT_INITIAL_MATRIX");
261:             MatConvert(T,MATAIJ,MAT_REUSE_MATRIX,&T2);
262:             CheckMat(T,T2,PETSC_TRUE,"MAT_REUSE_MATRIX");
263:             MatDestroy(&T2);
264:             MatDuplicate(T,MAT_COPY_VALUES,&T2);
265:             MatConvert(T2,MATAIJ,MAT_INPLACE_MATRIX,&T2);
266:             CheckMat(T,T2,PETSC_TRUE,"MAT_INPLACE_MATRIX");
267:             MatDestroy(&T);
268:             MatDestroy(&T2);
269:           }
270:           ISLocalToGlobalMappingDestroy(&trmap);
271:         }
272:       }
273:     }
274:   }

276:   /* test MatDiagonalScale */
277:   PetscPrintf(PETSC_COMM_WORLD,"Test MatDiagonalScale\n");
278:   MatDuplicate(A,MAT_COPY_VALUES,&A2);
279:   MatDuplicate(B,MAT_COPY_VALUES,&B2);
280:   MatCreateVecs(A,&x,&y);
281:   VecSetRandom(x,NULL);
282:   if (issymmetric) {
283:     VecCopy(x,y);
284:   } else {
285:     VecSetRandom(y,NULL);
286:     VecScale(y,8.);
287:   }
288:   MatDiagonalScale(A2,y,x);
289:   MatDiagonalScale(B2,y,x);
290:   CheckMat(A2,B2,PETSC_FALSE,"MatDiagonalScale");
291:   MatDestroy(&A2);
292:   MatDestroy(&B2);
293:   VecDestroy(&x);
294:   VecDestroy(&y);

296:   /* test MatPtAP (A IS and B AIJ) */
297:   if (isaij && m == n) {
298:     PetscPrintf(PETSC_COMM_WORLD,"Test MatPtAP\n");
299:     MatISStoreL2L(A,PETSC_TRUE);
300:     MatPtAP(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&A2);
301:     MatPtAP(B,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&B2);
302:     CheckMat(A2,B2,PETSC_FALSE,"MatPtAP MAT_INITIAL_MATRIX");
303:     MatPtAP(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&A2);
304:     CheckMat(A2,B2,PETSC_FALSE,"MatPtAP MAT_REUSE_MATRIX");
305:     MatDestroy(&A2);
306:     MatDestroy(&B2);
307:   }

309:   /* test MatGetLocalSubMatrix */
310:   PetscPrintf(PETSC_COMM_WORLD,"Test MatGetLocalSubMatrix\n");
311:   MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,&A2);
312:   ISCreateStride(PETSC_COMM_SELF,lm/2+lm%2,0,2,&reven);
313:   ISComplement(reven,0,lm,&rodd);
314:   ISCreateStride(PETSC_COMM_SELF,ln/2+ln%2,0,2,&ceven);
315:   ISComplement(ceven,0,ln,&codd);
316:   MatGetLocalSubMatrix(A2,reven,ceven,&Aee);
317:   MatGetLocalSubMatrix(A2,reven,codd,&Aeo);
318:   MatGetLocalSubMatrix(A2,rodd,ceven,&Aoe);
319:   MatGetLocalSubMatrix(A2,rodd,codd,&Aoo);
320:   for (i=0; i<lm; i++) {
321:     PetscInt    j,je,jo,colse[3], colso[3];
322:     PetscScalar ve[3], vo[3];
323:     PetscScalar v[3];
324:     PetscInt    cols[3];
325:     PetscInt    row = i/2;

327:     cols[0] = (i-1+n)%n;
328:     cols[1] = i%n;
329:     cols[2] = (i+1)%n;
330:     v[0] = -1.*(symmetric ? PetscMin(i+1,cols[0]+1) : i+1);
331:     v[1] =  2.*(symmetric ? PetscMin(i+1,cols[1]+1) : i+1);
332:     v[2] = -1.*(symmetric ? PetscMin(i+1,cols[2]+1) : i+1);
333:     ISGlobalToLocalMappingApply(cmap,IS_GTOLM_MASK,3,cols,NULL,cols);
334:     for (j=0,je=0,jo=0;j<3;j++) {
335:       if (cols[j]%2) {
336:         vo[jo] = v[j];
337:         colso[jo++] = cols[j]/2;
338:       } else {
339:         ve[je] = v[j];
340:         colse[je++] = cols[j]/2;
341:       }
342:     }
343:     if (i%2) {
344:       MatSetValuesLocal(Aoe,1,&row,je,colse,ve,ADD_VALUES);
345:       MatSetValuesLocal(Aoo,1,&row,jo,colso,vo,ADD_VALUES);
346:     } else {
347:       MatSetValuesLocal(Aee,1,&row,je,colse,ve,ADD_VALUES);
348:       MatSetValuesLocal(Aeo,1,&row,jo,colso,vo,ADD_VALUES);
349:     }
350:   }
351:   MatRestoreLocalSubMatrix(A2,reven,ceven,&Aee);
352:   MatRestoreLocalSubMatrix(A2,reven,codd,&Aeo);
353:   MatRestoreLocalSubMatrix(A2,rodd,ceven,&Aoe);
354:   MatRestoreLocalSubMatrix(A2,rodd,codd,&Aoo);
355:   ISDestroy(&reven);
356:   ISDestroy(&ceven);
357:   ISDestroy(&rodd);
358:   ISDestroy(&codd);
359:   MatAssemblyBegin(A2,MAT_FINAL_ASSEMBLY);
360:   MatAssemblyEnd(A2,MAT_FINAL_ASSEMBLY);
361:   MatAXPY(A2,-1.,A,SAME_NONZERO_PATTERN);
362:   CheckMat(A2,NULL,PETSC_FALSE,"MatGetLocalSubMatrix");
363:   MatDestroy(&A2);

365:   /* test MatConvert_Nest_IS */
366:   testT = PETSC_FALSE;
367:   PetscOptionsGetBool(NULL,NULL,"-test_trans",&testT,NULL);

369:   PetscPrintf(PETSC_COMM_WORLD,"Test MatConvert_Nest_IS\n");
370:   nr   = 2;
371:   nc   = 2;
372:   PetscOptionsGetInt(NULL,NULL,"-nr",&nr,NULL);
373:   PetscOptionsGetInt(NULL,NULL,"-nc",&nc,NULL);
374:   if (testT) {
375:     MatGetOwnershipRange(A,&cst,&cen);
376:     MatGetOwnershipRangeColumn(A,&rst,&ren);
377:   } else {
378:     MatGetOwnershipRange(A,&rst,&ren);
379:     MatGetOwnershipRangeColumn(A,&cst,&cen);
380:   }
381:   PetscMalloc3(nr,&rows,nc,&cols,2*nr*nc,&mats);
382:   for (i=0;i<nr*nc;i++) {
383:     if (testT) {
384:       MatCreateTranspose(A,&mats[i]);
385:       MatTranspose(B,MAT_INITIAL_MATRIX,&mats[i+nr*nc]);
386:     } else {
387:       MatDuplicate(A,MAT_COPY_VALUES,&mats[i]);
388:       MatDuplicate(B,MAT_COPY_VALUES,&mats[i+nr*nc]);
389:     }
390:   }
391:   for (i=0;i<nr;i++) {
392:     ISCreateStride(PETSC_COMM_WORLD,ren-rst,i+rst,nr,&rows[i]);
393:   }
394:   for (i=0;i<nc;i++) {
395:     ISCreateStride(PETSC_COMM_WORLD,cen-cst,i+cst,nc,&cols[i]);
396:   }
397:   MatCreateNest(PETSC_COMM_WORLD,nr,rows,nc,cols,mats,&A2);
398:   MatCreateNest(PETSC_COMM_WORLD,nr,rows,nc,cols,mats+nr*nc,&B2);
399:   for (i=0;i<nr;i++) {
400:     ISDestroy(&rows[i]);
401:   }
402:   for (i=0;i<nc;i++) {
403:     ISDestroy(&cols[i]);
404:   }
405:   for (i=0;i<2*nr*nc;i++) {
406:     MatDestroy(&mats[i]);
407:   }
408:   PetscFree3(rows,cols,mats);
409:   MatConvert(B2,MATAIJ,MAT_INITIAL_MATRIX,&T);
410:   MatDestroy(&B2);
411:   MatConvert(A2,MATIS,MAT_INITIAL_MATRIX,&B2);
412:   CheckMat(B2,T,PETSC_TRUE,"MatConvert_Nest_IS MAT_INITIAL_MATRIX");
413:   MatConvert(A2,MATIS,MAT_REUSE_MATRIX,&B2);
414:   CheckMat(B2,T,PETSC_TRUE,"MatConvert_Nest_IS MAT_REUSE_MATRIX");
415:   MatDestroy(&B2);
416:   MatConvert(A2,MATIS,MAT_INPLACE_MATRIX,&A2);
417:   CheckMat(A2,T,PETSC_TRUE,"MatConvert_Nest_IS MAT_INPLACE_MATRIX");
418:   MatDestroy(&T);
419:   MatDestroy(&A2);

421:   /* test MatCreateSubMatrix */
422:   PetscPrintf(PETSC_COMM_WORLD,"Test MatCreateSubMatrix\n");
423:   if (rank == 0) {
424:     ISCreateStride(PETSC_COMM_WORLD,1,1,1,&is);
425:     ISCreateStride(PETSC_COMM_WORLD,2,0,1,&is2);
426:   } else if (rank == 1) {
427:     ISCreateStride(PETSC_COMM_WORLD,1,0,1,&is);
428:     if (n > 3) {
429:       ISCreateStride(PETSC_COMM_WORLD,1,3,1,&is2);
430:     } else {
431:       ISCreateStride(PETSC_COMM_WORLD,0,0,1,&is2);
432:     }
433:   } else if (rank == 2 && n > 4) {
434:     ISCreateStride(PETSC_COMM_WORLD,0,0,1,&is);
435:     ISCreateStride(PETSC_COMM_WORLD,n-4,4,1,&is2);
436:   } else {
437:     ISCreateStride(PETSC_COMM_WORLD,0,0,1,&is);
438:     ISCreateStride(PETSC_COMM_WORLD,0,0,1,&is2);
439:   }
440:   MatCreateSubMatrix(A,is,is,MAT_INITIAL_MATRIX,&A2);
441:   MatCreateSubMatrix(B,is,is,MAT_INITIAL_MATRIX,&B2);
442:   CheckMat(A2,B2,PETSC_TRUE,"first MatCreateSubMatrix");

444:   MatCreateSubMatrix(A,is,is,MAT_REUSE_MATRIX,&A2);
445:   MatCreateSubMatrix(B,is,is,MAT_REUSE_MATRIX,&B2);
446:   CheckMat(A2,B2,PETSC_FALSE,"reuse MatCreateSubMatrix");
447:   MatDestroy(&A2);
448:   MatDestroy(&B2);

450:   if (!issymmetric) {
451:     MatCreateSubMatrix(A,is,is2,MAT_INITIAL_MATRIX,&A2);
452:     MatCreateSubMatrix(B,is,is2,MAT_INITIAL_MATRIX,&B2);
453:     MatCreateSubMatrix(A,is,is2,MAT_REUSE_MATRIX,&A2);
454:     MatCreateSubMatrix(B,is,is2,MAT_REUSE_MATRIX,&B2);
455:     CheckMat(A2,B2,PETSC_FALSE,"second MatCreateSubMatrix");
456:   }

458:   MatDestroy(&A2);
459:   MatDestroy(&B2);
460:   ISDestroy(&is);
461:   ISDestroy(&is2);

463:   /* Create an IS required by MatZeroRows(): just rank zero provides the rows to be eliminated */
464:   if (size > 1) {
465:     if (rank == 0) {
466:       PetscInt st,len;

468:       st   = (m+1)/2;
469:       len  = PetscMin(m/2,PetscMax(m-(m+1)/2-1,0));
470:       ISCreateStride(PETSC_COMM_WORLD,len,st,1,&is);
471:     } else {
472:       ISCreateStride(PETSC_COMM_WORLD,0,0,1,&is);
473:     }
474:   } else {
475:     ISCreateStride(PETSC_COMM_WORLD,1,0,1,&is);
476:   }

478:   if (squaretest) { /* tests for square matrices only, with same maps for rows and columns */
479:     /* test MatDiagonalSet */
480:     PetscPrintf(PETSC_COMM_WORLD,"Test MatDiagonalSet\n");
481:     MatDuplicate(A,MAT_COPY_VALUES,&A2);
482:     MatDuplicate(B,MAT_COPY_VALUES,&B2);
483:     MatCreateVecs(A,NULL,&x);
484:     VecSetRandom(x,NULL);
485:     MatDiagonalSet(A2,x,INSERT_VALUES);
486:     MatDiagonalSet(B2,x,INSERT_VALUES);
487:     CheckMat(A2,B2,PETSC_FALSE,"MatDiagonalSet");
488:     VecDestroy(&x);
489:     MatDestroy(&A2);
490:     MatDestroy(&B2);

492:     /* test MatShift (MatShift_IS internally uses MatDiagonalSet_IS with ADD_VALUES) */
493:     PetscPrintf(PETSC_COMM_WORLD,"Test MatShift\n");
494:     MatDuplicate(A,MAT_COPY_VALUES,&A2);
495:     MatDuplicate(B,MAT_COPY_VALUES,&B2);
496:     MatShift(A2,2.0);
497:     MatShift(B2,2.0);
498:     CheckMat(A2,B2,PETSC_FALSE,"MatShift");
499:     MatDestroy(&A2);
500:     MatDestroy(&B2);

502:     /* nonzero diag value is supported for square matrices only */
503:     TestMatZeroRows(A,B,PETSC_TRUE,is,diag);
504:   }
505:   TestMatZeroRows(A,B,squaretest,is,0.0);
506:   ISDestroy(&is);

508:   /* test MatTranspose */
509:   PetscPrintf(PETSC_COMM_WORLD,"Test MatTranspose\n");
510:   MatTranspose(A,MAT_INITIAL_MATRIX,&A2);
511:   MatTranspose(B,MAT_INITIAL_MATRIX,&B2);
512:   CheckMat(A2,B2,PETSC_FALSE,"initial matrix MatTranspose");

514:   MatTranspose(A,MAT_REUSE_MATRIX,&A2);
515:   CheckMat(A2,B2,PETSC_FALSE,"reuse matrix (not in place) MatTranspose");
516:   MatDestroy(&A2);

518:   MatDuplicate(A,MAT_COPY_VALUES,&A2);
519:   MatTranspose(A2,MAT_INPLACE_MATRIX,&A2);
520:   CheckMat(A2,B2,PETSC_FALSE,"reuse matrix (in place) MatTranspose");
521:   MatDestroy(&A2);

523:   MatTranspose(A,MAT_INITIAL_MATRIX,&A2);
524:   CheckMat(A2,B2,PETSC_FALSE,"reuse matrix (different type) MatTranspose");
525:   MatDestroy(&A2);
526:   MatDestroy(&B2);

528:   /* test MatISFixLocalEmpty */
529:   if (isaij) {
530:     PetscInt r[2];

532:     r[0] = 0;
533:     r[1] = PetscMin(m,n)-1;
534:     PetscPrintf(PETSC_COMM_WORLD,"Test MatISFixLocalEmpty\n");
535:     MatDuplicate(A,MAT_COPY_VALUES,&A2);

537:     MatISFixLocalEmpty(A2,PETSC_TRUE);
538:     MatAssemblyBegin(A2,MAT_FINAL_ASSEMBLY);
539:     MatAssemblyEnd(A2,MAT_FINAL_ASSEMBLY);
540:     CheckMat(A2,B,PETSC_FALSE,"MatISFixLocalEmpty (null)");

542:     MatZeroRows(A2,2,r,0.0,NULL,NULL);
543:     MatViewFromOptions(A2,NULL,"-fixempty_view");
544:     MatDuplicate(B,MAT_COPY_VALUES,&B2);
545:     MatZeroRows(B2,2,r,0.0,NULL,NULL);
546:     MatISFixLocalEmpty(A2,PETSC_TRUE);
547:     MatAssemblyBegin(A2,MAT_FINAL_ASSEMBLY);
548:     MatAssemblyEnd(A2,MAT_FINAL_ASSEMBLY);
549:     MatViewFromOptions(A2,NULL,"-fixempty_view");
550:     CheckMat(A2,B2,PETSC_FALSE,"MatISFixLocalEmpty (rows)");
551:     MatDestroy(&A2);

553:     MatDuplicate(A,MAT_COPY_VALUES,&A2);
554:     MatZeroRows(A2,2,r,0.0,NULL,NULL);
555:     MatTranspose(A2,MAT_INPLACE_MATRIX,&A2);
556:     MatTranspose(B2,MAT_INPLACE_MATRIX,&B2);
557:     MatViewFromOptions(A2,NULL,"-fixempty_view");
558:     MatISFixLocalEmpty(A2,PETSC_TRUE);
559:     MatAssemblyBegin(A2,MAT_FINAL_ASSEMBLY);
560:     MatAssemblyEnd(A2,MAT_FINAL_ASSEMBLY);
561:     MatViewFromOptions(A2,NULL,"-fixempty_view");
562:     CheckMat(A2,B2,PETSC_FALSE,"MatISFixLocalEmpty (cols)");

564:     MatDestroy(&A2);
565:     MatDestroy(&B2);

567:     if (squaretest) {
568:       MatDuplicate(A,MAT_COPY_VALUES,&A2);
569:       MatDuplicate(B,MAT_COPY_VALUES,&B2);
570:       MatZeroRowsColumns(A2,2,r,0.0,NULL,NULL);
571:       MatZeroRowsColumns(B2,2,r,0.0,NULL,NULL);
572:       MatViewFromOptions(A2,NULL,"-fixempty_view");
573:       MatISFixLocalEmpty(A2,PETSC_TRUE);
574:       MatAssemblyBegin(A2,MAT_FINAL_ASSEMBLY);
575:       MatAssemblyEnd(A2,MAT_FINAL_ASSEMBLY);
576:       MatViewFromOptions(A2,NULL,"-fixempty_view");
577:       CheckMat(A2,B2,PETSC_FALSE,"MatISFixLocalEmpty (rows+cols)");
578:       MatDestroy(&A2);
579:       MatDestroy(&B2);
580:     }

582:   }

584:   /* test MatInvertBlockDiagonal
585:        special cases for block-diagonal matrices */
586:   if (m == n) {
587:     ISLocalToGlobalMapping map;
588:     Mat                    Abd,Bbd;
589:     IS                     is,bis;
590:     const PetscScalar      *isbd,*aijbd;
591:     PetscScalar            *vals;
592:     const PetscInt         *sts,*idxs;
593:     PetscInt               *idxs2,diff,perm,nl,bs,st,en,in;
594:     PetscBool              ok;

596:     for (diff = 0; diff < 3; diff++) {
597:       for (perm = 0; perm < 3; perm++) {
598:         for (bs = 1; bs < 4; bs++) {
599:           PetscPrintf(PETSC_COMM_WORLD,"Test MatInvertBlockDiagonal blockdiag %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",n,diff,perm,bs);
600:           PetscMalloc1(bs*bs,&vals);
601:           MatGetOwnershipRanges(A,&sts);
602:           switch (diff) {
603:           case 1: /* inverted layout by processes */
604:             in = 1;
605:             st = sts[size - rank - 1];
606:             en = sts[size - rank];
607:             nl = en - st;
608:             break;
609:           case 2: /* round-robin layout */
610:             in = size;
611:             st = rank;
612:             nl = n/size;
613:             if (rank < n%size) nl++;
614:             break;
615:           default: /* same layout */
616:             in = 1;
617:             st = sts[rank];
618:             en = sts[rank + 1];
619:             nl = en - st;
620:             break;
621:           }
622:           ISCreateStride(PETSC_COMM_WORLD,nl,st,in,&is);
623:           ISGetLocalSize(is,&nl);
624:           ISGetIndices(is,&idxs);
625:           PetscMalloc1(nl,&idxs2);
626:           for (i=0;i<nl;i++) {
627:             switch (perm) { /* invert some of the indices */
628:             case 2:
629:               idxs2[i] = rank%2 ? idxs[i] : idxs[nl-i-1];
630:               break;
631:             case 1:
632:               idxs2[i] = rank%2 ? idxs[nl-i-1] : idxs[i];
633:               break;
634:             default:
635:               idxs2[i] = idxs[i];
636:               break;
637:             }
638:           }
639:           ISRestoreIndices(is,&idxs);
640:           ISCreateBlock(PETSC_COMM_WORLD,bs,nl,idxs2,PETSC_OWN_POINTER,&bis);
641:           ISLocalToGlobalMappingCreateIS(bis,&map);
642:           MatCreateIS(PETSC_COMM_WORLD,bs,PETSC_DECIDE,PETSC_DECIDE,bs*n,bs*n,map,map,&Abd);
643:           ISLocalToGlobalMappingDestroy(&map);
644:           MatISSetPreallocation(Abd,bs,NULL,0,NULL);
645:           for (i=0;i<nl;i++) {
646:             PetscInt b1,b2;

648:             for (b1=0;b1<bs;b1++) {
649:               for (b2=0;b2<bs;b2++) {
650:                 vals[b1*bs + b2] = i*bs*bs + b1*bs + b2 + 1 + (b1 == b2 ? 1.0 : 0);
651:               }
652:             }
653:             MatSetValuesBlockedLocal(Abd,1,&i,1,&i,vals,INSERT_VALUES);
654:           }
655:           MatAssemblyBegin(Abd,MAT_FINAL_ASSEMBLY);
656:           MatAssemblyEnd(Abd,MAT_FINAL_ASSEMBLY);
657:           MatConvert(Abd,MATAIJ,MAT_INITIAL_MATRIX,&Bbd);
658:           MatInvertBlockDiagonal(Abd,&isbd);
659:           MatInvertBlockDiagonal(Bbd,&aijbd);
660:           MatGetLocalSize(Bbd,&nl,NULL);
661:           ok   = PETSC_TRUE;
662:           for (i=0;i<nl/bs;i++) {
663:             PetscInt b1,b2;

665:             for (b1=0;b1<bs;b1++) {
666:               for (b2=0;b2<bs;b2++) {
667:                 if (PetscAbsScalar(isbd[i*bs*bs+b1*bs + b2]-aijbd[i*bs*bs+b1*bs + b2]) > PETSC_SMALL) ok = PETSC_FALSE;
668:                 if (!ok) {
669:                   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]));
670:                   break;
671:                 }
672:               }
673:               if (!ok) break;
674:             }
675:             if (!ok) break;
676:           }
677:           MatDestroy(&Abd);
678:           MatDestroy(&Bbd);
679:           PetscFree(vals);
680:           ISDestroy(&is);
681:           ISDestroy(&bis);
682:         }
683:       }
684:     }
685:   }
686:   /* free testing matrices */
687:   ISLocalToGlobalMappingDestroy(&cmap);
688:   ISLocalToGlobalMappingDestroy(&rmap);
689:   MatDestroy(&A);
690:   MatDestroy(&B);
691:   PetscFinalize();
692:   return 0;
693: }

695: PetscErrorCode CheckMat(Mat A, Mat B, PetscBool usemult, const char* func)
696: {
697:   Mat            Bcheck;
698:   PetscReal      error;

701:   if (!usemult && B) {
702:     PetscBool hasnorm;

704:     MatHasOperation(B,MATOP_NORM,&hasnorm);
705:     if (!hasnorm) usemult = PETSC_TRUE;
706:   }
707:   if (!usemult) {
708:     if (B) {
709:       MatType Btype;

711:       MatGetType(B,&Btype);
712:       MatConvert(A,Btype,MAT_INITIAL_MATRIX,&Bcheck);
713:     } else {
714:       MatConvert(A,MATAIJ,MAT_INITIAL_MATRIX,&Bcheck);
715:     }
716:     if (B) { /* if B is present, subtract it */
717:       MatAXPY(Bcheck,-1.,B,DIFFERENT_NONZERO_PATTERN);
718:     }
719:     MatNorm(Bcheck,NORM_INFINITY,&error);
720:     if (error > PETSC_SQRT_MACHINE_EPSILON) {
721:       ISLocalToGlobalMapping rl2g,cl2g;

723:       PetscObjectSetName((PetscObject)Bcheck,"Assembled Bcheck");
724:       MatView(Bcheck,NULL);
725:       if (B) {
726:         PetscObjectSetName((PetscObject)B,"Assembled AIJ");
727:         MatView(B,NULL);
728:         MatDestroy(&Bcheck);
729:         MatConvert(A,MATAIJ,MAT_INITIAL_MATRIX,&Bcheck);
730:         PetscObjectSetName((PetscObject)Bcheck,"Assembled IS");
731:         MatView(Bcheck,NULL);
732:       }
733:       MatDestroy(&Bcheck);
734:       PetscObjectSetName((PetscObject)A,"MatIS");
735:       MatView(A,NULL);
736:       MatGetLocalToGlobalMapping(A,&rl2g,&cl2g);
737:       ISLocalToGlobalMappingView(rl2g,NULL);
738:       ISLocalToGlobalMappingView(cl2g,NULL);
739:       SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"ERROR ON %s: %g",func,(double)error);
740:     }
741:     MatDestroy(&Bcheck);
742:   } else {
743:     PetscBool ok,okt;

745:     MatMultEqual(A,B,3,&ok);
746:     MatMultTransposeEqual(A,B,3,&okt);
748:   }
749:   return 0;
750: }

752: PetscErrorCode TestMatZeroRows(Mat A, Mat Afull, PetscBool squaretest, IS is, PetscScalar diag)
753: {
754:   Mat                    B,Bcheck,B2 = NULL,lB;
755:   Vec                    x = NULL, b = NULL, b2 = NULL;
756:   ISLocalToGlobalMapping l2gr,l2gc;
757:   PetscReal              error;
758:   char                   diagstr[16];
759:   const PetscInt         *idxs;
760:   PetscInt               rst,ren,i,n,N,d;
761:   PetscMPIInt            rank;
762:   PetscBool              miss,haszerorows;

765:   if (diag == 0.) {
766:     PetscStrcpy(diagstr,"zero");
767:   } else {
768:     PetscStrcpy(diagstr,"nonzero");
769:   }
770:   ISView(is,NULL);
771:   MatGetLocalToGlobalMapping(A,&l2gr,&l2gc);
772:   /* tests MatDuplicate and MatCopy */
773:   if (diag == 0.) {
774:     MatDuplicate(A,MAT_COPY_VALUES,&B);
775:   } else {
776:     MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,&B);
777:     MatCopy(A,B,SAME_NONZERO_PATTERN);
778:   }
779:   MatISGetLocalMat(B,&lB);
780:   MatHasOperation(lB,MATOP_ZERO_ROWS,&haszerorows);
781:   if (squaretest && haszerorows) {

783:     MatCreateVecs(B,&x,&b);
784:     MatDuplicate(B,MAT_COPY_VALUES,&B2);
785:     VecSetLocalToGlobalMapping(b,l2gr);
786:     VecSetLocalToGlobalMapping(x,l2gc);
787:     VecSetRandom(x,NULL);
788:     VecSetRandom(b,NULL);
789:     /* mimic b[is] = x[is] */
790:     VecDuplicate(b,&b2);
791:     VecSetLocalToGlobalMapping(b2,l2gr);
792:     VecCopy(b,b2);
793:     ISGetLocalSize(is,&n);
794:     ISGetIndices(is,&idxs);
795:     VecGetSize(x,&N);
796:     for (i=0;i<n;i++) {
797:       if (0 <= idxs[i] && idxs[i] < N) {
798:         VecSetValue(b2,idxs[i],diag,INSERT_VALUES);
799:         VecSetValue(x,idxs[i],1.,INSERT_VALUES);
800:       }
801:     }
802:     VecAssemblyBegin(b2);
803:     VecAssemblyEnd(b2);
804:     VecAssemblyBegin(x);
805:     VecAssemblyEnd(x);
806:     ISRestoreIndices(is,&idxs);
807:     /*  test ZeroRows on MATIS */
808:     PetscPrintf(PETSC_COMM_WORLD,"Test MatZeroRows (diag %s)\n",diagstr);
809:     MatZeroRowsIS(B,is,diag,x,b);
810:     PetscPrintf(PETSC_COMM_WORLD,"Test MatZeroRowsColumns (diag %s)\n",diagstr);
811:     MatZeroRowsColumnsIS(B2,is,diag,NULL,NULL);
812:   } else if (haszerorows) {
813:     /*  test ZeroRows on MATIS */
814:     PetscPrintf(PETSC_COMM_WORLD,"Test MatZeroRows (diag %s)\n",diagstr);
815:     MatZeroRowsIS(B,is,diag,NULL,NULL);
816:     b = b2 = x = NULL;
817:   } else {
818:     PetscPrintf(PETSC_COMM_WORLD,"Skipping MatZeroRows (diag %s)\n",diagstr);
819:     b = b2 = x = NULL;
820:   }

822:   if (squaretest && haszerorows) {
823:     VecAXPY(b2,-1.,b);
824:     VecNorm(b2,NORM_INFINITY,&error);
826:   }

828:   /* test MatMissingDiagonal */
829:   PetscPrintf(PETSC_COMM_WORLD,"Test MatMissingDiagonal\n");
830:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
831:   MatMissingDiagonal(B,&miss,&d);
832:   MatGetOwnershipRange(B,&rst,&ren);
833:   PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD);
834:   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);
835:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
836:   PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD);

838:   VecDestroy(&x);
839:   VecDestroy(&b);
840:   VecDestroy(&b2);

842:   /* check the result of ZeroRows with that from MPIAIJ routines
843:      assuming that MatConvert_IS_XAIJ and MatZeroRows_MPIAIJ work fine */
844:   if (haszerorows) {
845:     MatDuplicate(Afull,MAT_COPY_VALUES,&Bcheck);
846:     MatSetOption(Bcheck,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
847:     MatZeroRowsIS(Bcheck,is,diag,NULL,NULL);
848:     CheckMat(B,Bcheck,PETSC_FALSE,"Zerorows");
849:     MatDestroy(&Bcheck);
850:   }
851:   MatDestroy(&B);

853:   if (B2) { /* test MatZeroRowsColumns */
854:     MatDuplicate(Afull,MAT_COPY_VALUES,&B);
855:     MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
856:     MatZeroRowsColumnsIS(B,is,diag,NULL,NULL);
857:     CheckMat(B2,B,PETSC_FALSE,"MatZeroRowsColumns");
858:     MatDestroy(&B);
859:     MatDestroy(&B2);
860:   }
861:   return 0;
862: }

864: /*TEST

866:    test:
867:       args: -test_trans

869:    test:
870:       suffix: 2
871:       nsize: 4
872:       args: -matis_convert_local_nest -nr 3 -nc 4

874:    test:
875:       suffix: 3
876:       nsize: 5
877:       args: -m 11 -n 10 -matis_convert_local_nest -nr 2 -nc 1

879:    test:
880:       suffix: 4
881:       nsize: 6
882:       args: -m 9 -n 12 -test_trans -nr 2 -nc 7

884:    test:
885:       suffix: 5
886:       nsize: 6
887:       args: -m 12 -n 12 -test_trans -nr 3 -nc 1

889:    test:
890:       suffix: 6
891:       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap

893:    test:
894:       suffix: 7
895:       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap -permmap

897:    test:
898:       suffix: 8
899:       args: -m 12 -n 17 -test_trans -nr 2 -nc 3 -permmap

901:    test:
902:       suffix: 9
903:       nsize: 5
904:       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap

906:    test:
907:       suffix: 10
908:       nsize: 5
909:       args: -m 12 -n 12 -test_trans -nr 2 -nc 3 -diffmap -permmap

911:    test:
912:       suffix: vscat_default
913:       nsize: 5
914:       args: -m 12 -n 17 -test_trans -nr 2 -nc 3 -permmap
915:       output_file: output/ex23_11.out

917:    test:
918:       suffix: 12
919:       nsize: 3
920:       args: -m 12 -n 12 -symmetric -matis_localmat_type sbaij -test_trans -nr 2 -nc 3

922:    testset:
923:       output_file: output/ex23_13.out
924:       nsize: 3
925:       args: -m 12 -n 17 -test_trans -nr 2 -nc 3 -diffmap -permmap
926:       filter: grep -v "type:"
927:       test:
928:         suffix: baij
929:         args: -matis_localmat_type baij
930:       test:
931:         requires: viennacl
932:         suffix: viennacl
933:         args: -matis_localmat_type aijviennacl
934:       test:
935:         requires: cuda
936:         suffix: cusparse
937:         args: -matis_localmat_type aijcusparse

939:    test:
940:       suffix: negrep
941:       nsize: {{1 3}separate output}
942:       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}

944: TEST*/