In this work by using a numerical analysis, we investigate in a quantitative way the late-time dynamics of scalar coupled $f(R,mathcal{G})$ gravity. Particularly, we consider a Gauss-Bonnet term coupled to the scalar field coupling function $xi(phi)$, and we study three types of models, one with $f(R)$ terms that are known to provide a viable late-time phenomenology, and two Einstein-Gauss-Bonnet types of models. Our aim is to write the Friedmann equation in terms of appropriate statefinder quantities frequently used in the literature, and we numerically solve it by using physically motivated initial conditions. In the case that $f(R)$ gravity terms are present, the contribution of the Gauss-Bonnet related terms is minor, as we actually expected. This result is robust against changes in the initial conditions of the scalar field, and the reason is the dominating parts of the $f(R)$ gravity sector at late times. In the Einstein-Gauss-Bonnet type of models, we examine two distinct scenarios, firstly by choosing freely the scalar potential and the scalar Gauss-Bonnet coupling $xi(phi)$, in which case the resulting phenomenology is compatible with the latest Planck data and mimics the $Lambda$-Cold-Dark-Matter model. In the second case, since there is no fundamental particle physics reason for the graviton to change its mass, we assume that primordially the tensor perturbations propagate with the speed equal to that of lights, and thus this constraint restricts the functional form of the scalar coupling function $xi(phi)$, which must satisfy the differential equation $ddot{xi}=Hdot{xi}$.