  We present a general formalism for computing the matrix-element squared for
the emission of soft energy-ordered gluons beyond two loops in QCD perturbation
theory at finite $N_c$. Our formalism is valid in the eikonal approximation. A
Mathematica program has been developed for the automated calculation of all
real/virtual eikonal squared amplitudes needed at a given loop order. For the
purpose of illustration we show the explicit forms of the eikonal squared
amplitudes up to the fifth-loop order. In the large-$N_c$ limit our results
coincide with those previously reported in literature.
