xref: /petsc/src/mat/tests/ex115.c (revision f9dc6c1fa0305d45f2b8340c2fda6c66624b5bfb)
1c4762a1bSJed Brown 
2c4762a1bSJed Brown static char help[] = "Tests MatHYPRE\n";
3c4762a1bSJed Brown 
4c4762a1bSJed Brown #include <petscmathypre.h>
5c4762a1bSJed Brown 
6d71ae5a4SJacob Faibussowitsch int main(int argc, char **args)
7d71ae5a4SJacob Faibussowitsch {
8c4762a1bSJed Brown   Mat                 A, B, C, D;
9c4762a1bSJed Brown   Mat                 pAB, CD, CAB;
10c4762a1bSJed Brown   hypre_ParCSRMatrix *parcsr;
11c4762a1bSJed Brown   PetscReal           err;
12c4762a1bSJed Brown   PetscInt            i, j, N = 6, M = 6;
13c4762a1bSJed Brown   PetscBool           flg, testptap = PETSC_TRUE, testmatmatmult = PETSC_TRUE;
14c4762a1bSJed Brown   PetscReal           norm;
15c4762a1bSJed Brown   char                file[256];
16c4762a1bSJed Brown 
17327415f7SBarry Smith   PetscFunctionBeginUser;
189566063dSJacob Faibussowitsch   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
199566063dSJacob Faibussowitsch   PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
20c4762a1bSJed Brown #if defined(PETSC_USE_COMPLEX)
21c4762a1bSJed Brown   testptap       = PETSC_FALSE;
22c4762a1bSJed Brown   testmatmatmult = PETSC_FALSE;
239566063dSJacob Faibussowitsch   PetscCall(PetscOptionsInsertString(NULL, "-options_left 0"));
24c4762a1bSJed Brown #endif
259566063dSJacob Faibussowitsch   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ptap", &testptap, NULL));
269566063dSJacob Faibussowitsch   PetscCall(PetscOptionsGetBool(NULL, NULL, "-matmatmult", &testmatmatmult, NULL));
279566063dSJacob Faibussowitsch   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
28c4762a1bSJed Brown   if (!flg) { /* Create a matrix and test MatSetValues */
29c4762a1bSJed Brown     PetscMPIInt size;
30c4762a1bSJed Brown 
319566063dSJacob Faibussowitsch     PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
329566063dSJacob Faibussowitsch     PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL));
339566063dSJacob Faibussowitsch     PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL));
349566063dSJacob Faibussowitsch     PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, M, N));
359566063dSJacob Faibussowitsch     PetscCall(MatSetType(A, MATAIJ));
369566063dSJacob Faibussowitsch     PetscCall(MatSeqAIJSetPreallocation(A, 9, NULL));
379566063dSJacob Faibussowitsch     PetscCall(MatMPIAIJSetPreallocation(A, 9, NULL, 9, NULL));
389566063dSJacob Faibussowitsch     PetscCall(MatCreate(PETSC_COMM_WORLD, &B));
399566063dSJacob Faibussowitsch     PetscCall(MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, M, N));
409566063dSJacob Faibussowitsch     PetscCall(MatSetType(B, MATHYPRE));
41c4762a1bSJed Brown     if (M == N) {
429566063dSJacob Faibussowitsch       PetscCall(MatHYPRESetPreallocation(B, 9, NULL, 9, NULL));
43c4762a1bSJed Brown     } else {
449566063dSJacob Faibussowitsch       PetscCall(MatHYPRESetPreallocation(B, 6, NULL, 6, NULL));
45c4762a1bSJed Brown     }
46c4762a1bSJed Brown     if (M == N) {
47c4762a1bSJed Brown       for (i = 0; i < M; i++) {
48c4762a1bSJed Brown         PetscInt    cols[] = {0, 1, 2, 3, 4, 5};
49c4762a1bSJed Brown         PetscScalar vals[] = {0, 1. / size, 2. / size, 3. / size, 4. / size, 5. / size};
50c4762a1bSJed Brown         for (j = i - 2; j < i + 1; j++) {
51c4762a1bSJed Brown           if (j >= N) {
529566063dSJacob Faibussowitsch             PetscCall(MatSetValue(A, i, N - 1, (1. * j * N + i) / (3. * N * size), ADD_VALUES));
539566063dSJacob Faibussowitsch             PetscCall(MatSetValue(B, i, N - 1, (1. * j * N + i) / (3. * N * size), ADD_VALUES));
54c4762a1bSJed Brown           } else if (i > j) {
559566063dSJacob Faibussowitsch             PetscCall(MatSetValue(A, i, PetscMin(j, N - 1), (1. * j * N + i) / (2. * N * size), ADD_VALUES));
569566063dSJacob Faibussowitsch             PetscCall(MatSetValue(B, i, PetscMin(j, N - 1), (1. * j * N + i) / (2. * N * size), ADD_VALUES));
57c4762a1bSJed Brown           } else {
589566063dSJacob Faibussowitsch             PetscCall(MatSetValue(A, i, PetscMin(j, N - 1), -1. - (1. * j * N + i) / (4. * N * size), ADD_VALUES));
599566063dSJacob Faibussowitsch             PetscCall(MatSetValue(B, i, PetscMin(j, N - 1), -1. - (1. * j * N + i) / (4. * N * size), ADD_VALUES));
60c4762a1bSJed Brown           }
61c4762a1bSJed Brown         }
629566063dSJacob Faibussowitsch         PetscCall(MatSetValues(A, 1, &i, PetscMin(6, N), cols, vals, ADD_VALUES));
639566063dSJacob Faibussowitsch         PetscCall(MatSetValues(B, 1, &i, PetscMin(6, N), cols, vals, ADD_VALUES));
64c4762a1bSJed Brown       }
65c4762a1bSJed Brown     } else {
66c4762a1bSJed Brown       PetscInt  rows[2];
67c4762a1bSJed Brown       PetscBool test_offproc = PETSC_FALSE;
68c4762a1bSJed Brown 
699566063dSJacob Faibussowitsch       PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_offproc", &test_offproc, NULL));
70c4762a1bSJed Brown       if (test_offproc) {
71c4762a1bSJed Brown         const PetscInt *ranges;
72c4762a1bSJed Brown         PetscMPIInt     rank;
73c4762a1bSJed Brown 
749566063dSJacob Faibussowitsch         PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
759566063dSJacob Faibussowitsch         PetscCall(MatGetOwnershipRanges(A, &ranges));
76c4762a1bSJed Brown         rows[0] = ranges[(rank + 1) % size];
77c4762a1bSJed Brown         rows[1] = ranges[(rank + 1) % size + 1];
78c4762a1bSJed Brown       } else {
799566063dSJacob Faibussowitsch         PetscCall(MatGetOwnershipRange(A, &rows[0], &rows[1]));
80c4762a1bSJed Brown       }
81c4762a1bSJed Brown       for (i = rows[0]; i < rows[1]; i++) {
82c4762a1bSJed Brown         PetscInt    cols[] = {0, 1, 2, 3, 4, 5};
83c4762a1bSJed Brown         PetscScalar vals[] = {-1, 1, -2, 2, -3, 3};
84c4762a1bSJed Brown 
859566063dSJacob Faibussowitsch         PetscCall(MatSetValues(A, 1, &i, PetscMin(6, N), cols, vals, INSERT_VALUES));
869566063dSJacob Faibussowitsch         PetscCall(MatSetValues(B, 1, &i, PetscMin(6, N), cols, vals, INSERT_VALUES));
87c4762a1bSJed Brown       }
88c4762a1bSJed Brown     }
89c4762a1bSJed Brown     /* MAT_FLUSH_ASSEMBLY currently not supported */
909566063dSJacob Faibussowitsch     PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
919566063dSJacob Faibussowitsch     PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
929566063dSJacob Faibussowitsch     PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
939566063dSJacob Faibussowitsch     PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
94c4762a1bSJed Brown 
95c4762a1bSJed Brown #if defined(PETSC_USE_COMPLEX)
96c4762a1bSJed Brown     /* make the matrix imaginary */
979566063dSJacob Faibussowitsch     PetscCall(MatScale(A, PETSC_i));
989566063dSJacob Faibussowitsch     PetscCall(MatScale(B, PETSC_i));
99c4762a1bSJed Brown #endif
100c4762a1bSJed Brown 
101c4762a1bSJed Brown     /* MatAXPY further exercises MatSetValues_HYPRE */
1029566063dSJacob Faibussowitsch     PetscCall(MatAXPY(B, -1., A, DIFFERENT_NONZERO_PATTERN));
1039566063dSJacob Faibussowitsch     PetscCall(MatConvert(B, MATMPIAIJ, MAT_INITIAL_MATRIX, &C));
1049566063dSJacob Faibussowitsch     PetscCall(MatNorm(C, NORM_INFINITY, &err));
10554c59aa7SJacob Faibussowitsch     PetscCheck(err <= PETSC_SMALL, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatSetValues %g", err);
1069566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&B));
1079566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&C));
108c4762a1bSJed Brown   } else {
109c4762a1bSJed Brown     PetscViewer viewer;
110c4762a1bSJed Brown 
1119566063dSJacob Faibussowitsch     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &viewer));
1129566063dSJacob Faibussowitsch     PetscCall(MatSetFromOptions(A));
1139566063dSJacob Faibussowitsch     PetscCall(MatLoad(A, viewer));
1149566063dSJacob Faibussowitsch     PetscCall(PetscViewerDestroy(&viewer));
1159566063dSJacob Faibussowitsch     PetscCall(MatGetSize(A, &M, &N));
116c4762a1bSJed Brown   }
117c4762a1bSJed Brown 
118c4762a1bSJed Brown   /* check conversion routines */
119*f9dc6c1fSStefano Zampini   PetscCall(MatViewFromOptions(A, NULL, "-view_A"));
1209566063dSJacob Faibussowitsch   PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &B));
121*f9dc6c1fSStefano Zampini   PetscCall(MatViewFromOptions(B, NULL, "-view_convert"));
1229566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(B, A, 4, &flg));
123*f9dc6c1fSStefano Zampini   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error Mat HYPRE init");
124*f9dc6c1fSStefano Zampini   PetscCall(MatConvert(A, MATHYPRE, MAT_REUSE_MATRIX, &B));
125*f9dc6c1fSStefano Zampini   PetscCall(MatViewFromOptions(B, NULL, "-view_convert"));
126*f9dc6c1fSStefano Zampini   PetscCall(MatMultEqual(B, A, 4, &flg));
127*f9dc6c1fSStefano Zampini   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error Mat HYPRE reuse");
1289566063dSJacob Faibussowitsch   PetscCall(MatConvert(B, MATIS, MAT_INITIAL_MATRIX, &D));
1299566063dSJacob Faibussowitsch   PetscCall(MatConvert(B, MATIS, MAT_REUSE_MATRIX, &D));
1309566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(D, A, 4, &flg));
13154c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error Mat IS");
1329566063dSJacob Faibussowitsch   PetscCall(MatConvert(B, MATAIJ, MAT_INITIAL_MATRIX, &C));
1339566063dSJacob Faibussowitsch   PetscCall(MatConvert(B, MATAIJ, MAT_REUSE_MATRIX, &C));
1349566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(C, A, 4, &flg));
13554c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error Mat AIJ");
1369566063dSJacob Faibussowitsch   PetscCall(MatAXPY(C, -1., A, SAME_NONZERO_PATTERN));
1379566063dSJacob Faibussowitsch   PetscCall(MatNorm(C, NORM_INFINITY, &err));
13854c59aa7SJacob Faibussowitsch   PetscCheck(err <= PETSC_SMALL, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error Mat AIJ %g", err);
1399566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&C));
1409566063dSJacob Faibussowitsch   PetscCall(MatConvert(D, MATAIJ, MAT_INITIAL_MATRIX, &C));
1419566063dSJacob Faibussowitsch   PetscCall(MatAXPY(C, -1., A, SAME_NONZERO_PATTERN));
1429566063dSJacob Faibussowitsch   PetscCall(MatNorm(C, NORM_INFINITY, &err));
14354c59aa7SJacob Faibussowitsch   PetscCheck(err <= PETSC_SMALL, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error Mat IS %g", err);
1449566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&C));
1459566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&D));
146c4762a1bSJed Brown 
147c4762a1bSJed Brown   /* check MatCreateFromParCSR */
1489566063dSJacob Faibussowitsch   PetscCall(MatHYPREGetParCSR(B, &parcsr));
1499566063dSJacob Faibussowitsch   PetscCall(MatCreateFromParCSR(parcsr, MATAIJ, PETSC_COPY_VALUES, &D));
1509566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&D));
1519566063dSJacob Faibussowitsch   PetscCall(MatCreateFromParCSR(parcsr, MATHYPRE, PETSC_USE_POINTER, &C));
152c4762a1bSJed Brown 
153c4762a1bSJed Brown   /* check MatMult operations */
1549566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(A, B, 4, &flg));
15554c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMult B");
1569566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(A, C, 4, &flg));
15754c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMult C");
1589566063dSJacob Faibussowitsch   PetscCall(MatMultAddEqual(A, B, 4, &flg));
15954c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMultAdd B");
1609566063dSJacob Faibussowitsch   PetscCall(MatMultAddEqual(A, C, 4, &flg));
16154c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMultAdd C");
1629566063dSJacob Faibussowitsch   PetscCall(MatMultTransposeEqual(A, B, 4, &flg));
16354c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMultTranspose B");
1649566063dSJacob Faibussowitsch   PetscCall(MatMultTransposeEqual(A, C, 4, &flg));
16554c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMultTranspose C");
1669566063dSJacob Faibussowitsch   PetscCall(MatMultTransposeAddEqual(A, B, 4, &flg));
16754c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMultTransposeAdd B");
1689566063dSJacob Faibussowitsch   PetscCall(MatMultTransposeAddEqual(A, C, 4, &flg));
16954c59aa7SJacob Faibussowitsch   PetscCheck(flg, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMultTransposeAdd C");
170c4762a1bSJed Brown 
171c4762a1bSJed Brown   /* check PtAP */
172c4762a1bSJed Brown   if (testptap && M == N) {
173c4762a1bSJed Brown     Mat pP, hP;
174c4762a1bSJed Brown 
175c4762a1bSJed Brown     /* PETSc MatPtAP -> output is a MatAIJ
176c4762a1bSJed Brown        It uses HYPRE functions when -matptap_via hypre is specified at command line */
1779566063dSJacob Faibussowitsch     PetscCall(MatPtAP(A, A, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &pP));
1789566063dSJacob Faibussowitsch     PetscCall(MatPtAP(A, A, MAT_REUSE_MATRIX, PETSC_DEFAULT, &pP));
1799566063dSJacob Faibussowitsch     PetscCall(MatNorm(pP, NORM_INFINITY, &norm));
1809566063dSJacob Faibussowitsch     PetscCall(MatPtAPMultEqual(A, A, pP, 10, &flg));
18154c59aa7SJacob Faibussowitsch     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatPtAP_MatAIJ");
182c4762a1bSJed Brown 
183c4762a1bSJed Brown     /* MatPtAP_HYPRE_HYPRE -> output is a MatHYPRE */
1849566063dSJacob Faibussowitsch     PetscCall(MatPtAP(C, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &hP));
1859566063dSJacob Faibussowitsch     PetscCall(MatPtAP(C, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &hP));
1869566063dSJacob Faibussowitsch     PetscCall(MatPtAPMultEqual(C, B, hP, 10, &flg));
18754c59aa7SJacob Faibussowitsch     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatPtAP_HYPRE_HYPRE");
188c20d7725SJed Brown 
189c20d7725SJed Brown     /* Test MatAXPY_Basic() */
1909566063dSJacob Faibussowitsch     PetscCall(MatAXPY(hP, -1., pP, DIFFERENT_NONZERO_PATTERN));
1919566063dSJacob Faibussowitsch     PetscCall(MatHasOperation(hP, MATOP_NORM, &flg));
192c4762a1bSJed Brown     if (!flg) { /* TODO add MatNorm_HYPRE */
1939566063dSJacob Faibussowitsch       PetscCall(MatConvert(hP, MATAIJ, MAT_INPLACE_MATRIX, &hP));
194c4762a1bSJed Brown     }
1959566063dSJacob Faibussowitsch     PetscCall(MatNorm(hP, NORM_INFINITY, &err));
196e00437b9SBarry Smith     PetscCheck(err / norm <= PETSC_SMALL, PetscObjectComm((PetscObject)hP), PETSC_ERR_PLIB, "Error MatPtAP %g %g", err, norm);
1979566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&pP));
1989566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&hP));
199c4762a1bSJed Brown 
200c4762a1bSJed Brown     /* MatPtAP_AIJ_HYPRE -> output can be decided at runtime with -matptap_hypre_outtype */
2019566063dSJacob Faibussowitsch     PetscCall(MatPtAP(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &hP));
2029566063dSJacob Faibussowitsch     PetscCall(MatPtAP(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &hP));
2039566063dSJacob Faibussowitsch     PetscCall(MatPtAPMultEqual(A, B, hP, 10, &flg));
20454c59aa7SJacob Faibussowitsch     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatPtAP_AIJ_HYPRE");
2059566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&hP));
206c4762a1bSJed Brown   }
2079566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&C));
2089566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&B));
209c4762a1bSJed Brown 
210c4762a1bSJed Brown   /* check MatMatMult */
211c4762a1bSJed Brown   if (testmatmatmult) {
2129566063dSJacob Faibussowitsch     PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B));
2139566063dSJacob Faibussowitsch     PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &C));
2149566063dSJacob Faibussowitsch     PetscCall(MatConvert(B, MATHYPRE, MAT_INITIAL_MATRIX, &D));
215c4762a1bSJed Brown 
216c4762a1bSJed Brown     /* PETSc MatMatMult -> output is a MatAIJ
217c4762a1bSJed Brown        It uses HYPRE functions when -matmatmult_via hypre is specified at command line */
2189566063dSJacob Faibussowitsch     PetscCall(MatMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &pAB));
2199566063dSJacob Faibussowitsch     PetscCall(MatMatMult(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &pAB));
2209566063dSJacob Faibussowitsch     PetscCall(MatNorm(pAB, NORM_INFINITY, &norm));
2219566063dSJacob Faibussowitsch     PetscCall(MatMatMultEqual(A, B, pAB, 10, &flg));
22254c59aa7SJacob Faibussowitsch     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatMatMult_AIJ_AIJ");
223c4762a1bSJed Brown 
224c4762a1bSJed Brown     /* MatMatMult_HYPRE_HYPRE -> output is a MatHYPRE */
2259566063dSJacob Faibussowitsch     PetscCall(MatMatMult(C, D, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &CD));
2269566063dSJacob Faibussowitsch     PetscCall(MatMatMult(C, D, MAT_REUSE_MATRIX, PETSC_DEFAULT, &CD));
2279566063dSJacob Faibussowitsch     PetscCall(MatMatMultEqual(C, D, CD, 10, &flg));
22854c59aa7SJacob Faibussowitsch     PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatMatMult_HYPRE_HYPRE");
229c20d7725SJed Brown 
230c20d7725SJed Brown     /* Test MatAXPY_Basic() */
2319566063dSJacob Faibussowitsch     PetscCall(MatAXPY(CD, -1., pAB, DIFFERENT_NONZERO_PATTERN));
232c20d7725SJed Brown 
2339566063dSJacob Faibussowitsch     PetscCall(MatHasOperation(CD, MATOP_NORM, &flg));
234c4762a1bSJed Brown     if (!flg) { /* TODO add MatNorm_HYPRE */
2359566063dSJacob Faibussowitsch       PetscCall(MatConvert(CD, MATAIJ, MAT_INPLACE_MATRIX, &CD));
236c4762a1bSJed Brown     }
2379566063dSJacob Faibussowitsch     PetscCall(MatNorm(CD, NORM_INFINITY, &err));
23854c59aa7SJacob Faibussowitsch     PetscCheck((err / norm) <= PETSC_SMALL, PetscObjectComm((PetscObject)CD), PETSC_ERR_PLIB, "Error MatMatMult %g %g", err, norm);
239c20d7725SJed Brown 
2409566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&C));
2419566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&D));
2429566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&pAB));
2439566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&CD));
244c4762a1bSJed Brown 
245c4762a1bSJed Brown     /* When configured with HYPRE, MatMatMatMult is available for the triplet transpose(aij)-aij-aij */
2469566063dSJacob Faibussowitsch     PetscCall(MatCreateTranspose(A, &C));
2479566063dSJacob Faibussowitsch     PetscCall(MatMatMatMult(C, A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &CAB));
2489566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&C));
2499566063dSJacob Faibussowitsch     PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &C));
2509566063dSJacob Faibussowitsch     PetscCall(MatMatMult(C, A, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &D));
2519566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&C));
2529566063dSJacob Faibussowitsch     PetscCall(MatMatMult(D, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C));
2539566063dSJacob Faibussowitsch     PetscCall(MatNorm(C, NORM_INFINITY, &norm));
2549566063dSJacob Faibussowitsch     PetscCall(MatAXPY(C, -1., CAB, DIFFERENT_NONZERO_PATTERN));
2559566063dSJacob Faibussowitsch     PetscCall(MatNorm(C, NORM_INFINITY, &err));
25654c59aa7SJacob Faibussowitsch     PetscCheck((err / norm) <= PETSC_SMALL, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatMatMatMult %g %g", err, norm);
2579566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&C));
2589566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&D));
2599566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&CAB));
2609566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&B));
261c4762a1bSJed Brown   }
262c4762a1bSJed Brown 
263c4762a1bSJed Brown   /* Check MatView */
2649566063dSJacob Faibussowitsch   PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &B));
2659566063dSJacob Faibussowitsch   PetscCall(MatViewFromOptions(B, NULL, "-view_B"));
266c4762a1bSJed Brown 
267c4762a1bSJed Brown   /* Check MatDuplicate/MatCopy */
268c4762a1bSJed Brown   for (j = 0; j < 3; j++) {
269c4762a1bSJed Brown     MatDuplicateOption dop;
270c4762a1bSJed Brown 
271c4762a1bSJed Brown     dop = MAT_COPY_VALUES;
272c4762a1bSJed Brown     if (j == 1) dop = MAT_DO_NOT_COPY_VALUES;
273c4762a1bSJed Brown     if (j == 2) dop = MAT_SHARE_NONZERO_PATTERN;
274c4762a1bSJed Brown 
275c4762a1bSJed Brown     for (i = 0; i < 3; i++) {
276c4762a1bSJed Brown       MatStructure str;
277c4762a1bSJed Brown 
2789566063dSJacob Faibussowitsch       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Dup/Copy tests: %" PetscInt_FMT " %" PetscInt_FMT "\n", j, i));
279c4762a1bSJed Brown 
280c4762a1bSJed Brown       str = DIFFERENT_NONZERO_PATTERN;
281c4762a1bSJed Brown       if (i == 1) str = SAME_NONZERO_PATTERN;
282c4762a1bSJed Brown       if (i == 2) str = SUBSET_NONZERO_PATTERN;
283c4762a1bSJed Brown 
2849566063dSJacob Faibussowitsch       PetscCall(MatDuplicate(A, dop, &C));
2859566063dSJacob Faibussowitsch       PetscCall(MatDuplicate(B, dop, &D));
286c4762a1bSJed Brown       if (dop != MAT_COPY_VALUES) {
2879566063dSJacob Faibussowitsch         PetscCall(MatCopy(A, C, str));
2889566063dSJacob Faibussowitsch         PetscCall(MatCopy(B, D, str));
289c4762a1bSJed Brown       }
290c4762a1bSJed Brown       /* AXPY with AIJ and HYPRE */
2919566063dSJacob Faibussowitsch       PetscCall(MatAXPY(C, -1.0, D, str));
2929566063dSJacob Faibussowitsch       PetscCall(MatNorm(C, NORM_INFINITY, &err));
293c4762a1bSJed Brown       if (err > PETSC_SMALL) {
2949566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(A, NULL, "-view_duplicate_diff"));
2959566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(B, NULL, "-view_duplicate_diff"));
2969566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(C, NULL, "-view_duplicate_diff"));
297*f9dc6c1fSStefano Zampini         PetscCall(MatViewFromOptions(D, NULL, "-view_duplicate_diff"));
29898921bdaSJacob Faibussowitsch         SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error test 1 MatDuplicate/MatCopy %g (%" PetscInt_FMT ",%" PetscInt_FMT ")", err, j, i);
299c4762a1bSJed Brown       }
300c4762a1bSJed Brown       /* AXPY with HYPRE and HYPRE */
3019566063dSJacob Faibussowitsch       PetscCall(MatAXPY(D, -1.0, B, str));
302c4762a1bSJed Brown       if (err > PETSC_SMALL) {
3039566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(A, NULL, "-view_duplicate_diff"));
3049566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(B, NULL, "-view_duplicate_diff"));
3059566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(D, NULL, "-view_duplicate_diff"));
30698921bdaSJacob Faibussowitsch         SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error test 2 MatDuplicate/MatCopy %g (%" PetscInt_FMT ",%" PetscInt_FMT ")", err, j, i);
307c4762a1bSJed Brown       }
308c4762a1bSJed Brown       /* Copy from HYPRE to AIJ */
3099566063dSJacob Faibussowitsch       PetscCall(MatCopy(B, C, str));
310c4762a1bSJed Brown       /* Copy from AIJ to HYPRE */
3119566063dSJacob Faibussowitsch       PetscCall(MatCopy(A, D, str));
312c4762a1bSJed Brown       /* AXPY with HYPRE and AIJ */
3139566063dSJacob Faibussowitsch       PetscCall(MatAXPY(D, -1.0, C, str));
3149566063dSJacob Faibussowitsch       PetscCall(MatHasOperation(D, MATOP_NORM, &flg));
315c4762a1bSJed Brown       if (!flg) { /* TODO add MatNorm_HYPRE */
3169566063dSJacob Faibussowitsch         PetscCall(MatConvert(D, MATAIJ, MAT_INPLACE_MATRIX, &D));
317c4762a1bSJed Brown       }
3189566063dSJacob Faibussowitsch       PetscCall(MatNorm(D, NORM_INFINITY, &err));
319c4762a1bSJed Brown       if (err > PETSC_SMALL) {
3209566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(A, NULL, "-view_duplicate_diff"));
3219566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(C, NULL, "-view_duplicate_diff"));
3229566063dSJacob Faibussowitsch         PetscCall(MatViewFromOptions(D, NULL, "-view_duplicate_diff"));
32398921bdaSJacob Faibussowitsch         SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error test 3 MatDuplicate/MatCopy %g (%" PetscInt_FMT ",%" PetscInt_FMT ")", err, j, i);
324c4762a1bSJed Brown       }
3259566063dSJacob Faibussowitsch       PetscCall(MatDestroy(&C));
3269566063dSJacob Faibussowitsch       PetscCall(MatDestroy(&D));
327c4762a1bSJed Brown     }
328c4762a1bSJed Brown   }
3299566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&B));
330c4762a1bSJed Brown 
3319566063dSJacob Faibussowitsch   PetscCall(MatHasCongruentLayouts(A, &flg));
332c4762a1bSJed Brown   if (flg) {
333c4762a1bSJed Brown     Vec y, y2;
334c4762a1bSJed Brown 
3359566063dSJacob Faibussowitsch     PetscCall(MatConvert(A, MATHYPRE, MAT_INITIAL_MATRIX, &B));
3369566063dSJacob Faibussowitsch     PetscCall(MatCreateVecs(A, NULL, &y));
3379566063dSJacob Faibussowitsch     PetscCall(MatCreateVecs(B, NULL, &y2));
3389566063dSJacob Faibussowitsch     PetscCall(MatGetDiagonal(A, y));
3399566063dSJacob Faibussowitsch     PetscCall(MatGetDiagonal(B, y2));
3409566063dSJacob Faibussowitsch     PetscCall(VecAXPY(y2, -1.0, y));
3419566063dSJacob Faibussowitsch     PetscCall(VecNorm(y2, NORM_INFINITY, &err));
342c4762a1bSJed Brown     if (err > PETSC_SMALL) {
3439566063dSJacob Faibussowitsch       PetscCall(VecViewFromOptions(y, NULL, "-view_diagonal_diff"));
3449566063dSJacob Faibussowitsch       PetscCall(VecViewFromOptions(y2, NULL, "-view_diagonal_diff"));
34598921bdaSJacob Faibussowitsch       SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Error MatGetDiagonal %g", err);
346c4762a1bSJed Brown     }
3479566063dSJacob Faibussowitsch     PetscCall(MatDestroy(&B));
3489566063dSJacob Faibussowitsch     PetscCall(VecDestroy(&y));
3499566063dSJacob Faibussowitsch     PetscCall(VecDestroy(&y2));
350c4762a1bSJed Brown   }
351c4762a1bSJed Brown 
3529566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&A));
353c4762a1bSJed Brown 
3549566063dSJacob Faibussowitsch   PetscCall(PetscFinalize());
355b122ec5aSJacob Faibussowitsch   return 0;
356c4762a1bSJed Brown }
357c4762a1bSJed Brown 
358c4762a1bSJed Brown /*TEST
359c4762a1bSJed Brown 
360c4762a1bSJed Brown    build:
361c4762a1bSJed Brown       requires: hypre
362c4762a1bSJed Brown 
363c4762a1bSJed Brown    test:
364263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
365c4762a1bSJed Brown       suffix: 1
366c4762a1bSJed Brown       args: -N 11 -M 11
367c4762a1bSJed Brown       output_file: output/ex115_1.out
368c4762a1bSJed Brown 
369c4762a1bSJed Brown    test:
370263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
371c4762a1bSJed Brown       suffix: 2
372c4762a1bSJed Brown       nsize: 3
373c4762a1bSJed Brown       args: -N 13 -M 13 -matmatmult_via hypre
374c4762a1bSJed Brown       output_file: output/ex115_1.out
375c4762a1bSJed Brown 
376c4762a1bSJed Brown    test:
377263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
378c4762a1bSJed Brown       suffix: 3
379c4762a1bSJed Brown       nsize: 4
380c4762a1bSJed Brown       args: -M 13 -N 7 -matmatmult_via hypre
381c4762a1bSJed Brown       output_file: output/ex115_1.out
382c4762a1bSJed Brown 
383c4762a1bSJed Brown    test:
384263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
385c4762a1bSJed Brown       suffix: 4
386c4762a1bSJed Brown       nsize: 2
387c4762a1bSJed Brown       args: -M 12 -N 19
388c4762a1bSJed Brown       output_file: output/ex115_1.out
389c4762a1bSJed Brown 
390c4762a1bSJed Brown    test:
391263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
392c4762a1bSJed Brown       suffix: 5
393c4762a1bSJed Brown       nsize: 3
394c4762a1bSJed Brown       args: -M 13 -N 13 -matptap_via hypre -matptap_hypre_outtype hypre
395c4762a1bSJed Brown       output_file: output/ex115_1.out
396c4762a1bSJed Brown 
397c4762a1bSJed Brown    test:
398263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
399c4762a1bSJed Brown       suffix: 6
400c4762a1bSJed Brown       nsize: 3
401c4762a1bSJed Brown       args: -M 12 -N 19 -test_offproc
402c4762a1bSJed Brown       output_file: output/ex115_1.out
403c4762a1bSJed Brown 
404c4762a1bSJed Brown    test:
405263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
406c4762a1bSJed Brown       suffix: 7
407c4762a1bSJed Brown       nsize: 3
408c4762a1bSJed Brown       args: -M 19 -N 12 -test_offproc -view_B ::ascii_info_detail
409c4762a1bSJed Brown       output_file: output/ex115_7.out
410c4762a1bSJed Brown 
411c4762a1bSJed Brown    test:
412263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
413c4762a1bSJed Brown       suffix: 8
414c4762a1bSJed Brown       nsize: 3
415c4762a1bSJed Brown       args: -M 1 -N 12 -test_offproc
416c4762a1bSJed Brown       output_file: output/ex115_1.out
417c4762a1bSJed Brown 
418c4762a1bSJed Brown    test:
419263f2b91SStefano Zampini       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
420c4762a1bSJed Brown       suffix: 9
421c4762a1bSJed Brown       nsize: 3
422c4762a1bSJed Brown       args: -M 1 -N 2 -test_offproc
423c4762a1bSJed Brown       output_file: output/ex115_1.out
424c4762a1bSJed Brown 
425c4762a1bSJed Brown TEST*/
426