# Compton scattering in QCD

### Theory

Let us consider Compton scattering in QCD. There are three Feynman diagrams:

Since the last diagram contains 3-gluon vertex and in order to avoid ghosts contributions, we will use QCD axial gauge:

$\mathcal L_{g.f.} = -\frac{1}{2\xi} \left( n_\mu A^{a \mu}\right)^2,$ where $n_\mu n^\mu = 1$. Further, we'll take $\xi = 0$. The dependence on the auxiliary unit vector $n_\mu$ must cancel in all observables.

The Feynman rules for QCD in axial gauge are: \begin{eqnarray} \mbox{quark propagator:}&\qquad& D(k) = \frac{-i\,\left(m + k_\mu \gamma^\mu \right)}{m^2 - k^2}, \\ \mbox{gluon propagator:}&\qquad& G_{\mu\nu}^{ab}(p) = -I\, \delta^{ab}\,\frac{P_{\mu\nu}(p)}{p^2}, \\ \mbox{quark-gluon vertex:} &\qquad& V_\mu^a = -i\, Q\, \gamma_ \mu\, T^a, \\ \mbox{3-gluon vertex:} &\qquad& V_{\mu\nu\rho}^{abc} (p_1, p_2, p_3) = \\ && = \, Q\, f^{abc}\,\left( g_{\mu\nu}\,(p_2 - p_1)_\rho + g_{\nu\rho}\,(p_3 - p_2)_\mu + g_{\mu\rho}\,(p_1-p_3)_\nu\right), \end{eqnarray} where $Q$ is quark charge, $\gamma_\mu$ is Dirac matrix, $T^a$ is SU(3) matrix and all momentums in 3-gluon vertex are incoming. Tensor $P_{\mu\nu}$ that appears in gluon propagator is equal to the sum of gluon polarzations: $\sum\, \epsilon_\mu(p)\, \epsilon_\nu(p) \,=\, P_{\mu\nu}(p) = -g_{\mu\nu} + \frac{p_\mu \, n_\nu + p_\nu \, n_\mu}{(p\,n)} - \frac{p_\mu \, p_\nu}{p^2}$

Using these Feynman rules it is easy to write amplitudes corresponding to the above Feynman diagrams:

\begin{eqnarray} \mbox{diagram a):}&\qquad& (\mathcal M_a)^{ab} =\bar u(k_2)\, V_\mu^a \, \epsilon^\mu (p_2) \, D(k_1 + p_1) \, V_\nu^b \, \epsilon^\nu (p_1) \, u(k_1), \\ \mbox{diagram b):}&\qquad& (\mathcal M_b)^{ab} =\bar u(k_2)\, V_\mu^b \, \epsilon^\mu (p_1) \, D(k_1 - p_2) \, V_\nu^a \, \epsilon^\nu (p_2) \, u(k_1), \\ \mbox{diagram c):}&\qquad& (\mathcal M_c)^{ab} =\bar u(k_2)\, V_\mu^c \, u(k_1) G_{\mu\nu}^{cd}(k_1 - k_2)\,V_{\nu\alpha\beta}^{dba}(k_1-k_2, p_1, -p_2) \epsilon^\alpha (p_1) \, \epsilon^\beta (p_2), \end{eqnarray} where $u$ and $\bar u$ are quark wave function and its conjugation respectively, $\epsilon$ is a polarisation vector of gluon and indices $a,\, b$ are colour indices of gluons. The total matrix element is $\mathcal M^{ab} = (\mathcal M_a)^{ab} \,+\, (\mathcal M_b)^{ab} \,+\, (\mathcal M_c)^{ab}$

The final goal is squared matrix element summed over final and averaged over initial quantum numbers:

$\frac{1}{32}\sum |\mathcal M|^2 = \frac{1}{32}\sum \mathcal M^{ab}\mathcal M^*_{ab}$

Sum over quark polarizations: \begin{gather} \sum\, u(k) \, \bar u(k)\, = \,m\, +\, k^\mu\,\gamma_\mu \end{gather}

