  In this thesis, we use numerical relativity to investigate gravitational
waves from binary black holes in extensions of GR. We first study spherically
symmetric gravitational collapse in cubic Horndeski theories of gravity. By
varying the coupling constants and the initial amplitude of the scalar field,
we determine the region in the space of couplings and amplitudes for which it
is possible to construct global solutions to the Horndeski theories.
Furthermore, we identify the regime of validity of effective field theory (EFT)
as the sub-region for which a certain weak coupling condition remains small at
all times. We study black hole binary mergers in these cubic Horndeski theories
of gravity, treating them fully non-linearly. In the regime of validity of EFT,
the mismatch of the gravitational wave strain between Horndeski and GR (coupled
to a scalar field) can be larger than $30\%$ in the Advanced LIGO mass range.
Initial data and coupling constants are chosen so the theory always remains in
the weakly coupled regime. We observe that the waveform in Horndeski theories
is shifted by an amount much larger than the smallness parameter that controls
initial data. This effect is generic and may be present in other theories of
gravity involving higher derivatives. We explore a higher-order curvature
correction of GR. Guided by toy models, we develop systems capable of
reproducing the low energy behaviour of many such theories with a fully
nonlinear/non-perturbative approach. We evolve binary black holes, observing a
shift in phase accumulated over time which is not statistically significant
when compared to GR, for the methods and coupling used. Finally, we present
AHFinder, a flexible multi-purpose tool to find apparent horizons in the
open-source numerical relativity code GRChombo.
