order quantities, equal to a quadratic “source” term formed with the first order perturbations [@GlNiPrPucqg]. The source term is complicated and is delicate to handle numerically when evolving the perturbations. A place where this was explicitly done was for instance in the evolution of boosted black hole collisions to second order [@GlNiPrPuboost2]. The calculations are lengthy and complicated. Fortunately, in our case one can proceed in a different way. We have just shown that there is a gauge in which the first order perturbations (the $\ell=1$ ones) vanish. Therefore in that gauge one can write a second order Zerilli equation that is source-free. Moreover, the linear portion of the second order perturbative equation is exactly the same as the first order perturbative equation (the Zerilli equation), and this equation can be written in terms of quantities that are gauge invariant. We notice, however, that eliminating the first order $\ell=1$ terms through a first order gauge transformation introduces second-order changes in the metric that are not a pure second order gauge transformation, and must be taken into account. One must then be careful in making sure that the initial data for the $\ell=2$ perturbations used in the Zerilli equation corresponds precisely 