We analyze the efficiency of available algorithms for the simulation of classical fidelity and show that their computational costs increase exponentially with the number of degrees of freedom. Then we present an algorithm for which the number of trajectories needed for convergence is independent of the system's dimensionality and show that, within a continuous family of algorithms, this algorithm is the only one with this property. Simultaneously we propose a general analytical approach to estimate efficiency of trajectory-based methods and suggest how to use it to accelerate calculations of other classical correlation functions. Converged numerical results are provided for systems with phase space volume 2(1700) times larger than the volume of the initial state.