The final standard thing is complex conjugation of matrix element. This can be done using the following obvious formula: \begin{multline} \left( \bar u(p_2)\, T_{a_1} \, T_{a_2}\,\dots \,T_{a_m} \, \gamma_{\alpha_1} \,\gamma_{\alpha_2}\,\dots \gamma_{\alpha_n} \, u(p_1) \right)^* = \\ = \bar u(p_1)\, T_{a_m} \, T_{a_{m-1}}\,\dots \,T_{a_1} \, \gamma_{\alpha_n} \,\gamma_{\alpha_{n-1}}\,\dots \gamma_{\alpha_1} \, u(p_2) \end{multline}

After summing squared matrix element over polarisations, all combinations of gamma and SU(N) matrices will automatically transform to combinations traces.

The dependence on the auxiliary unit vector $n_\mu$ must cancel in the final result for squared matrix element.

### Code

The following code reproduces exactly same steps as one need to perform with paper and pencil to calculate squared matrix element of Compton scattering in QCD (although it is hardly possible to solve this problem with paper and pencil):

//setting up matrix objects
//unitary matrices
defineMatrices 'T_A', Matrix2.matrix,
//gamma matrices
'G_a', Matrix1.matrix,
//quark wave function
'u[p_a]', Matrix1.vector, Matrix2.vector,
//its conjugation
'cu[p_a]', Matrix1.covector, Matrix2.covector,
//quark-gluon vertex
'V_iA', Matrix1.matrix, Matrix2.matrix,
//quark propagator
'D[p_m]', Matrix1.matrix

//setting up symmetries of tensors
//SU(N) symmetric constants
setSymmetric 'd_ABC'
//SU(N) structure constants
setAntiSymmetric 'f_ABC'

//Feynman rules
//quark-gluon vertex
def V = 'V_mA = -I*Q*G_m*T_A'.t
//quark propagator
def D = 'D[p_m] = -I*(m + p_m*G^m)/(m**2 - p_m*p^m)'.t
//sum over gluon polarizations in axial gauge
def P = '''P_mn[k_a] =
- g_mn + 1/(k_a*n^a)*(k_m*n_n + k_n*n_m)
- 1/(k_a*n^a)**2 * k_m*k_n'''.t
//auxiliary vector is unit
def n2 = 'n_a*n^a = 1'.t
//gluon propagator
def G = P >> 'G_mnAB[k_a] = I*g_AB*P_mn[k_a]/(k_a*k^a)'.t
//3-gluon vertex
def V3 = '''V_{mnr}^{ABC}[k1_m, k2_m, k3_m] = Q*f^{ABC}*(
g_mn*(k2_r - k1_r)
+ g_nr*(k3_m - k2_m)
+ g_mr*(k1_n-k3_n))'''.t

//Matrix element
//diagram a)
def Ma = 'cu[k2_m]*V_mA*e^m[p2_m]*D[p1_m+k1_m]*V_nB*e^n[p1_m]*u[k1_m]'.t
//diagram b)
def Mb = 'cu[k2_m]*V_mB*e^m[p1_m]*D[k1_m-p2_m]*V_nA*e^n[p2_m]*u[k1_m]'.t
//diagram c)
def Mc = '''cu[k2_m]*V_mC*u[k1_m]*G^mnCD[k1_a-k2_a]*
V_{nabDBA}[k1_a-k2_a, p1_m, -p2_m]*e^a[p1_m]*e^b[p2_m]'''.t
//matrix element
def M = Ma + Mb + Mc
//substitute Feynman rules
M = (D & G & V & V3) >> M
//Mandelstam and mass shell substitutions
def mandelstam = setMandelstam(
['k1_m': 'm', 'p1_m': '0', 'k2_m': 'm', 'p2_m': '0'])
//simplify matrix element
M = (ExpandAll & EliminateMetrics & n2 & mandelstam) >> M

