I have quickly looked through the paper by Elber et al. (J. Phys. Chem. B 103, 899 (1999)) and the earlier relevant paper by Olender and Elber (J. Chem. Phys. 105, 9299 (1996)).

Two remarks. Firstly, as the authors indicate, the approach is applicable only when 'the reactants and the products are known'. This is a major limitation to a full-blown molecular dynamics (MD) simulation. Secondly, while theoretically the time steps in the method by the authors can be chosen to be arbitrary long, this does not mean that arbitrary long time steps necessarily result in reasonable trajectories. To appreciate this, note that the method proposed by the authors solves the finite-difference version of the underlying MD equations in the norm, with the proviso that unless the time-step is infinitesimally small, even the global minimum of the adopted norm need not be close to the exact trajectory/trajectories corresponding to the given reactants and products (I note in passing that the norm adopted by the authors is the sum of the squares of the local errors. This norm is one of an uncountable set of norms that one can adopt and it is well-known that minimisation of different norms yields in general different solutions.) The situation can be likened with that of solving the simple ordinary differential equation d2f(x)/dx2 + f(x) = 0 over [0,π] subject to the conditions f(0) = 0, f(π) = 0. One has f(x) = sin(x) and f(x) = -sin(x). Now express the differential equation for f(x) in some finite-difference form and demand that a similar error as introduced by Olender and Elber, and Elber et al. for this finite-difference equation (subject to f(0) = 0, 'reactant; and f(π) = 0, 'product') be minimised for the subdivision of the interval [0,π] into an 'arbitrary' number of subintervals. Clearly (and you can readily verify this numerically), unless the lengths of the subintervals are sufficiently small, the calculated solutions that minimise the norm do not resemble sin(x) and -sin(x).   □

The MD timesteps depends on the interaction within the system and the system topology. Generally, it must be smaller than any period of the vibration modes in the system. In order to increase the timestep, which allows a longer simulation, some of the fast freedom degrees must be filtered out or frozen. In practice the length of the timestep must be iteratively decided to obtain the total energy convergence for a NVT ensemble, which indicate the stability of the motion integration. Once the system is equilibrated a longer timestep can be used. 

Thanks for your replies. I would just try to point out my view : There is a difference between an accurate time-integration related to the trajectory space (vibrations) and related to the energy space. The method corrects the deviations in energy space and obtains accurate time-averages, while the authors do not claim that the trajectory space is identical with other spaces with a smaller timestep. They also do not make any assumptions on the kinetics, which underlines their view. In more detail, they look at the path which the system performs in the space of the potential energy. In the case for an arbitrary differential equation, that means that one conserved quantity has to be defined for that case (MD - Energy). If this quantity is defined, its deviation at larger integration steps can be defined, which means that integration at a larger timestep still conserves this quantity. That is my view. 

Dear Emanuel, by minimising the sum of the squares of the errors, the method forces the finite-difference errors, in particular the error corresponding to the finite-difference approximation of the time derivative, to be compensated by the trajectory. This means, as I have indicated earlier, that while theoretically one can choose the time steps to be arbitrarily long, the minimising path/paths corresponding to such long time steps have no relation to the actual path/paths. You could test this point on the simple ordinary differential equation that I have given in my earlier comment. Over and above this aspect, there are issues that I discussed earlier. 

I have now written a short Mathematica® program that highlights some of the points I have raised earlier on this page.