  Dynamical black holes in the non-perturbative regime are not mathematically
well understood. Studying approximate symmetries of spacetimes describing
dynamical black holes gives an insight into their structure. Utilising the
property that approximate symmetries coincide with actual symmetries when they
are present allows one to construct geometric invariants characterising the
symmetry. In this paper, we extend Dain's construction of geometric invariants
characterising stationarity to the case of initial data sets for the Einstein
equations corresponding to black hole spacetimes. We prove the existence and
uniqueness of solutions to a boundary value problem showing that one can always
find approximate Killing vectors in black hole spacetimes and these coincide
with actual Killing vectors when they are present. In the time-symmetric
setting we make use of a 2+1 decomposition to construct a geometric invariant
on a MOTS that vanishes if and only if the Killing initial data equations are
locally satisfied.