//Complex conjugation of matrix element
def MC = Conjugate >> M
//exchange spinor momentums
MC = 'u[k1_m]*cu[k2_m] = u[k2_m]*cu[k1_m]'.t >> MC
//reorder gamma and SU(N) matrices
MC = (Reverse[Matrix1, Matrix2]) >> MC

//Squared matrix element
//set M_AB = matrix element
M = 'M_AB'.eq M
//set MC_AB = conjugated matrix element
MC = 'MC_AB'.eq MC
//squared matrix element
def M2 = (M & MC) >> 'M_AB*MC^AB/32'.t
//expand and eliminate contractions with metrics and deltas
M2 = (ExpandAndEliminate & n2) >> M2

//Sum over polarizations
//sum over gluon polarizations
M2 = ('e_m[p1_a]*e_n[p1_a] = P_mn[p1_a]'.t & P) >> M2
M2 = ('e_m[p2_a]*e_n[p2_a] = P_mn[p2_a]'.t & P) >> M2
//sum over quark polarizations
M2 = 'u[k2_m]*cu[k2_m] = m + k2^m*G_m'.t >> M2
M2 = 'u[k1_m]*cu[k1_m] = m + k1^m*G_m'.t >> M2

//trace of gamma matrices
M2 = DiracTrace['G_a'] >> M2
//trace of unitary matrices
M2 = UnitaryTrace >> M2
//simplifications and substitutions
M2 = (ExpandAndEliminate & n2 & mandelstam) >> M2
//simplify combinations of unitary constants
M2 = UnitarySimplify >> M2

//performing further simplifications
M2 = (ExpandAndEliminate & mandelstam) >> M2
//mandelstam u
M2 = 'u = 2*m**2-s-t'.t >> M2
//momentum conservation:
// using to reduce the number of scalar combinations
M2 = 'p2_a = k1_a + p1_a - k2_a'.t >> M2
M2 = ExpandAll >> M2
// replace scalar contractions with auxiliary vector n_m
// (like p1_m*n^m) with some symbols (e.g. p1_m*n^m = p1n)
['k1', 'k2', 'p1'].each {
M2 = "${it}_a *n^a =${it}n".t >> M2
}
println M2.toString(OutputFormat.WolframMathematica)


Let us consider the above code in more detail. In the first few lines we specify the matrix objects that will be used in further calculations. In contrast to Compton scattering in QED, quark wave function has both spinor and SU(N) matrix indices (with types Matrix1 and Matrix2 respectively). In the code we used Latin lowercase indices as Lorentz indices and uppercase indices as SU(N) indices.

Instead of the last few lines (100-105) we coul simply use Factor >> M2 to simplify and output the result. But, currently, the Factor transformation can work slowly when working with huge expressions. So, in this particular case it is better to output the result and paste it to some other system and perform simplifications there. For this, we replaced all remaining tensorial scalars (i.e. contractions of momentums with auxiliary vector $n_\mu$) with symbols and output the result in Wolfram Mathematica format. These symbols must cancel in the final result.

The factorized result is: \begin{multline} \frac{1}{32}\sum |\mathcal M|^2 =\frac{(-1 + N^2) Q^4 }{ 16\, N\, (m^2 - s)^2\, (m^2 - s - t)^2\, t^2} (2 m^4 N^2 - 4 m^2 N^2 s + 2 N^2 s^2 - 2 m^2 N^2 t + \\ + 2 N^2 s t - t^2 + N^2 t^2) \, \times \, (2 m^8 - 8 m^6 s + 12 m^4 s^2 - 8 m^2 s^3 + 2 s^4 + \\ + 4 m^4 s t - 8 m^2 s^2 t + 4 s^3 t + 3 m^4 t^2 - 2 m^2 s t^2 + 3 s^2 t^2 - m^2 t^3 + s t^3) \end{multline} where $N$ is a dimension of SU(N) matrices. As we see, the dependence on the auxiliary vector $n_\mu$ cancelled in the physical result.