In the context of non-destructive evaluation, the modelling of Ultrasonic Testing (UT) is popular across numerous industrial fields. The construction of efficient numerical schemes for high-frequency transient wave propagation to model UT experiments is an active research topic. In fact, combining accuracy of the simulation with performance can be challenging, in particular when complex configurations are considered. Combining high-order lumped finite elements with explicit time schemes is a popular approach. Such schemes have proven to be efficient when considering quasi-uniform meshes with reasonable variations of the material properties. However, when dealing with non-uniform geometries -- e.g. in the presence of thin layers -- or strong contrasts in the wave speeds the scheme stability condition, a.k.a. CFL, can deteriorate quickly and become very restrictive. One way to efficiently address these issues, when modelling thin layers for example, is to use effective transmission conditions. In this context, we can avoid meshing the intermediate layer. It can be proven that the CFL condition is independent of the effective properties of the layer, including its width. Another way to alleviate the stability condition is to use locally implicit time schemes. In this context, the discriminating elements are discretized using an implicit scheme, while the rest of the domain are discretized with an explicit one. It can be shown through energy arguments that the formulation is stable, and the CFL condition is not constrained by the problematic cells anymore. However, this necessitates solving a linear system at each time step, which can severely affect the performance. In this work, we present another approach to address such issues: the stabilization of leapfrog schemes, that we extend to visco-elastic models. They consist in using, in the problematic subdomains, a polynomial of the stiffness operators so that the maximal stable time step is increased, depending on the order of the polynomial. This approach results in a fully explicit scheme parameterized by the polynomial function. Combined with the mortar element method, this method leads to an improvement on the global CFL condition, as well as on the global performances. We give illustrative results in 2D and 3D in a non-destructive testing context.