<p>Liquid bridge force plays an important role in geology and behavior of wet granular materials. At present, many discrete-element-method (DEM) numerical models involving liquid bridges have been developed, but few of them were capable of characterizing liquid bridges comprehensively. In addition, these models usually neglected surface roughness of spheres, leading to overestimation of liquid bridge force. In this study, an analytical expression of liquid bridge force between two unequal-sized spheres in contact when the liquid bridge volume tends to zero was derived. It showed that a large liquid bridge force could be induced by an extremely small liquid volume, resulting in overestimation of liquid bridge force under nearly dry condition. To solve this problem, the scheme of an improved liquid bridge model (ILBM) incorporating roughness coefficient was proposed first. Subsequently, a numerical method incorporating ILBM into DEM was developed. The developed method can calculate liquid bridge forces accurately and encompasses all the parameters of liquid bridges. Furthermore, it avoids unreasonable liquid bridge force under extremely small liquid bridge volume and considers the formation and rupture of liquid bridges. The accuracy and reliability of the developed method were validated by axial compression tests. Finally, the evolution of number of liquid bridges, coordination number, and liquid bridge force during the axial compression tests were analyzed and the discrepancy between numerical and experimental results was discussed.</p> Graphical Abstract <p></p>