(MG) models [@Koyama:2015vza; @Ishak:2018his].

Modifying the Einstein-Hilbert action, however, introduces in principle an additional degree of freedom that is conformally coupled to matter and can produce significant deviations from the predictions of GR, which have passed a wide array of precise observational tests, especially in the Solar System [@Will:2005va]. Furthermore, the recent simultaneous detection of gravitational waves and EM counterparts by the LIGO/Virgo collaboration [@TheLIGOScientific:2017qsa; @Goldstein:2017mmi; @Savchenko:2017ffs; @Monitor:2017mdv; @GBM:2017lvd], has placed additional constraints [@Sakstein:2017xjx; @Ezquiaga:2017ekz; @Creminelli:2017sry; @Baker:2017hug] into the form of the most general expression of a scalar-tensor theory that produces second order equations of motion, described by the Horndeski action [@Horndeski1974; @PhysRevD.84.064039].

In order to be able to confront such tight constraints successfully, while at the same time provide a stable self-accelerative cosmic mechanism, viable candidates contain a restoring property, called “screening" [@Khoury:2010xi; @Khoury:2013tda], which is a dynamical mechanism that weakens the additional fifth forces in high-density environments through the corresponding scalar field self-interactions. In the Vainshtein mechanism [@VAINSHTEIN1972393; @Babichev:2013usa], GR is recovered thanks to the second derivative terms in the scalar field Lagrangian, that become large in high density environments and effectively weaken the coupling to the matter sources. The Vainshtein mechanism is very efficient in the vicinity 