-- Verify Casimir trick for electron positron annihilation E = sqrt(p^2 + m^2) p1 = (E, 0, 0, p) p2 = (E, 0, 0, -p) p3 = (E, E expsin(theta) expcos(phi), E expsin(theta) expsin(phi), E expcos(theta)) p4 = (E, -E expsin(theta) expcos(phi), -E expsin(theta) expsin(phi), -E expcos(theta)) u11 = (E + m, 0, p1[4], p1[2] + i p1[3]) u12 = (0, E + m, p1[2] - i p1[3], -p1[4]) v21 = (p2[4], p2[2] + i p2[3], E + m, 0) v22 = (p2[2] - i p2[3], -p2[4], 0, E + m) u1 = (u11,u12) v2 = (v21,v22) N = (E + m)^2 I = ((1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1)) gmunu = ((1,0,0,0),(0,-1,0,0),(0,0,-1,0),(0,0,0,-1)) gamma0 = ((1,0,0,0),(0,1,0,0),(0,0,-1,0),(0,0,0,-1)) gamma1 = ((0,0,0,1),(0,0,1,0),(0,-1,0,0),(-1,0,0,0)) gamma2 = ((0,0,0,-i),(0,0,i,0),(0,i,0,0),(-i,0,0,0)) gamma3 = ((0,0,1,0),(0,0,0,-1),(-1,0,0,0),(0,1,0,0)) gamma = (gamma0,gamma1,gamma2,gamma3) gammaT = transpose(gamma) gammaL = transpose(dot(gmunu,gamma)) v2bar = dot(conj(v2),gamma0) -- adjoint of v2 q1 = p1 - p3 q2 = p1 - p4 qslash1 = dot(q1,gmunu,gamma) qslash2 = dot(q2,gmunu,gamma) t = dot(q1,gmunu,q1) u = dot(q2,gmunu,q2) "Sum over spin states" a11 = 0 a12 = 0 a22 = 0 for(a,1,2,for(b,1,2, a1 = dot(v2bar[b],gammaT,qslash1 + m I,gammaT,u1[a]), a2 = dot(v2bar[b],gammaT,qslash2 + m I,gammaT,u1[a]), a11 = a11 + contract(dot(a1,gmunu,transpose(conj(a1)),gmunu)), a12 = a12 + contract(dot(a1,gmunu,conj(a2),gmunu)), a22 = a22 + contract(dot(a2,gmunu,transpose(conj(a2)),gmunu)) )) a11 a12 a22 "Casimir trick" pslash1 = dot(p1,gmunu,gamma) pslash2 = dot(p2,gmunu,gamma) P1 = pslash1 + m I P2 = pslash2 - m I Q1 = qslash1 + m I Q2 = qslash2 + m I T = dot(P1,gammaT,Q1,gammaT,P2,gammaL,Q1,gammaL) f11 = contract(contract(contract(T,3,4),2,3)) T = dot(P1,gammaT,Q2,gammaT,P2,gammaL,Q1,gammaL) f12 = contract(contract(contract(T,3,5),2,3)) T = dot(P1,gammaT,Q2,gammaT,P2,gammaL,Q2,gammaL) f22 = contract(contract(contract(T,3,4),2,3)) f11 f12 f22 "Verify Casimir trick (1=ok)" f11 == 1/N a11 f12 == 1/N a12 f22 == 1/N a22