Actual source code: ex110.c

  1: static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **argv)
  6: {
  7:   Mat             A, B;
  8:   PetscInt        i, j, column, *ooj;
  9:   PetscInt       *di, *dj, *oi, *oj, nd;
 10:   const PetscInt *garray;
 11:   PetscScalar    *oa, *da;
 12:   PetscScalar     value;
 13:   PetscRandom     rctx;
 14:   PetscBool       equal, done;
 15:   Mat             AA, AB;
 16:   PetscMPIInt     size, rank;

 19:   PetscInitialize(&argc, &argv, (char *)0, help);
 20:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 22:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);

 24:   /* Create a mpiaij matrix for checking */
 25:   MatCreateAIJ(PETSC_COMM_WORLD, 5, 5, PETSC_DECIDE, PETSC_DECIDE, 0, NULL, 0, NULL, &A);
 26:   MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE);
 27:   MatSetUp(A);
 28:   PetscRandomCreate(PETSC_COMM_WORLD, &rctx);
 29:   PetscRandomSetFromOptions(rctx);

 31:   for (i = 5 * rank; i < 5 * rank + 5; i++) {
 32:     for (j = 0; j < 5 * size; j++) {
 33:       PetscRandomGetValue(rctx, &value);
 34:       column = (PetscInt)(5 * size * PetscRealPart(value));
 35:       PetscRandomGetValue(rctx, &value);
 36:       MatSetValues(A, 1, &i, 1, &column, &value, INSERT_VALUES);
 37:     }
 38:   }
 39:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 40:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);

 42:   MatMPIAIJGetSeqAIJ(A, &AA, &AB, &garray);
 43:   MatGetRowIJ(AA, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&di, (const PetscInt **)&dj, &done);
 44:   MatSeqAIJGetArray(AA, &da);
 45:   MatGetRowIJ(AB, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&oi, (const PetscInt **)&oj, &done);
 46:   MatSeqAIJGetArray(AB, &oa);

 48:   PetscMalloc1(oi[5], &ooj);
 49:   PetscArraycpy(ooj, oj, oi[5]);
 50:   /* modify the column entries in the non-diagonal portion back to global numbering */
 51:   for (i = 0; i < oi[5]; i++) ooj[i] = garray[ooj[i]];

 53:   MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD, 5, 5, PETSC_DETERMINE, PETSC_DETERMINE, di, dj, da, oi, ooj, oa, &B);
 54:   MatSetUp(B);
 55:   MatEqual(A, B, &equal);

 57:   MatRestoreRowIJ(AA, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&di, (const PetscInt **)&dj, &done);
 58:   MatSeqAIJRestoreArray(AA, &da);
 59:   MatRestoreRowIJ(AB, 0, PETSC_FALSE, PETSC_FALSE, &nd, (const PetscInt **)&oi, (const PetscInt **)&oj, &done);
 60:   MatSeqAIJRestoreArray(AB, &oa);


 64:   /* Free spaces */
 65:   PetscFree(ooj);
 66:   PetscRandomDestroy(&rctx);
 67:   MatDestroy(&A);
 68:   MatDestroy(&B);
 69:   PetscFree(oj);
 70:   PetscFinalize();
 71:   return 0;
 72: }

 74: /*TEST

 76:    test:
 77:       nsize: 3

 79: TEST*/