We present a method to compute efficiently solutions of systems of ordinary differential equations that possess highly oscillatory forcing terms. This approach is based on asymptotic expansions in inverse powers of the oscillatory parameter,and features two fundamental advantages with respect to standard ODE solvers: rstly, the construction of the numerical solution is more efficient when the system is highly oscillatory, and secondly, the cost of the computation is essentially independent of the oscillatory parameter. Numerical examples are provided, motivated by problems in electronic engineering.