Many studies have used fractal theory to characterize pore structure distribution heterogeneity through mercury intake curves. However, there is relatively little research on the fractal model calculation of mercury removal curves. In this study, a high-pressure mercury intrusion test is used to describe the pore and fracture distribution heterogeneity (PFDH). The fractal physical meaning of the mercury removal curve was determined by calculating the change in the curve’s fractal dimension value. The results are as follows. (1) According to mercury removal efficiency and porosity, samples can be divided into types A (mercury removal efficiency above 35%) and B (mercury removal efficiency below 35%). In general, type A sample belongs to micro-pore-developed types, and type B samples belong to the macro-pore-developed type. (2) The Menger model (M) represents the complexity of a specific surface area, while the Sierpinski model (S) represents the roughness of the pore volume. Among all the samples, the lower-pore-volume region controls PFDH. (3) According to the calculation results of the single fractal model, it can be seen that the PFDH of type B is stronger than that of type A, which is similar to the results of mercury intrusion. According to the calculation structure of the multifractal model, it can be seen that the volume distribution heterogeneity of type B under various pores is significantly stronger than that of type A. This is opposite to the result of mercury injection. (4) DM has a relationship with the pore volume percentage at different stages, so the M model at the mercury inlet stage can better characterize PFDH at the mercury inlet stage.