The shape of fractures in an elastic medium under different stress distributions have been well studied in the literature, however, the fracture closure process during unloading has not been investigated thoroughly. The fracture surface is normally assumed to be perfectly smooth so that only the width and stress intensity at the fracture tip are critical in the analysis. In reality, the creation of a fracture in a rock seldom produces smooth surfaces and the resulting asperities on the fracture surfaces can impact fracture closure. Correctly modeling this fracture closure behavior has numerous applications in structural engineering and earth sciences. In this article, we present an approach to model fracture closure behavior in a 2D and 3D elastic media. The fracture surface displacements under arbitrary normal load are derived using superposition method. The contact stress, deformation, and fracture volume evolution can be estimated in a computationally efficient manner for various fracture surface properties. When compared with the traditional integral transform methods used to model fracture closure, the superposition method presented in this study produces comparable results with significant less computation time. In addition, with the aid of parallel computation, large scale fracture closure and contact problems can be successfully simulated using our proposed dynamic fracture closure model (DFCM) with very modest computation times.