We investigate the possibility of $qq^\prime \bar Q \bar Q^\prime$ tetraquark bound states using $n_f=2+1$ lattice QCD with pion masses $\simeq 164$, $299$ and $415$ MeV. Two types of lattice interpolating operator are chosen, reflecting first diquark-antidiquark and second meson-meson structure. Performing variational analysis using these operators and their mixings, we determine the ground and first excited states from the lattice hadron correlators. Using non-relativistic QCD to simulate the bottom quarks and the Tsukuba formulation of relativistic heavy quarks for charm quarks, we study the $ud\bar b \bar b$, $\ell s\bar b \bar b$ as well as $ud\bar c \bar b$ channels, with $\ell=u,d$. In the case of the $ud\bar b \bar b$ and $\ell s\bar b \bar b$ channels unambiguous signals for $J^P=1^+$ tetraquarks are found with binding energies $189(10)$ and $98(7)$ MeV below the corresponding free two-meson thresholds at the physical point. These tetraquarks are therefore strong-interaction stable, implying they are stable under strong as well as electromagnetic interactions while they can decay weakly. So far these are the first exotic hadrons predicted to have this feature. Further evidence for binding is found in the $ud\bar c \bar b$ channel, whereby the binding energy broadly straddles the electromagnetic stability threshold.
Studying further the quark mass dependence, we vary the heavy quark mass in $ud\bar Q \bar Q$, $\ell s\bar Q \bar Q$ as well as $ud\bar Q \bar b$, $\ell s\bar Q \bar b$ between $\sim0.7$ and 6.3 times the bottom quark mass. The observed mass dependence closely follows a behaviour argued from phenomenological considerations of the heavy quark potential.