The contact stiffness of the mechanical joint usually becomes the weakest part of the stiffness for the whole machinery equipment, which is one of the important parameters affecting the dynamic characteristics of the engineering machinery. Based on the three-dimensional Weierstrass–Mandelbrot (WM) function, the novel normal contact stiffness model of the joint with the bi-fractal surface is proposed, which comprehensively considers the effects of elastoplastic deformation of asperity and friction factor. The effect of various parameters (fractal dimension, scaling parameter, material parameter, friction factor) on the normal contact stiffness of the joint is analyzed by numerical simulation. The normal contact stiffness of the joint increases with an increase in the fractal dimension, normal load, and material properties and decreases with an increase in the scaling parameter. Meanwhile, the fractal parameters of the equivalent rough surface of the joint are calculated by the structural function method. The experimental results show that when the load is between 14 and 38 N∙m, the error of the model is within 20%. The normal contact stiffness model of the bi-fractal surface joint can provide a theoretical basis for the analysis of the dynamic characteristics of the whole machine at the design stage.