/* [wxMaxima batch file version 1] [ DO NOT EDIT BY HAND! ]*/
/* [ Created with wxMaxima version 17.10.0 ] */
/* [wxMaxima: title   start ]
compton0.wxm   

    π(+, p1) +  γ(k1, λ1)   -->  π(+, p2)  + γ(k2, λ2)
   [wxMaxima: title   end   ] */


/* [wxMaxima: comment start ]
compton0.wxm,   Jan. 10, 2018
Edwin (Ted) Woollett, Maxima by Example,  
ch. 12,   Dirac Algebra and Quantum Electrodynamics, ver. 3, "Dirac3"

http://web.csulb.edu/~woollett/ 
woollett@charter.net

compton0.wxm evaluates Compton scattering in the rest frame
of the initial pion (the "lab frame").

 We use Heaviside-Lorentz electromagnetic units.
In general, we follow the conventions in Peskin & Schroeder: Quantum Field Theory.
In particular, the diagonal elements of the metric tensor g^{μ, ν} are diag (g)  = (1, -1, -1, -1).

Following Griffiths and Halzen/Martin, we interpret the Feynman rules
as giving an expression for - i M.
   [wxMaxima: comment end   ] */


/* [wxMaxima: comment start ]
First load the dirac3 package by loading dirac3.mac.

Your maxima-init.mac file should have the lines (if you are working
in c:/work5/ and the dirac3 package files are present in that folder)
( See Ch. 1 of Maxima by Example for a more leisurely discussion)

maxima_userdir: "c:/work5" $
maxima_tempdir : "c:/work5"$
file_search_maxima : append(["c:/work5/###.{mac,mc}"],file_search_maxima )$
file_search_lisp : append(["c:/work5/###.lisp"],file_search_lisp )$

   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
load(dirac3);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: section start ]
Amplitude for Specified Photon Polarizations λ1 and λ2
   [wxMaxima: section end   ] */


/* [wxMaxima: comment start ]
We write the amplitude in terms of the symbolic 4-vector product symbol D.
D( pa, pb )  = D ( pb, pa ) is the symbolic 4-vector product of pa and pb.

[As a reference, see Peter Renton,  Electroweak Interactions,  Fig. 4.10 (page 181)
and Eq (4.144) (page 185); translation:  p --> p1,  k --> k1,  p' --> p2,  k' --> k2,
  ε  -->  ε1   ε'  -->  ε2,  q'  --> qp. ]

With our conventions, the amplitude for specific  photon polarizations  λ1, λ2  is

M( λ1, λ2) = -e^2 Mr,
 
and the "reduced amplitude"  Mr = M1 + M2 + M3, where
-------------------------------------------------------------------------
M1 =   D ( ε1, p1 + q)  * D( ε2*, p2 + q ) / (s - m^2),

   where  
   m is the mass of the pion,  
   s = (p1 + k1)^2, 
   q = p1 + k1 = p2 + k2  by conservation of 4-momentum,  
   ε1 = ε (k1, λ1), 
   ε2 = ε (k2, λ2) 
   ε2* is the complex conjugate of  ε2.

M2 =   D ( ε2*, p1 + qp)  * D( ε1, p2 + qp ) / (u - m^2),

    where u = (k2 - p1)^2 = (p1 - k2)^2, qp (q-prime) = p1 - k2 = p2 - k1.

M3 =  - 2 * D ( ε2*, ε1 ). 

The photon polarization 4-vectors ε1 and ε2 obey the Lorentz
 gauge conditions

    D ( k1, ε1 ) = 0,
    D ( k2, ε2 ) = 0,
    D ( k2, ε2* ) = 0.

which imply that M1 and M2 can be simplified to

  M1 = 4 * D ( ε1, p1 )  *  D ( ε2*, p2 ) / ( s - m^2)

   M2 = 4 * D ( ε2*, p1 )  *  D ( ε1, p2 ) / ( u - m^2)

We use gauge freedom within the Lorentz gauge to use 
  a "transverse gauge", such that  ε1^{0} = 0 and ε2^{0} = 0.

We now specialize to the lab frame of the initial pion, in which
p1 = [ m, 0, 0, 0],  so D ( p1, ε1 ) = 0, and  D ( p1, ε2* ) = 0. 

Thus in the lab frame of the initial pion, we have

  Mr ( λ1, λ2)   =   M3   =  - 2 * D ( ε2*, ε1 ). 

Using explicit circular polarization vectors, we can use this to
evaluate the possible (reduced) helicity amplitudes for R --> R,  R --> L,
L --> R, and L --> L  (vide infra).

   [wxMaxima: comment end   ] */


