Interpolating precise ephemerides of GPS satellites

My objective is to obtain at least a sub-meter position accuracy by interpolation of the GPS ephemerides.

There is conflicting research outside, where some authors state that a simple Lagrange or least-squares best fit polynomial interpolation between 15-minute GPS ephemeris intervals is sufficient to provide orbit position accuracy (geocentric Earth-fixed frame) in the decimeter to millimeter level (M. Horemuz, "Polynomial interpolation of GPS satellite coordinates", GPS Solutions, Feb 2006).

On the other hand, other sources (for example, the Bernese GNSS v5.2 manual, the ESA Navigation Guide Book etc) state that it is necessary to involve the dynamics too, and integrate the equations of motion of the GPS satellites in order to achieve the desired accuracy, at the expense of computational effort, which other researchers claim it is too much of an 'overkill'.

Basically, one side shows positive results for a mathematically straightforward interpolation procedure, whereas the other side vouches for a complex interpolation that involves the integration (in the calculus sense) of vehicle dynamics.

I would like to get some professional advice on people who have interpolated the orbits of GPS satellites before, on which method is actually necessary to obtain at least a sub-meter accuracy?

Thanks! -Sam