Actual source code: ex35.c

  1: static char help[] = "Exhaustive memory tracking for DMPlex.\n\n\n";

  3: #include <petscdmplex.h>

  5: static PetscErrorCode EstimateMemory(DM dm, PetscLogDouble *est)
  6: {
  7:   DMLabel  marker;
  8:   PetscInt cdim, depth, d, pStart, pEnd, p, Nd[4] = {0, 0, 0, 0}, lsize = 0, rmem = 0, imem = 0;
  9:   PetscInt coneSecMem = 0, coneMem = 0, supportSecMem = 0, supportMem = 0, labelMem = 0;

 12:   PetscPrintf(PETSC_COMM_SELF, "Memory Estimates\n");
 13:   DMGetCoordinateDim(dm, &cdim);
 14:   DMPlexGetDepth(dm, &depth);
 15:   DMPlexGetChart(dm, &pStart, &pEnd);
 16:   for (d = 0; d <= depth; ++d) {
 17:     PetscInt start, end;

 19:     DMPlexGetDepthStratum(dm, d, &start, &end);
 20:     Nd[d] = end - start;
 21:   }
 22:   /* Coordinates: 3 Nv reals + 2*Nv + 2*Nv ints */
 23:   rmem += cdim * Nd[0];
 24:   imem += 2 * Nd[0] + 2 * Nd[0];
 25:   PetscPrintf(PETSC_COMM_SELF, "  Coordinate mem:  %" PetscInt_FMT " %" PetscInt_FMT "\n", (PetscInt)(cdim * Nd[0] * sizeof(PetscReal)), (PetscInt)(4 * Nd[0] * sizeof(PetscInt)));
 26:   /* Depth:       Nc+Nf+Ne+Nv ints */
 27:   for (d = 0; d <= depth; ++d) labelMem += Nd[d];
 28:   /* Cell Type:   Nc+Nf+Ne+Nv ints */
 29:   for (d = 0; d <= depth; ++d) labelMem += Nd[d];
 30:   /* Marker */
 31:   DMGetLabel(dm, "marker", &marker);
 32:   if (marker) DMLabelGetStratumSize(marker, 1, &lsize);
 33:   labelMem += lsize;
 34:   PetscPrintf(PETSC_COMM_SELF, "  Label mem:       %" PetscInt_FMT "\n", (PetscInt)(labelMem * sizeof(PetscInt)));
 35:   //imem += labelMem;
 36:   /* Cones and Orientations:       4 Nc + 3 Nf + 2 Ne ints + (Nc+Nf+Ne) ints no separate orientation section */
 37:   for (d = 0; d <= depth; ++d) coneSecMem += 2 * Nd[d];
 38:   for (p = pStart; p < pEnd; ++p) {
 39:     PetscInt csize;

 41:     DMPlexGetConeSize(dm, p, &csize);
 42:     coneMem += csize;
 43:   }
 44:   PetscPrintf(PETSC_COMM_SELF, "  Cone mem:        %" PetscInt_FMT " %" PetscInt_FMT " (%" PetscInt_FMT ")\n", (PetscInt)(coneMem * sizeof(PetscInt)), (PetscInt)(coneSecMem * sizeof(PetscInt)), (PetscInt)(coneMem * sizeof(PetscInt)));
 45:   imem += 2 * coneMem + coneSecMem;
 46:   /* Supports:       4 Nc + 3 Nf + 2 Ne ints + Nc+Nf+Ne ints */
 47:   for (d = 0; d <= depth; ++d) supportSecMem += 2 * Nd[d];
 48:   for (p = pStart; p < pEnd; ++p) {
 49:     PetscInt ssize;

 51:     DMPlexGetSupportSize(dm, p, &ssize);
 52:     supportMem += ssize;
 53:   }
 54:   PetscPrintf(PETSC_COMM_SELF, "  Support mem:     %" PetscInt_FMT " %" PetscInt_FMT "\n", (PetscInt)(supportMem * sizeof(PetscInt)), (PetscInt)(supportSecMem * sizeof(PetscInt)));
 55:   imem += supportMem + supportSecMem;
 56:   *est = ((PetscLogDouble)imem) * sizeof(PetscInt) + ((PetscLogDouble)rmem) * sizeof(PetscReal);
 57:   PetscPrintf(PETSC_COMM_WORLD, "  Estimated memory %" PetscInt_FMT "\n", (PetscInt)*est);
 58:   return 0;
 59: }

 61: int main(int argc, char **argv)
 62: {
 63:   DM             dm;
 64:   PetscBool      trace = PETSC_FALSE, checkMemory = PETSC_TRUE, auxMemory = PETSC_FALSE;
 65:   PetscLogDouble before, after, est                                       = 0, clean, max;

 68:   PetscInitialize(&argc, &argv, NULL, help);
 69:   PetscOptionsGetBool(NULL, NULL, "-trace", &trace, NULL);
 70:   PetscOptionsGetBool(NULL, NULL, "-check_memory", &checkMemory, NULL);
 71:   PetscOptionsGetBool(NULL, NULL, "-aux_memory", &auxMemory, NULL);
 72:   PetscMemorySetGetMaximumUsage();
 73:   PetscMallocGetCurrentUsage(&before);
 74:   if (trace) PetscMallocTraceSet(NULL, PETSC_TRUE, 5000.);
 75:   DMCreate(PETSC_COMM_WORLD, &dm);
 76:   DMSetType(dm, DMPLEX);
 77:   DMSetFromOptions(dm);
 78:   if (trace) PetscMallocTraceSet(NULL, PETSC_FALSE, 5000);
 79:   DMViewFromOptions(dm, NULL, "-dm_view");
 80:   PetscMallocGetCurrentUsage(&after);
 81:   PetscMemoryGetMaximumUsage(&max);
 82:   EstimateMemory(dm, &est);
 83:   DMDestroy(&dm);
 84:   PetscMallocGetCurrentUsage(&clean);
 85:   PetscPrintf(PETSC_COMM_WORLD, "Measured Memory\n");
 86:   if (auxMemory) {
 87:     PetscPrintf(PETSC_COMM_WORLD, "  Initial memory         %" PetscInt_FMT "\n  Extra memory for build %" PetscInt_FMT "\n  Memory after destroy   %" PetscInt_FMT "\n", (PetscInt)before, (PetscInt)(max - after), (PetscInt)clean);
 88:   }
 89:   if (checkMemory) {
 90:     PetscPrintf(PETSC_COMM_WORLD, "  Memory for mesh  %" PetscInt_FMT "\n", (PetscInt)(after - before));
 91:     PetscPrintf(PETSC_COMM_WORLD, "Discrepancy %" PetscInt_FMT "\n", (PetscInt)PetscAbsReal(after - before - est));
 92:   }
 93:   PetscFinalize();
 94:   return 0;
 95: }

 97: /*TEST
 98:   build:
 99:     requires: !defined(PETSC_USE_64BIT_INDICES) double !complex !defined(PETSCTEST_VALGRIND)

101:   # Memory checks cannot be included in tests because the allocated memory differs among environments
102:   testset:
103:     args: -malloc_requested_size -dm_plex_box_faces 5,5 -check_memory 0
104:     test:
105:       suffix: tri
106:       requires: triangle
107:       args: -dm_plex_simplex 1 -dm_plex_interpolate 0

109:     test:
110:       suffix: tri_interp
111:       requires: triangle
112:       args: -dm_plex_simplex 1 -dm_plex_interpolate 1

114:     test:
115:       suffix: quad
116:       args: -dm_plex_simplex 0 -dm_plex_interpolate 0

118:     test:
119:       suffix: quad_interp
120:       args: -dm_plex_simplex 0 -dm_plex_interpolate 1

122:   # Memory checks cannot be included in tests because the allocated memory differs among environments
123:   testset:
124:     args: -malloc_requested_size -dm_plex_dim 3 -dm_plex_box_faces 5,5,5 -check_memory 0

126:     # Filter out label memory because tet mesher produce different surface meshes for different compilers
127:     test:
128:       suffix: tet
129:       requires: ctetgen
130:       filter: grep -v "Label mem:"
131:       args: -dm_plex_simplex 1 -dm_plex_interpolate 0

133:     # Filter out label memory because tet mesher produce different surface meshes for different compilers
134:     test:
135:       suffix: tet_interp
136:       requires: ctetgen
137:       filter: grep -v "Label mem:"
138:       args: -dm_plex_simplex 1 -dm_plex_interpolate 1

140:     test:
141:       suffix: hex
142:       args: -dm_plex_simplex 0 -dm_plex_interpolate 0

144:     test:
145:       suffix: hex_interp
146:       args: -dm_plex_simplex 0 -dm_plex_interpolate 1
147: TEST*/