The integral model of dry friction components is built with assumption of classical Coulomb friction law and with specially developed model of normal stress distribution coupled with rolling resistance for elliptic contact shape. In order to avoid a necessity of numerical integration over the contact area at each the numerical simulation step, few versions of approximate model are developed and then tested numerically. In the numerical experiments the simulation results of the Celtic stone with the friction forces modelled by the use of approximants of different complexity (from no coupling between friction force and torque to the second order Padé approximation) are compared to results obtained from model with friction approximated in the form of piecewise polynomial functions (based on the Taylor series with hertzian stress distribution). The coefficients of the corresponding approximate models are found by the use of optimization methods, like as in identification process using the real experiment data.