1 static char help[] = "Exhaustive memory tracking for DMPlex.\n\n\n"; 2 3 #include <petscdmplex.h> 4 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; 10 PetscErrorCode ierr; 11 12 PetscFunctionBeginUser; 13 ierr = PetscPrintf(PETSC_COMM_SELF, "Memory Estimates\n");CHKERRQ(ierr); 14 ierr = DMGetCoordinateDim(dm, &cdim);CHKERRQ(ierr); 15 ierr = DMPlexGetDepth(dm, &depth);CHKERRQ(ierr); 16 ierr = DMPlexGetChart(dm, &pStart, &pEnd);CHKERRQ(ierr); 17 for (d = 0; d <= depth; ++d) { 18 PetscInt start, end; 19 20 ierr = DMPlexGetDepthStratum(dm, d, &start, &end);CHKERRQ(ierr); 21 Nd[d] = end - start; 22 } 23 /* Coordinates: 3 Nv reals + 2*Nv + 2*Nv ints */ 24 rmem += cdim*Nd[0]; 25 imem += 2*Nd[0] + 2*Nd[0]; 26 ierr = PetscPrintf(PETSC_COMM_SELF, " Coordinate mem: %D %D\n", cdim*Nd[0]*sizeof(PetscReal), 4*Nd[0]*sizeof(PetscInt));CHKERRQ(ierr); 27 /* Depth: Nc+Nf+Ne+Nv ints */ 28 for (d = 0; d <= depth; ++d) labelMem += Nd[d]; 29 /* Cell Type: Nc+Nf+Ne+Nv ints */ 30 for (d = 0; d <= depth; ++d) labelMem += Nd[d]; 31 /* Marker */ 32 ierr = DMGetLabel(dm, "marker", &marker);CHKERRQ(ierr); 33 if (marker) {ierr = DMLabelGetStratumSize(marker, 1, &lsize);CHKERRQ(ierr);} 34 labelMem += lsize; 35 ierr = PetscPrintf(PETSC_COMM_SELF, " Label mem: %D\n", labelMem*sizeof(PetscInt));CHKERRQ(ierr); 36 //imem += labelMem; 37 /* Cones and Orientations: 4 Nc + 3 Nf + 2 Ne ints + (Nc+Nf+Ne) ints no separate orientation section */ 38 for (d = 0; d <= depth; ++d) coneSecMem += 2*Nd[d]; 39 for (p = pStart; p < pEnd; ++p) { 40 PetscInt csize; 41 42 ierr = DMPlexGetConeSize(dm, p, &csize);CHKERRQ(ierr); 43 coneMem += csize; 44 } 45 ierr = PetscPrintf(PETSC_COMM_SELF, " Cone mem: %D %D (%D)\n", coneMem*sizeof(PetscInt), coneSecMem*sizeof(PetscInt), coneMem*sizeof(PetscInt));CHKERRQ(ierr); 46 imem += 2*coneMem + coneSecMem; 47 /* Supports: 4 Nc + 3 Nf + 2 Ne ints + Nc+Nf+Ne ints */ 48 for (d = 0; d <= depth; ++d) supportSecMem += 2*Nd[d]; 49 for (p = pStart; p < pEnd; ++p) { 50 PetscInt ssize; 51 52 ierr = DMPlexGetSupportSize(dm, p, &ssize);CHKERRQ(ierr); 53 supportMem += ssize; 54 } 55 ierr = PetscPrintf(PETSC_COMM_SELF, " Support mem: %D %D\n", supportMem*sizeof(PetscInt), supportSecMem*sizeof(PetscInt));CHKERRQ(ierr); 56 imem += supportMem + supportSecMem; 57 *est = ((PetscLogDouble) imem)*sizeof(PetscInt) + ((PetscLogDouble) rmem)*sizeof(PetscReal); 58 ierr = PetscPrintf(PETSC_COMM_WORLD, " Estimated memory %D\n", (PetscInt) *est);CHKERRQ(ierr); 59 PetscFunctionReturn(0); 60 } 61 62 int main(int argc, char **argv) 63 { 64 DM dm; 65 PetscBool trace = PETSC_FALSE, checkMemory = PETSC_TRUE, auxMemory = PETSC_FALSE; 66 PetscLogDouble before, after, est = 0, clean, max; 67 PetscErrorCode ierr; 68 69 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 70 ierr = PetscOptionsGetBool(NULL, NULL, "-trace", &trace, NULL);CHKERRQ(ierr); 71 ierr = PetscOptionsGetBool(NULL, NULL, "-check_memory", &checkMemory, NULL);CHKERRQ(ierr); 72 ierr = PetscOptionsGetBool(NULL, NULL, "-aux_memory", &auxMemory, NULL);CHKERRQ(ierr); 73 ierr = PetscMemorySetGetMaximumUsage();CHKERRQ(ierr); 74 ierr = PetscMallocGetCurrentUsage(&before);CHKERRQ(ierr); 75 if (trace) {ierr = PetscMallocTraceSet(NULL, PETSC_TRUE, 5000.);CHKERRQ(ierr);} 76 ierr = DMCreate(PETSC_COMM_WORLD, &dm);CHKERRQ(ierr); 77 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 78 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 79 if (trace) {ierr = PetscMallocTraceSet(NULL, PETSC_FALSE, 5000);CHKERRQ(ierr);} 80 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 81 ierr = PetscMallocGetCurrentUsage(&after);CHKERRQ(ierr); 82 ierr = PetscMemoryGetMaximumUsage(&max);CHKERRQ(ierr); 83 ierr = EstimateMemory(dm, &est);CHKERRQ(ierr); 84 ierr = DMDestroy(&dm);CHKERRQ(ierr); 85 ierr = PetscMallocGetCurrentUsage(&clean);CHKERRQ(ierr); 86 ierr = PetscPrintf(PETSC_COMM_WORLD, "Measured Memory\n");CHKERRQ(ierr); 87 if (auxMemory) { 88 ierr = PetscPrintf(PETSC_COMM_WORLD, " Initial memory %D\n Extra memory for build %D\n Memory after destroy %D\n", 89 (PetscInt) before, (PetscInt) (max-after), (PetscInt) clean);CHKERRQ(ierr); 90 } 91 if (checkMemory) { 92 ierr = PetscPrintf(PETSC_COMM_WORLD, " Memory for mesh %D\n", (PetscInt) (after-before));CHKERRQ(ierr); 93 ierr = PetscPrintf(PETSC_COMM_WORLD, "Discrepancy %D\n", (PetscInt) PetscAbsReal(after-before-est));CHKERRQ(ierr); 94 } 95 ierr = PetscFinalize(); 96 return ierr; 97 } 98 99 /*TEST 100 build: 101 requires: !defined(PETSC_USE_64BIT_INDICES) double !complex !defined(PETSCTEST_VALGRIND) 102 103 # Memory checks cannot be included in tests because the allocated memory differs among environments 104 testset: 105 args: -malloc_requested_size -dm_plex_box_faces 5,5 -check_memory 0 106 test: 107 suffix: tri 108 requires: triangle 109 args: -dm_plex_simplex 1 -dm_plex_interpolate 0 110 111 test: 112 suffix: tri_interp 113 requires: triangle 114 args: -dm_plex_simplex 1 -dm_plex_interpolate 1 115 116 test: 117 suffix: quad 118 args: -dm_plex_simplex 0 -dm_plex_interpolate 0 119 120 test: 121 suffix: quad_interp 122 args: -dm_plex_simplex 0 -dm_plex_interpolate 1 123 124 # Memory checks cannot be included in tests because the allocated memory differs among environments 125 testset: 126 args: -malloc_requested_size -dm_plex_dim 3 -dm_plex_box_faces 5,5,5 -check_memory 0 127 test: 128 suffix: tet 129 requires: ctetgen 130 args: -dm_plex_simplex 1 -dm_plex_interpolate 0 131 132 test: 133 suffix: tet_interp 134 requires: ctetgen 135 args: -dm_plex_simplex 1 -dm_plex_interpolate 1 136 137 test: 138 suffix: hex 139 args: -dm_plex_simplex 0 -dm_plex_interpolate 0 140 141 test: 142 suffix: hex_interp 143 args: -dm_plex_simplex 0 -dm_plex_interpolate 1 144 TEST*/ 145