Actual source code: ex169.c


  2: static char help[] = "Test memory leak when duplicating a redundant matrix.\n\n";

  4: /*
  5:   Include "petscmat.h" so that we can use matrices.
  6:   automatically includes:
  7:      petscsys.h    - base PETSc routines   petscvec.h    - vectors
  8:      petscmat.h    - matrices
  9:      petscis.h     - index sets            petscviewer.h - viewers
 10: */
 11: #include <petscmat.h>

 13: int main(int argc, char **args)
 14: {
 15:   Mat          A, Ar, C;
 16:   PetscViewer  fd;                       /* viewer */
 17:   char         file[PETSC_MAX_PATH_LEN]; /* input file name */
 18:   PetscInt     ns = 2;
 19:   PetscMPIInt  size;
 20:   PetscSubcomm subc;
 21:   PetscBool    flg;

 23:   PetscFunctionBeginUser;
 24:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 25:   /*
 26:      Determine files from which we read the two linear systems
 27:      (matrix and right-hand-side vector).
 28:   */
 29:   PetscCall(PetscOptionsGetString(NULL, NULL, "-f0", file, sizeof(file), &flg));
 30:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER, "Must indicate binary file with the -f0 option");
 31:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd));
 32:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 33:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Reading matrix with %d processors\n", size));
 34:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 35:   PetscCall(MatLoad(A, fd));
 36:   PetscCall(PetscViewerDestroy(&fd));
 37:   /*
 38:      Determines amount of subcomunicators
 39:   */
 40:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nsub", &ns, NULL));
 41:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Splitting in %" PetscInt_FMT " subcommunicators\n", ns));
 42:   PetscCall(PetscSubcommCreate(PetscObjectComm((PetscObject)A), &subc));
 43:   PetscCall(PetscSubcommSetNumber(subc, ns));
 44:   PetscCall(PetscSubcommSetType(subc, PETSC_SUBCOMM_CONTIGUOUS));
 45:   PetscCall(PetscSubcommSetFromOptions(subc));
 46:   PetscCall(MatCreateRedundantMatrix(A, 0, PetscSubcommChild(subc), MAT_INITIAL_MATRIX, &Ar));
 47:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Copying matrix\n"));
 48:   PetscCall(MatDuplicate(Ar, MAT_COPY_VALUES, &C));
 49:   PetscCall(MatAXPY(Ar, 0.1, C, DIFFERENT_NONZERO_PATTERN));
 50:   PetscCall(PetscSubcommDestroy(&subc));

 52:   /*
 53:      Free memory
 54:   */
 55:   PetscCall(MatDestroy(&A));
 56:   PetscCall(MatDestroy(&Ar));
 57:   PetscCall(MatDestroy(&C));
 58:   PetscCall(PetscFinalize());
 59:   return 0;
 60: }

 62: /*TEST

 64:    test:
 65:       nsize: 4
 66:       requires: !complex double !defined(PETSC_USE_64BIT_INDICES)
 67:       args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/ns-real-int32-float64 -malloc_dump

 69: TEST*/