xref: /petsc/src/mat/tests/ex115.c (revision c4762a1b19cd2af06abeed90e8f9d34fb975dd94)
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,256,&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);CHKERRQ(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);
70       if (test_offproc) {
71         const PetscInt *ranges;
72         PetscMPIInt    rank;
73 
74         ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(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 
170     /* MatPtAP_HYPRE_HYPRE -> output is a MatHYPRE */
171     ierr = MatPtAP(C,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
172     ierr = MatPtAP(C,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
173     ierr = MatAXPY(hP,-1.,pP,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
174     ierr = MatHasOperation(hP,MATOP_NORM,&flg);CHKERRQ(ierr);
175     if (!flg) { /* TODO add MatNorm_HYPRE */
176       ierr = MatConvert(hP,MATAIJ,MAT_INPLACE_MATRIX,&hP);CHKERRQ(ierr);
177     }
178     ierr = MatNorm(hP,NORM_INFINITY,&err);CHKERRQ(ierr);
179     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatPtAP %g %g",err,norm);
180     ierr = MatDestroy(&hP);CHKERRQ(ierr);
181 
182     /* MatPtAP_AIJ_HYPRE -> output can be decided at runtime with -matptap_hypre_outtype */
183     ierr = MatPtAP(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
184     ierr = MatPtAP(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&hP);CHKERRQ(ierr);
185     ierr = MatAXPY(hP,-1.,pP,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
186     ierr = MatHasOperation(hP,MATOP_NORM,&flg);CHKERRQ(ierr);
187     if (!flg) { /* TODO add MatNorm_HYPRE */
188       ierr = MatConvert(hP,MATAIJ,MAT_INPLACE_MATRIX,&hP);CHKERRQ(ierr);
189     }
190     ierr = MatNorm(hP,NORM_INFINITY,&err);CHKERRQ(ierr);
191     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatPtAP mixed %g %g",err,norm);
192     ierr = MatDestroy(&hP);CHKERRQ(ierr);
193 
194     ierr = MatDestroy(&pP);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 
211     /* MatMatMult_HYPRE_HYPRE -> output is a MatHYPRE */
212     ierr = MatMatMult(C,D,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&CD);CHKERRQ(ierr);
213     ierr = MatMatMult(C,D,MAT_REUSE_MATRIX,PETSC_DEFAULT,&CD);CHKERRQ(ierr);
214     ierr = MatAXPY(CD,-1.,pAB,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
215     ierr = MatHasOperation(CD,MATOP_NORM,&flg);CHKERRQ(ierr);
216     if (!flg) { /* TODO add MatNorm_HYPRE */
217       ierr = MatConvert(CD,MATAIJ,MAT_INPLACE_MATRIX,&CD);CHKERRQ(ierr);
218     }
219     ierr = MatNorm(CD,NORM_INFINITY,&err);CHKERRQ(ierr);
220     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMatMult %g %g",err,norm);
221     ierr = MatDestroy(&C);CHKERRQ(ierr);
222     ierr = MatDestroy(&D);CHKERRQ(ierr);
223     ierr = MatDestroy(&CD);CHKERRQ(ierr);
224     ierr = MatDestroy(&pAB);CHKERRQ(ierr);
225 
226     /* When configured with HYPRE, MatMatMatMult is available for the triplet transpose(aij)-aij-aij */
227     ierr = MatCreateTranspose(A,&C);CHKERRQ(ierr);
228     ierr = MatMatMatMult(C,A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&CAB);CHKERRQ(ierr);
229     ierr = MatDestroy(&C);CHKERRQ(ierr);
230     ierr = MatTranspose(A,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);
231     ierr = MatMatMult(C,A,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&D);CHKERRQ(ierr);
232     ierr = MatDestroy(&C);CHKERRQ(ierr);
233     ierr = MatMatMult(D,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C);CHKERRQ(ierr);
234     ierr = MatNorm(C,NORM_INFINITY,&norm);CHKERRQ(ierr);
235     ierr = MatAXPY(C,-1.,CAB,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
236     ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
237     if (err/norm > PETSC_SMALL) SETERRQ2(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatMatMatMult %g %g",err,norm);
238     ierr = MatDestroy(&C);CHKERRQ(ierr);
239     ierr = MatDestroy(&D);CHKERRQ(ierr);
240     ierr = MatDestroy(&CAB);CHKERRQ(ierr);
241     ierr = MatDestroy(&B);CHKERRQ(ierr);
242   }
243 
244   /* Check MatView */
245   ierr = MatViewFromOptions(A,NULL,"-view_A");CHKERRQ(ierr);
246   ierr = MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);
247   ierr = MatViewFromOptions(B,NULL,"-view_B");CHKERRQ(ierr);
248 
249   /* Check MatDuplicate/MatCopy */
250   for (j=0;j<3;j++) {
251     MatDuplicateOption dop;
252 
253     dop = MAT_COPY_VALUES;
254     if (j==1) dop = MAT_DO_NOT_COPY_VALUES;
255     if (j==2) dop = MAT_SHARE_NONZERO_PATTERN;
256 
257     for (i=0;i<3;i++) {
258       MatStructure str;
259 
260       ierr = PetscPrintf(PETSC_COMM_WORLD,"Dup/Copy tests: %D %D\n",j,i);CHKERRQ(ierr);
261 
262       str = DIFFERENT_NONZERO_PATTERN;
263       if (i==1) str = SAME_NONZERO_PATTERN;
264       if (i==2) str = SUBSET_NONZERO_PATTERN;
265 
266       ierr = MatDuplicate(A,dop,&C);CHKERRQ(ierr);
267       ierr = MatDuplicate(B,dop,&D);CHKERRQ(ierr);
268       if (dop != MAT_COPY_VALUES) {
269         ierr = MatCopy(A,C,str);CHKERRQ(ierr);
270         ierr = MatCopy(B,D,str);CHKERRQ(ierr);
271       }
272       /* AXPY with AIJ and HYPRE */
273       ierr = MatAXPY(C,-1.0,D,str);CHKERRQ(ierr);
274       ierr = MatNorm(C,NORM_INFINITY,&err);CHKERRQ(ierr);
275       if (err > PETSC_SMALL) {
276         ierr = MatViewFromOptions(A,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
277         ierr = MatViewFromOptions(B,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
278         ierr = MatViewFromOptions(C,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
279         SETERRQ3(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error test 1 MatDuplicate/MatCopy %g (%D,%D)",err,j,i);
280       }
281       /* AXPY with HYPRE and HYPRE */
282       ierr = MatAXPY(D,-1.0,B,str);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(D,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
287         SETERRQ3(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error test 2 MatDuplicate/MatCopy %g (%D,%D)",err,j,i);
288       }
289       /* Copy from HYPRE to AIJ */
290       ierr = MatCopy(B,C,str);CHKERRQ(ierr);
291       /* Copy from AIJ to HYPRE */
292       ierr = MatCopy(A,D,str);CHKERRQ(ierr);
293       /* AXPY with HYPRE and AIJ */
294       ierr = MatAXPY(D,-1.0,C,str);CHKERRQ(ierr);
295       ierr = MatHasOperation(D,MATOP_NORM,&flg);CHKERRQ(ierr);
296       if (!flg) { /* TODO add MatNorm_HYPRE */
297         ierr = MatConvert(D,MATAIJ,MAT_INPLACE_MATRIX,&D);CHKERRQ(ierr);
298       }
299       ierr = MatNorm(D,NORM_INFINITY,&err);CHKERRQ(ierr);
300       if (err > PETSC_SMALL) {
301         ierr = MatViewFromOptions(A,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
302         ierr = MatViewFromOptions(C,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
303         ierr = MatViewFromOptions(D,NULL,"-view_duplicate_diff");CHKERRQ(ierr);
304         SETERRQ3(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error test 3 MatDuplicate/MatCopy %g (%D,%D)",err,j,i);
305       }
306       ierr = MatDestroy(&C);CHKERRQ(ierr);
307       ierr = MatDestroy(&D);CHKERRQ(ierr);
308     }
309   }
310   ierr = MatDestroy(&B);CHKERRQ(ierr);
311 
312   ierr = MatHasCongruentLayouts(A,&flg);CHKERRQ(ierr);
313   if (flg) {
314     Vec y,y2;
315 
316     ierr = MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);
317     ierr = MatCreateVecs(A,NULL,&y);CHKERRQ(ierr);
318     ierr = MatCreateVecs(B,NULL,&y2);CHKERRQ(ierr);
319     ierr = MatGetDiagonal(A,y);CHKERRQ(ierr);
320     ierr = MatGetDiagonal(B,y2);CHKERRQ(ierr);
321     ierr = VecAXPY(y2,-1.0,y);CHKERRQ(ierr);
322     ierr = VecNorm(y2,NORM_INFINITY,&err);CHKERRQ(ierr);
323     if (err > PETSC_SMALL) {
324       ierr = VecViewFromOptions(y,NULL,"-view_diagonal_diff");CHKERRQ(ierr);
325       ierr = VecViewFromOptions(y2,NULL,"-view_diagonal_diff");CHKERRQ(ierr);
326       SETERRQ1(PetscObjectComm((PetscObject)A),PETSC_ERR_PLIB,"Error MatGetDiagonal %g",err);
327     }
328     ierr = MatDestroy(&B);CHKERRQ(ierr);
329     ierr = VecDestroy(&y);CHKERRQ(ierr);
330     ierr = VecDestroy(&y2);CHKERRQ(ierr);
331   }
332 
333   ierr = MatDestroy(&A);CHKERRQ(ierr);
334 
335   ierr = PetscFinalize();
336   return ierr;
337 }
338 
339 
340 /*TEST
341 
342    build:
343       requires: hypre
344 
345    test:
346       suffix: 1
347       args: -N 11 -M 11
348       output_file: output/ex115_1.out
349 
350    test:
351       suffix: 2
352       nsize: 3
353       args: -N 13 -M 13 -matmatmult_via hypre
354       output_file: output/ex115_1.out
355 
356    test:
357       suffix: 3
358       nsize: 4
359       args: -M 13 -N 7 -matmatmult_via hypre
360       output_file: output/ex115_1.out
361 
362    test:
363       suffix: 4
364       nsize: 2
365       args: -M 12 -N 19
366       output_file: output/ex115_1.out
367 
368    test:
369       suffix: 5
370       nsize: 3
371       args: -M 13 -N 13 -matptap_via hypre -matptap_hypre_outtype hypre
372       output_file: output/ex115_1.out
373 
374    test:
375       suffix: 6
376       nsize: 3
377       args: -M 12 -N 19 -test_offproc
378       output_file: output/ex115_1.out
379 
380    test:
381       suffix: 7
382       nsize: 3
383       args: -M 19 -N 12 -test_offproc -view_B ::ascii_info_detail
384       output_file: output/ex115_7.out
385 
386    test:
387       suffix: 8
388       nsize: 3
389       args: -M 1 -N 12 -test_offproc
390       output_file: output/ex115_1.out
391 
392    test:
393       suffix: 9
394       nsize: 3
395       args: -M 1 -N 2 -test_offproc
396       output_file: output/ex115_1.out
397 
398 TEST*/
399