/* [wxMaxima: comment start ]
For specified initial and final polarization states of the photons, we can then write

dsigdo = A_lab * | Mr ( λ1, λ2) |^2  = 4 * A_lab * | D ( ε2*, ε1 ) |^2

where 
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
A_lab :  (alpha^2 / (4*m^2) ) * (kp/k)^2;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Let ek1 stand for ε1 and ek1c stand for ε1*
Likewise, let ek2 stand for ε2 and ek2c stand for ε2*
then,
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
dsigdo_lab_pol : A_lab * 4 * D(ek2c,ek1) * D(ek2, ek1c);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: section start ]
Use comp_def(...) to define the 4-momenta components
   [wxMaxima: section end   ] */


/* [wxMaxima: comment start ]
These define the 4-momenta components in the initial pion rest frame.

All the particle 3-momentum vectors lie in the z-x plane.

comp_def ( p1 (E1, px, py, pz), ...) for example assigns to the particle
with 4-momentum p1 the energy E1, and the 3-momentum components the
values px, py, and pz.

th, standing for θ, is the angle between the emergent photon k2's
3-momentum vector and the positive z-axis.

We write the components of the emerging pion's 4-momentum  p2 using
conservation of 4-momentum: p1 + k1 = p2 + k2.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
assume ( k > 0, kp > 0, th >= 0, th <= %pi )$

comp_def (
            p1( m, 0, 0, 0),
            k1 (k, 0, 0, k),
            k2 (kp, kp*sin(th), 0, kp*cos(th)),
            p2( m + k - kp, - kp*sin(th), 0, k - kp*cos(th)) )$
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
comp_def creates a separate (zero based)  array for each 4-vector.
Use listarray to see the components
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
listarray (k2);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
k2 [0];
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
VP(...) is our explicit (rather than symbolic)  "4-vector product"  function.  For example
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
VP(pa,pa);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
grind (%)$
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
From the definition of s and conservation of 4-momentum, we have

s = (p1 + k1)^2 = m^2 + 2 * D(p1, k1)

s = (p2 + k2)^2 = m^2 + 2 * D(p2, k2),

so D(p1,k1) = D(p2,k2).

We exploit this equality to solve for kp in terms of k, using Maxima's solve function.

   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
eqn : VP(p1,k1) = VP(p2,k2);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
solns : solve (eqn,kp);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
kp_soln : rhs (solns[1]);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
 So we can write the ratio  (kp/k) as   F (m, k, th) in the form:

  (kp/k) =  1 / [ 1 + (k/m) * ( 1 - cos(th)) ]
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
kp_k_ratio : kp/k, solns;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
The minus sign out front goes with the symbol m on top, if we
play with denom and num functions. Let's try to force the minus
sign downstairs and collect terms proportional to k/m.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
rd :  -denom(kp_k_ratio);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
expand (  %/m);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
The dirac3 package function collect_terms (expr, fac) is used here.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
collect_terms (%, k/m);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
kp_k_ratio :  1 / %;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
which achieves our goals for this expression, the ratio of k'/k == kp/k as
a function of k, m, and cos(th)
   [wxMaxima: comment end   ] */


/* [wxMaxima: section start ]
dσ/dΩ_lab for specific helicities
   [wxMaxima: section end   ] */


/* [wxMaxima: comment start ]
Above, we have the expression

dsigdo_lab_pol : A_lab * 4 * D(ek2c,ek1) * D(ek2, ek1c);


   [wxMaxima: comment end   ] */


