xref: /petsc/src/mat/tests/ex110.c (revision 58d68138c660dfb4e9f5b03334792cd4f2ffd7cc)
1 static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n";
2 
3 #include <petscmat.h>
4 
5 int main(int argc, char **argv) {
6   Mat             A, B;
7   PetscInt        i, j, column, *ooj;
8   PetscInt       *di, *dj, *oi, *oj, nd;
9   const PetscInt *garray;
10   PetscScalar    *oa, *da;
11   PetscScalar     value;
12   PetscRandom     rctx;
13   PetscBool       equal, done;
14   Mat             AA, AB;
15   PetscMPIInt     size, rank;
16 
17   PetscFunctionBeginUser;
18   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
19   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
20   PetscCheck(size > 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run with 2 or more processes");
21   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
22 
23   /* Create a mpiaij matrix for checking */
24   PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, 5, 5, PETSC_DECIDE, PETSC_DECIDE, 0, NULL, 0, NULL, &A));
25   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE));
26   PetscCall(MatSetUp(A));
27   PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
28   PetscCall(PetscRandomSetFromOptions(rctx));
29 
30   for (i = 5 * rank; i < 5 * rank + 5; i++) {
31     for (j = 0; j < 5 * size; j++) {
32       PetscCall(PetscRandomGetValue(rctx, &value));
33       column = (PetscInt)(5 * size * PetscRealPart(value));
34       PetscCall(PetscRandomGetValue(rctx, &value));
35       PetscCall(MatSetValues(A, 1, &i, 1, &column, &value, INSERT_VALUES));
36     }
37   }
38   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
39   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
40 
41   PetscCall(MatMPIAIJGetSeqAIJ(A, &AA, &AB, &garray));
42   PetscCall(MatGetRowIJ(AA, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&di, (const PetscInt **)&dj, &done));
43   PetscCall(MatSeqAIJGetArray(AA, &da));
44   PetscCall(MatGetRowIJ(AB, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&oi, (const PetscInt **)&oj, &done));
45   PetscCall(MatSeqAIJGetArray(AB, &oa));
46 
47   PetscCall(PetscMalloc1(oi[5], &ooj));
48   PetscCall(PetscArraycpy(ooj, oj, oi[5]));
49   /* modify the column entries in the non-diagonal portion back to global numbering */
50   for (i = 0; i < oi[5]; i++) { ooj[i] = garray[ooj[i]]; }
51 
52   PetscCall(MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD, 5, 5, PETSC_DETERMINE, PETSC_DETERMINE, di, dj, da, oi, ooj, oa, &B));
53   PetscCall(MatSetUp(B));
54   PetscCall(MatEqual(A, B, &equal));
55 
56   PetscCall(MatRestoreRowIJ(AA, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&di, (const PetscInt **)&dj, &done));
57   PetscCall(MatSeqAIJRestoreArray(AA, &da));
58   PetscCall(MatRestoreRowIJ(AB, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&oi, (const PetscInt **)&oj, &done));
59   PetscCall(MatSeqAIJRestoreArray(AB, &oa));
60 
61   PetscCheck(equal, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Likely a bug in MatCreateMPIAIJWithSplitArrays()");
62 
63   /* Free spaces */
64   PetscCall(PetscFree(ooj));
65   PetscCall(PetscRandomDestroy(&rctx));
66   PetscCall(MatDestroy(&A));
67   PetscCall(MatDestroy(&B));
68   PetscCall(PetscFree(oj));
69   PetscCall(PetscFinalize());
70   return 0;
71 }
72 
73 /*TEST
74 
75    test:
76       nsize: 3
77 
78 TEST*/
79