  The correlator of two gluonic operators plays an important role for example
in transport properties of a Quark Gluon Plasma (QGP) or in sum rules for
glueballs.
  In [1] an operator product expansion (OPE) at zero temperature was performed
for the correlators of two scalar operators $O_1=-\frac{1}{4} G^{\mu \nu}G_{\mu
\nu}$ and two QCD energy-momentum tensors $T^{\mu\nu}$. There we presented
analytical two-loop results for the Wilson coefficient $C_1$ in front of the
gluon condensate operator $O_1$. In this paper these results are extended to
three-loop order.
  The three-loop Wilson coefficient $C_0$ in front of the unity operator
$O_0=1$ was already presented in [1] for the $T^{\mu\nu}$-correlator. For the
$O_1$-correlator the coefficient $C_0$ is known to four loop order from [2].
For the correlator of two pseudoscalar operators
$\tilde{O}_1=\varepsilon_{\mu\nu\rho\sigma} G^{\mu \nu} G^{\rho \sigma}$ both
coefficients $C_0$ and $C_1$ were computed in [3] to three-loop order. At zero
temperature $C_0$ and $C_1$ are the leading Wilson coefficients in massless
QCD.
