  We perform a non-perturbative study of the scale-dependent renormalisation
factors of a complete set of dimension-six four-fermion operators. The
renormalisation-group (RG) running is determined in the continuum limit for a
specific Schr\"dinger Functional (SF) renormalisation scheme in the framework
of lattice QCD with two dynamical flavours ( $N_f = 2$ ). The theory is
regularised on a lattice with a plaquette Wilson action and
$\mathcal{O}(a)$-improved Wilson fermions. For one of these operators, the
computation had been performed in ref. [1]; the present work completes the
study for the rest of the operator basis, on the same simulations
(configuration ensembles). The related weak matrix elements arise in several
operator product expansions; in $\Delta F = 2$ transitions they contain the QCD
long-distance effects, including contributions from beyond-Standard Model (BSM)
processes. Some of these operators mix under renormalisation and their
RG-running is governed by anomalous dimension matrices. In ref. [2] the RG
formalism for the operator basis has been worked out in full generality and the
anomalous dimension matrix has been calculated in NLO perturbation theory. Here
the discussion is extended to the matrix step-scaling functions (matrix-SSFs),
which are used in finite-size recursive techniques. We rely on these
matrix-SSFs to obtain non-perturbative estimates of the operator anomalous
dimensions for scales ranging from $\mathcal{O}(\Lambda_{\rm QCD})$ to
$\mathcal{O}(M_W)$.
