  A search for supersymmetry involving the pair production of gluinos decaying
via off-shell third-generation squarks into the lightest neutralino
($\tilde\chi^0_1$) is reported. It exploits LHC proton$-$proton collision data
at a centre-of-mass energy $\sqrt{s} = 13$ TeV with an integrated luminosity of
139 fb$^{-1}$ collected with the ATLAS detector from 2015 to 2018. The search
uses events containing large missing transverse momentum, up to one electron or
muon, and several energetic jets, at least three of which must be identified as
containing $b$-hadrons. Both a simple kinematic event selection and an event
selection based upon a deep neural-network are used. No significant excess
above the predicted background is found. In simplified models involving the
pair production of gluinos that decay via off-shell top (bottom) squarks,
gluino masses less than 2.44 TeV (2.35 TeV) are excluded at 95% CL for a
massless $\tilde\chi^0_1$. Limits are also set on the gluino mass in models
with variable branching ratios for gluino decays to $b\bar{b}\tilde\chi^0_1$,
$t\bar{t}\tilde\chi^0_1$ and $t\bar{b}\tilde\chi^-_1$ /
$\bar{t}b\tilde\chi^+_1$.