/* [wxMaxima: comment start ]
In photon2.wxm we derived circular polarization 4-vectors in the "transverse gauge", 
for various cases, and we use those results here to explicitly evaluate the
factor  | D ( ε2*, ε1 ) |^2 which appears above in dsigdo_lab_pol.

Let the circular polarization 4-vector ek1R  describe a right-handed photon whose
3-momentum vector is aligned along the positive z-axis. Likewise, let the circular polarization 
4-vector ek1L  describe a left-handed photon whose 3-momentum vector is aligned along the 
positive z-axis.  Also ek1Rc is the complex-conjugate of ek1R, etc.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
arrays;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
comp_def (
     ek1R (0, -1/sqrt(2),  -%i/sqrt(2), 0 ),
     ek1L (0,  1/sqrt(2),  -%i/sqrt(2),  0),
     ek1Rc (0, -1/sqrt(2),  %i/sqrt(2), 0 ),
     ek1Lc  (0,  1/sqrt(2),  %i/sqrt(2),  0) )$
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
arrays;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Let the circular polarization 4-vector ek2R  describe a right-handed photon whose
3-momentum vector is parallel to  [sin(th), 0, cos(th)]. Likewise, let the circular polarization 
4-vector ek2L  describe a left-handed photon whose 3-momentum vector is parallel to  
[sin(th), 0, cos(th)].  Also ek2Rc is the complex-conjugate of ek2R, etc.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
comp_def ( ek2R (0, - cos (th)/sqrt(2), - %i/sqrt(2), sin(th)/sqrt(2) ),
                   ek2L (0, cos (th)/sqrt(2), - %i/sqrt(2), - sin(th)/sqrt(2)),
                   ek2Rc (0, - cos (th)/sqrt(2), %i/sqrt(2), sin(th)/sqrt(2) ),
                   ek2Lc (0, cos (th)/sqrt(2), %i/sqrt(2), - sin(th)/sqrt(2)) )$
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
We can then use VP(a,b) to evaluate the 4-vector dot product (with our metric conventions).
D ( ε2*, ε1 ) is then given by VP (ek2c, ek1) for each of the four helicity cases. These will all turn
out to be real, so we can then square the result to get | D ( ε2*, ε1 ) |^2, which we call Dsq.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
VP (ek2Rc, ek1R);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
ratsimp(%);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Case R --> R
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
Dsq_R_R : ratsimp ( VP (ek2Rc, ek1R) )^2;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
dsig_do_lab_RR : A_lab * 4 * Dsq_R_R;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Case R --> L
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
Dsq_R_L : ratsimp ( VP (ek2Lc, ek1R) )^2;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
dsig_do_lab_RL : A_lab * 4 * Dsq_R_L;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Case L --> R
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
Dsq_L_R : ratsimp ( VP (ek2Rc, ek1L) )^2;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
dsig_do_lab_LR : A_lab * 4 * Dsq_L_R;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Case L --> L
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
Dsq_L_L : ratsimp ( VP (ek2Lc, ek1L) )^2;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
dsig_do_lab_LL : A_lab * 4 * Dsq_L_L;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: section start ]
Unpolarized dσ/dΩ_lab 
   [wxMaxima: section end   ] */


/* [wxMaxima: comment start ]
Since the reduced amplitude for the case  λ1 -->  λ2  is
Mr ( λ1, λ2)   =   M3   =  - 2 * D ( ε2*, ε1 ) we get the unpolarized
differential scattering cross section by averaging over the 
initial polarization λ1 and summing over the final polarization  λ2  the
quantity 4 times the product of D ( ε2*, ε1 ) with D ( ε2, ε1* )

We can then use 

Σ(λ1)  ε1*^{nu}  ε1^{mu} = - g^{nu, mu} = - g^{mu,nu}

and likewise

Σ(λ2)  ε2*_{mu}  ε2_{nu} =  - g_{mu,nu}

and then Σ (mu,nu)  g^{mu,nu} g_{mu,nu} = 4.

   [wxMaxima: comment end   ] */


/* [wxMaxima: comment start ]
Our symbolic metric tensor is Gm (mu, nu); explicit values can
be obtained using noncov.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
Gm(0,0);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: input   start ] */
noncov(%);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Symbolic contraction of repeated indices can be accomplished
using either Con or scon.  
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
Con ( Gm(mu,nu)*Gm(mu,nu));
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Our explicit metric tensor is the array gmet[mu,nu].
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
gmet[0,0];
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
Contractions involving explicit metric tensors and/or matrices is
carried out using mcon, and we must supply the contraction indices.
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
mcon (gmet[mu,nu] * gmet[mu,nu], mu, nu);
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
So Σ (λ1, λ2)  | Mr (λ1, λ2) |^2 = 16, and 
   [wxMaxima: comment end   ] */


/* [wxMaxima: input   start ] */
dsigdo_lab_unpol : (A_lab/2) * 16;
/* [wxMaxima: input   end   ] */


/* [wxMaxima: comment start ]
which is proportional to the angle dependent factor (kp/k)^2.
   [wxMaxima: comment end   ] */



/* Maxima can't load/batch files which end with a comment! */
"Created with wxMaxima"$