xref: /petsc/src/mat/tests/ex115.c (revision 030f984af8d8bb4c203755d35bded3c05b3d83ce)
1 
2 static char help[] = "Tests MatHYPRE\n";
3 
4 #include <petscmathypre.h>
5 
6 int main(int argc,char **args)
7 {
8   Mat                A,B,C,D;
9   Mat                pAB,CD,CAB;
10   hypre_ParCSRMatrix *parcsr;
11   PetscReal          err;
12   PetscInt           i,j,N = 6, M = 6;
13   PetscErrorCode     ierr;
14   PetscBool          flg,testptap = PETSC_TRUE,testmatmatmult = PETSC_TRUE;
15   PetscReal          norm;
16   char               file[256];
17 
18   ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
19   ierr = PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg);CHKERRQ(ierr);
20 #if defined(PETSC_USE_COMPLEX)
21   testptap = PETSC_FALSE;
22   testmatmatmult = PETSC_FALSE;
23   ierr = PetscOptionsInsertString(NULL,"-options_left 0");CHKERRQ(ierr);
24 #endif
25   ierr = PetscOptionsGetBool(NULL,NULL,"-ptap",&testptap,NULL);CHKERRQ(ierr);
26   ierr = PetscOptionsGetBool(NULL,NULL,"-matmatmult",&testmatmatmult,NULL);CHKERRQ(ierr);
27   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
28   if (!flg) { /* Create a matrix and test MatSetValues */
29     PetscMPIInt size;
30 
31     ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRMPI(ierr);
32     ierr = PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);CHKERRQ(ierr);
33     ierr = PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);CHKERRQ(ierr);
34     ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr);
35     ierr = MatSetType(A,MATAIJ);CHKERRQ(ierr);
36     ierr = MatSeqAIJSetPreallocation(A,9,NULL);CHKERRQ(ierr);
37     ierr = MatMPIAIJSetPreallocation(A,9,NULL,9,NULL);CHKERRQ(ierr);
38     ierr = MatCreate(PETSC_COMM_WORLD,&B);CHKERRQ(ierr);
39     ierr = MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr);
40     ierr = MatSetType(B,MATHYPRE);CHKERRQ(ierr);
41     if (M == N) {
42       ierr = MatHYPRESetPreallocation(B,9,NULL,9,NULL);CHKERRQ(ierr);
43     } else {
44       ierr = MatHYPRESetPreallocation(B,6,NULL,6,NULL);CHKERRQ(ierr);
45     }
46     if (M == N) {
47       for (i=0; i<M; i++) {
48         PetscInt    cols[] = {0,1,2,3,4,5};
49         PetscScalar vals[] = {0,1./size,2./size,3./size,4./size,5./size};
50         for (j=i-2; j<i+1; j++) {
51           if (j >= N) {
52             ierr = MatSetValue(A,i,N-1,(1.*j*N+i)/(3.*N*size),ADD_VALUES);CHKERRQ(ierr);
53             ierr = MatSetValue(B,i,N-1,(1.*j*N+i)/(3.*N*size),ADD_VALUES);CHKERRQ(ierr);
54           } else if (i > j) {
55             ierr = MatSetValue(A,i,PetscMin(j,N-1),(1.*j*N+i)/(2.*N*size),ADD_VALUES);CHKERRQ(ierr);
56             ierr = MatSetValue(B,i,PetscMin(j,N-1),(1.*j*N+i)/(2.*N*size),ADD_VALUES);CHKERRQ(ierr);
57           } else {
58             ierr = MatSetValue(A,i,PetscMin(j,N-1),-1.-(1.*j*N+i)/(4.*N*size),ADD_VALUES);CHKERRQ(ierr);
59             ierr = MatSetValue(B,i,PetscMin(j,N-1),-1.-(1.*j*N+i)/(4.*N*size),ADD_VALUES);CHKERRQ(ierr);
60           }
61         }
62         ierr = MatSetValues(A,1,&i,PetscMin(6,N),cols,vals,ADD_VALUES);CHKERRQ(ierr);
63         ierr = MatSetValues(B,1,&i,PetscMin(6,N),cols,vals,ADD_VALUES);CHKERRQ(ierr);
64       }
65     } else {
66       PetscInt  rows[2];
67       PetscBool test_offproc = PETSC_FALSE;
68 
69       ierr = PetscOptionsGetBool(NULL,NULL,"-test_offproc",&test_offproc,NULL);CHKERRQ(ierr);
70       if (test_offproc) {
71         const PetscInt *ranges;
72         PetscMPIInt    rank;
73 
74         ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr);
75         ierr = MatGetOwnershipRanges(A,&ranges);CHKERRQ(ierr);
76         rows[0] = ranges[(rank+1)%size];
77         rows[1] = ranges[(rank+1)%size + 1];
78       } else {
79         ierr = MatGetOwnershipRange(A,&rows[0],&rows[1]);CHKERRQ(ierr);
80       }
81       for (i=rows[0];i<rows[1];i++) {
82         PetscInt    cols[] = {0,1,2,3,4,5};
83         PetscScalar vals[] = {-1,1,-2,2,-3,3};
84 
85         ierr = MatSetValues(A,1,&i,PetscMin(6,N),cols,vals,INSERT_VALUES);CHKERRQ(ierr);
86         ierr = MatSetValues(B,1,&i,PetscMin(6,N),cols,vals,INSERT_VALUES);CHKERRQ(ierr);
87       }
88     }
89     /* MAT_FLUSH_ASSEMBLY currently not supported */
90     ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
91     ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
92     ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
93     ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
94 
95 #if defined(PETSC_USE_COMPLEX)
96     /* make the matrix imaginary */
97     ierr = MatScale(A,PETSC_i);CHKERRQ(ierr);
98     ierr = MatScale(B,PETSC_i);CHKERRQ(ierr);
99 #endif
100 
101     /* MatAXPY further exercises MatSetValues_HYPRE */
102     ierr = MatAXPY(B,-1.,A,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
103     ierr = MatConvert(B,MATMPIAIJ,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);
104     ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
105     if (err > PETSC_SMALL) SETERRQ1(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatSetValues %g",err);
106     ierr = MatDestroy(&B);CHKERRQ(ierr);
107     ierr = MatDestroy(&C);CHKERRQ(ierr);
108   } else {
109     PetscViewer viewer;
110 
111     ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer);CHKERRQ(ierr);
112     ierr = MatSetFromOptions(A);CHKERRQ(ierr);
113     ierr = MatLoad(A,viewer);CHKERRQ(ierr);
114     ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
115     ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
116   }
117 
118   /* check conversion routines */
119   ierr = MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);
120   ierr = MatConvert(A,MATHYPRE,MAT_REUSE_MATRIX,&B);CHKERRQ(ierr);
121   ierr = MatConvert(B,MATIS,MAT_INITIAL_MATRIX,&D);CHKERRQ(ierr);
122   ierr = MatConvert(B,MATIS,MAT_REUSE_MATRIX,&D);CHKERRQ(ierr);
123   ierr = MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);
124   ierr = MatConvert(B,MATAIJ,MAT_REUSE_MATRIX,&C);CHKERRQ(ierr);
125   ierr = MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN);CHKERRQ(ierr);
126   ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
127   if (err > PETSC_SMALL) SETERRQ1(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error Mat AIJ %g",err);
128   ierr = MatDestroy(&C);CHKERRQ(ierr);
129   ierr = MatConvert(D,MATAIJ,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);
130   ierr = MatAXPY(C,-1.,A,SAME_NONZERO_PATTERN);CHKERRQ(ierr);
131   ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
132   if (err > PETSC_SMALL) SETERRQ1(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error Mat IS %g",err);
133   ierr = MatDestroy(&C);CHKERRQ(ierr);
134   ierr = MatDestroy(&D);CHKERRQ(ierr);
135 
136   /* check MatCreateFromParCSR */
137   ierr = MatHYPREGetParCSR(B,&parcsr);CHKERRQ(ierr);
138   ierr = MatCreateFromParCSR(parcsr,MATAIJ,PETSC_COPY_VALUES,&D);CHKERRQ(ierr);
139   ierr = MatDestroy(&D);CHKERRQ(ierr);
140   ierr = MatCreateFromParCSR(parcsr,MATHYPRE,PETSC_USE_POINTER,&C);CHKERRQ(ierr);
141 
142   /* check MatMult operations */
143   ierr = MatMultEqual(A,B,4,&flg);CHKERRQ(ierr);
144   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMult B");
145   ierr = MatMultEqual(A,C,4,&flg);CHKERRQ(ierr);
146   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMult C");
147   ierr = MatMultAddEqual(A,B,4,&flg);CHKERRQ(ierr);
148   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMultAdd B");
149   ierr = MatMultAddEqual(A,C,4,&flg);CHKERRQ(ierr);
150   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMultAdd C");
151   ierr = MatMultTransposeEqual(A,B,4,&flg);CHKERRQ(ierr);
152   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMultTranspose B");
153   ierr = MatMultTransposeEqual(A,C,4,&flg);CHKERRQ(ierr);
154   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMultTranspose C");
155   ierr = MatMultTransposeAddEqual(A,B,4,&flg);CHKERRQ(ierr);
156   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMultTransposeAdd B");
157   ierr = MatMultTransposeAddEqual(A,C,4,&flg);CHKERRQ(ierr);
158   if (!flg) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMultTransposeAdd C");
159 
160   /* check PtAP */
161   if (testptap && M == N) {
162     Mat pP,hP;
163 
164     /* PETSc MatPtAP -> output is a MatAIJ
165        It uses HYPRE functions when -matptap_via hypre is specified at command line */
166     ierr = MatPtAP(A,A,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&pP);CHKERRQ(ierr);
167     ierr = MatPtAP(A,A,MAT_REUSE_MATRIX,PETSC_DEFAULT,&pP);CHKERRQ(ierr);
168     ierr = MatNorm(pP,NORM_INFINITY,&norm);CHKERRQ(ierr);
169     ierr = MatPtAPMultEqual(A,A,pP,10,&flg);CHKERRQ(ierr);
170     if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatPtAP_MatAIJ");
171 
172     /* MatPtAP_HYPRE_HYPRE -> output is a MatHYPRE */
173     ierr = MatPtAP(C,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
174     ierr = MatPtAP(C,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
175     ierr = MatPtAPMultEqual(C,B,hP,10,&flg);CHKERRQ(ierr);
176     if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatPtAP_HYPRE_HYPRE");
177 
178     /* Test MatAXPY_Basic() */
179     ierr = MatAXPY(hP,-1.,pP,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
180     ierr = MatHasOperation(hP,MATOP_NORM,&flg);CHKERRQ(ierr);
181     if (!flg) { /* TODO add MatNorm_HYPRE */
182       ierr = MatConvert(hP,MATAIJ,MAT_INPLACE_MATRIX,&hP);CHKERRQ(ierr);
183     }
184     ierr = MatNorm(hP,NORM_INFINITY,&err);CHKERRQ(ierr);
185     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)hP),PETSC_ERR_PLIB,"Error MatPtAP %g %g",err,norm);
186     ierr = MatDestroy(&pP);CHKERRQ(ierr);
187     ierr = MatDestroy(&hP);CHKERRQ(ierr);
188 
189     /* MatPtAP_AIJ_HYPRE -> output can be decided at runtime with -matptap_hypre_outtype */
190     ierr = MatPtAP(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
191     ierr = MatPtAP(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
192     ierr = MatPtAPMultEqual(A,B,hP,10,&flg);CHKERRQ(ierr);
193     if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatPtAP_AIJ_HYPRE");
194     ierr = MatDestroy(&hP);CHKERRQ(ierr);
195   }
196   ierr = MatDestroy(&C);CHKERRQ(ierr);
197   ierr = MatDestroy(&B);CHKERRQ(ierr);
198 
199   /* check MatMatMult */
200   if (testmatmatmult) {
201     ierr = MatTranspose(A,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);
202     ierr = MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);
203     ierr = MatConvert(B,MATHYPRE,MAT_INITIAL_MATRIX,&D);CHKERRQ(ierr);
204 
205     /* PETSc MatMatMult -> output is a MatAIJ
206        It uses HYPRE functions when -matmatmult_via hypre is specified at command line */
207     ierr = MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&pAB);CHKERRQ(ierr);
208     ierr = MatMatMult(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&pAB);CHKERRQ(ierr);
209     ierr = MatNorm(pAB,NORM_INFINITY,&norm);CHKERRQ(ierr);
210     ierr = MatMatMultEqual(A,B,pAB,10,&flg);CHKERRQ(ierr);
211     if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatMatMult_AIJ_AIJ");
212 
213     /* MatMatMult_HYPRE_HYPRE -> output is a MatHYPRE */
214     ierr = MatMatMult(C,D,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&CD);CHKERRQ(ierr);
215     ierr = MatMatMult(C,D,MAT_REUSE_MATRIX,PETSC_DEFAULT,&CD);CHKERRQ(ierr);
216     ierr = MatMatMultEqual(C,D,CD,10,&flg);CHKERRQ(ierr);
217     if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatMatMult_HYPRE_HYPRE");
218 
219     /* Test MatAXPY_Basic() */
220     ierr = MatAXPY(CD,-1.,pAB,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
221 
222     ierr = MatHasOperation(CD,MATOP_NORM,&flg);CHKERRQ(ierr);
223     if (!flg) { /* TODO add MatNorm_HYPRE */
224       ierr = MatConvert(CD,MATAIJ,MAT_INPLACE_MATRIX,&CD);CHKERRQ(ierr);
225     }
226     ierr = MatNorm(CD,NORM_INFINITY,&err);CHKERRQ(ierr);
227     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)CD),PETSC_ERR_PLIB,"Error MatMatMult %g %g",err,norm);
228 
229     ierr = MatDestroy(&C);CHKERRQ(ierr);
230     ierr = MatDestroy(&D);CHKERRQ(ierr);
231     ierr = MatDestroy(&pAB);CHKERRQ(ierr);
232     ierr = MatDestroy(&CD);CHKERRQ(ierr);
233 
234     /* When configured with HYPRE, MatMatMatMult is available for the triplet transpose(aij)-aij-aij */
235     ierr = MatCreateTranspose(A,&C);CHKERRQ(ierr);
236     ierr = MatMatMatMult(C,A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&CAB);CHKERRQ(ierr);
237     ierr = MatDestroy(&C);CHKERRQ(ierr);
238     ierr = MatTranspose(A,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);
239     ierr = MatMatMult(C,A,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&D);CHKERRQ(ierr);
240     ierr = MatDestroy(&C);CHKERRQ(ierr);
241     ierr = MatMatMult(D,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C);CHKERRQ(ierr);
242     ierr = MatNorm(C,NORM_INFINITY,&norm);CHKERRQ(ierr);
243     ierr = MatAXPY(C,-1.,CAB,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
244     ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
245     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMatMatMult %g %g",err,norm);
246     ierr = MatDestroy(&C);CHKERRQ(ierr);
247     ierr = MatDestroy(&D);CHKERRQ(ierr);
248     ierr = MatDestroy(&CAB);CHKERRQ(ierr);
249     ierr = MatDestroy(&B);CHKERRQ(ierr);
250   }
251 
252   /* Check MatView */
253   ierr = MatViewFromOptions(A,NULL,"-view_A");CHKERRQ(ierr);
254   ierr = MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);
255   ierr = MatViewFromOptions(B,NULL,"-view_B");CHKERRQ(ierr);
256 
257   /* Check MatDuplicate/MatCopy */
258   for (j=0;j<3;j++) {
259     MatDuplicateOption dop;
260 
261     dop = MAT_COPY_VALUES;
262     if (j==1) dop = MAT_DO_NOT_COPY_VALUES;
263     if (j==2) dop = MAT_SHARE_NONZERO_PATTERN;
264 
265     for (i=0;i<3;i++) {
266       MatStructure str;
267 
268       ierr = PetscPrintf(PETSC_COMM_WORLD,"Dup/Copy tests: %D %D\n",j,i);CHKERRQ(ierr);
269 
270       str = DIFFERENT_NONZERO_PATTERN;
271       if (i==1) str = SAME_NONZERO_PATTERN;
272       if (i==2) str = SUBSET_NONZERO_PATTERN;
273 
274       ierr = MatDuplicate(A,dop,&C);CHKERRQ(ierr);
275       ierr = MatDuplicate(B,dop,&D);CHKERRQ(ierr);
276       if (dop != MAT_COPY_VALUES) {
277         ierr = MatCopy(A,C,str);CHKERRQ(ierr);
278         ierr = MatCopy(B,D,str);CHKERRQ(ierr);
279       }
280       /* AXPY with AIJ and HYPRE */
281       ierr = MatAXPY(C,-1.0,D,str);CHKERRQ(ierr);
282       ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
283       if (err > PETSC_SMALL) {
284         ierr = MatViewFromOptions(A,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
285         ierr = MatViewFromOptions(B,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
286         ierr = MatViewFromOptions(C,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
287         SETERRQ3(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error test 1 MatDuplicate/MatCopy %g (%D,%D)",err,j,i);
288       }
289       /* AXPY with HYPRE and HYPRE */
290       ierr = MatAXPY(D,-1.0,B,str);CHKERRQ(ierr);
291       if (err > PETSC_SMALL) {
292         ierr = MatViewFromOptions(A,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
293         ierr = MatViewFromOptions(B,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
294         ierr = MatViewFromOptions(D,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
295         SETERRQ3(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error test 2 MatDuplicate/MatCopy %g (%D,%D)",err,j,i);
296       }
297       /* Copy from HYPRE to AIJ */
298       ierr = MatCopy(B,C,str);CHKERRQ(ierr);
299       /* Copy from AIJ to HYPRE */
300       ierr = MatCopy(A,D,str);CHKERRQ(ierr);
301       /* AXPY with HYPRE and AIJ */
302       ierr = MatAXPY(D,-1.0,C,str);CHKERRQ(ierr);
303       ierr = MatHasOperation(D,MATOP_NORM,&flg);CHKERRQ(ierr);
304       if (!flg) { /* TODO add MatNorm_HYPRE */
305         ierr = MatConvert(D,MATAIJ,MAT_INPLACE_MATRIX,&D);CHKERRQ(ierr);
306       }
307       ierr = MatNorm(D,NORM_INFINITY,&err);CHKERRQ(ierr);
308       if (err > PETSC_SMALL) {
309         ierr = MatViewFromOptions(A,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
310         ierr = MatViewFromOptions(C,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
311         ierr = MatViewFromOptions(D,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
312         SETERRQ3(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error test 3 MatDuplicate/MatCopy %g (%D,%D)",err,j,i);
313       }
314       ierr = MatDestroy(&C);CHKERRQ(ierr);
315       ierr = MatDestroy(&D);CHKERRQ(ierr);
316     }
317   }
318   ierr = MatDestroy(&B);CHKERRQ(ierr);
319 
320   ierr = MatHasCongruentLayouts(A,&flg);CHKERRQ(ierr);
321   if (flg) {
322     Vec y,y2;
323 
324     ierr = MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);
325     ierr = MatCreateVecs(A,NULL,&y);CHKERRQ(ierr);
326     ierr = MatCreateVecs(B,NULL,&y2);CHKERRQ(ierr);
327     ierr = MatGetDiagonal(A,y);CHKERRQ(ierr);
328     ierr = MatGetDiagonal(B,y2);CHKERRQ(ierr);
329     ierr = VecAXPY(y2,-1.0,y);CHKERRQ(ierr);
330     ierr = VecNorm(y2,NORM_INFINITY,&err);CHKERRQ(ierr);
331     if (err > PETSC_SMALL) {
332       ierr = VecViewFromOptions(y,NULL,"-view_diagonal_diff");CHKERRQ(ierr);
333       ierr = VecViewFromOptions(y2,NULL,"-view_diagonal_diff");CHKERRQ(ierr);
334       SETERRQ1(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatGetDiagonal %g",err);
335     }
336     ierr = MatDestroy(&B);CHKERRQ(ierr);
337     ierr = VecDestroy(&y);CHKERRQ(ierr);
338     ierr = VecDestroy(&y2);CHKERRQ(ierr);
339   }
340 
341   ierr = MatDestroy(&A);CHKERRQ(ierr);
342 
343   ierr = PetscFinalize();
344   return ierr;
345 }
346 
347 /*TEST
348 
349    build:
350       requires: hypre
351 
352    test:
353       suffix: 1
354       args: -N 11 -M 11
355       output_file: output/ex115_1.out
356 
357    test:
358       suffix: 2
359       nsize: 3
360       args: -N 13 -M 13 -matmatmult_via hypre
361       output_file: output/ex115_1.out
362 
363    test:
364       suffix: 3
365       nsize: 4
366       args: -M 13 -N 7 -matmatmult_via hypre
367       output_file: output/ex115_1.out
368 
369    test:
370       suffix: 4
371       nsize: 2
372       args: -M 12 -N 19
373       output_file: output/ex115_1.out
374 
375    test:
376       suffix: 5
377       nsize: 3
378       args: -M 13 -N 13 -matptap_via hypre -matptap_hypre_outtype hypre
379       output_file: output/ex115_1.out
380 
381    test:
382       suffix: 6
383       nsize: 3
384       args: -M 12 -N 19 -test_offproc
385       output_file: output/ex115_1.out
386 
387    test:
388       suffix: 7
389       nsize: 3
390       args: -M 19 -N 12 -test_offproc -view_B ::ascii_info_detail
391       output_file: output/ex115_7.out
392 
393    test:
394       suffix: 8
395       nsize: 3
396       args: -M 1 -N 12 -test_offproc
397       output_file: output/ex115_1.out
398 
399    test:
400       suffix: 9
401       nsize: 3
402       args: -M 1 -N 2 -test_offproc
403       output_file: output/ex115_1.out
404 
405 TEST*/
406