BGap.h

   1 * This file implements the Berends-Giele recursion to compute the
   2 * alpha' expansion of disk integrals using the method described in the
   3 * paper arXiv:1609.XXXX
   4 *
   5 * Non-abelian Z-theory: Berends-Giele recursion for the alpha' expansion
   6 * of disk integrals
   7 *
   8 * in collaboration with Oliver Schlotterer.
   9 *
  10 * You can use this file as you wish, as long as you cite the
  11 * paper 1609.XXXX afterwards. In addition, consider citing
  12 * the paper 1603.09731 as well, as it is very closely related
  13 *
  14 * Copyright (c) 2016 Carlos R. Mafra
  15
  16         Indices a,b,c,d,e,f,g,h,i,j,k,l;
  17         Indices m,n,p,q,r,ss,t,u,v,x;
  18         Indices <m1,n1>,...,<m20,n20>;
  19         Indices i1,...,i20;
  20         Indices j1,...,j20;
  21
  22         Symbols [p],[c];
  23         Symbol x1;
  24         Symbols <zeta2>,...,<zeta10>;
  25
  26 * Mandelstam s(i,j), Is(i,j)=1/s(i,j), KK(i,[c],j) = ki.kj
  27         CFunction s(s), Is(s), KK;
  28
  29         CFunction phi, Zint;
  30         CFunction nargs;
  31
  32         CFunction ctmpF1;
  33         Function tmpF1;
  34
  35 * phi(i1,i2,i3|j1,j2,j3) = ML(i1)*ML(i2)*ML(i3)*MR(j1)*MR(j2)*MR(j3)
  36         Function ML,MR,L,R;
  37         Function KKL,KKR;
  38
  39 #procedure DeconcatenateN(F,n)
  40 *
  41 * Deconcatenate function F into n-partitions or all partitions when n==0
  42 *
  43 * n==0: F(1,2,3) --> F(1)*F(2)*F(3) + F(1,2)*F(3) + F(1)*F(2,3)
  44 * n==2: F(1,2,3) --> F(1,2)*F(3) + F(1)*F(2,3)
  45 *
  46         if (match(tmpF1(?m)));
  47                 print "Error: temp function already present: %t";
  48                 exit;
  49         endif;
  50
  51         if ((count('F',1) <= -1) || (count('F',1) >= 2));
  52                 print "Error: Deconc does not work for : %t";
  53                 exit;
  54         endif;
  55
  56         chainout 'F';
  57
  58         repeat id once 'F'(?m)*'F'(?n) =
  59         + 'F'(?m,?n)
  60         + tmpF1(?m)*'F'(?n)
  61         ;
  62
  63         id tmpF1(?m) = 'F'(?m);
  64
  65         #if 'n' == 0
  66 * do nothing, keep all partitions
  67         #else
  68 * keep only n-partitions (or up to n-partitions)
  69 *        if ((match('F'(?m))) && (count('F',1) != 'n'));
  70         if ((match('F'(?m))) && (count('F',1) > 'n'));
  71                 discard;
  72         endif;
  73         #endif
  74
  75 * remove trivial partition F(1,2,3) -> F(1,2,3)
  76         if ((count('F',1) == 1));
  77                 discard;
  78         endif;
  79
  80 #endprocedure
  81
  82
  83 #procedure BiDeconcatenate(phi,n)
  84 *
  85 * phi(A,[p],B) --> all deconcatenations of A and B, separated by [c]
  86 *
  87 * s(1,2,3)*phi(1,2,3,[p],4,5,6) = phi(1,[c],2,3,[p],4,[c],5,6) + phi(1,2,[c],3,[p],4,[c],5,6) + ...
  88
  89 * A and B cannot be composed of smaller words, so no [c]s allowed
  90         if (match(phi(?m,[c],?n)));
  91                 print "Bug: cannot deconcatenate multiple slots %t";
  92                 exit;
  93         endif;
  94
  95         id once 'phi'(?m1,[p],?m2) = Is(?m1)*ML(?m1)*MR(?m2);
  96         #call DeconcatenateN(ML,'n')
  97         #call DeconcatenateN(MR,'n')
  98
  99 * keep only the same depth of deconcatenation on both sides
 100         if (count(ML,1) != count(MR,1));
 101         discard;
 102         endif;
 103
 104         .sort:bideconc'n';
 105 #endprocedure
 106
 107
 108 #procedure tLalphaExpand(n,w)
 109
 110         #if {'$Npts'-1-'n'} >= 9
 111         #include- BGphiN/BGphi9-ordered.h
 112         bracket L,R;
 113         .sort:brack9;
 114         keep brackets;
 115         id L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*L(?a8)*L(?a9)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)*R(?b8)*R(?b9) =
 116         + phi(?a9,[p],?b9)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*L(?a8)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)*R(?b8)
 117         - phi(?a9,[p],?b1)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*L(?a8)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)*R(?b8)*R(?b9)
 118         ;
 119         #call MultiparticleConstraint(9)
 120         #endif
 121
 122         #if {'$Npts'-1-'n'} >= 8
 123         #include- BGphiN/BGphi8-ordered.h
 124         bracket L,R;
 125         .sort:brack8;
 126         keep brackets;
 127         id L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*L(?a8)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)*R(?b8) =
 128         + phi(?a8,[p],?b8)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)
 129         - phi(?a8,[p],?b1)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)*R(?b8)
 130         ;
 131         #call MultiparticleConstraint(8)
 132         #endif
 133
 134         #if {'$Npts'-1-'n'} >= 7
 135         #include- BGphiN/BGphi7-ordered.h
 136         bracket L,R;
 137         .sort:brack7;
 138         keep brackets;
 139         id L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*L(?a7)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7) =
 140         + phi(?a7,[p],?b7)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)
 141         - phi(?a7,[p],?b1)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)*R(?b7)
 142         ;
 143         #call MultiparticleConstraint(7)
 144         #endif
 145
 146         #if {'$Npts'-1-'n'} >= 6
 147         #include- BGphiN/BGphi6-ordered.h
 148         bracket L,R;
 149         .sort:brack6;
 150         keep brackets;
 151         id L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*L(?a6)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6) =
 152         + phi(?a6,[p],?b6)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5)
 153         - phi(?a6,[p],?b1)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*R(?b2)*R(?b3)*R(?b4)*R(?b5)*R(?b6)
 154         ;
 155         #call MultiparticleConstraint(6)
 156         #endif
 157
 158         #if {'$Npts'-1-'n'} >= 5
 159         #include- BGphiN/BGphi5-ordered.h
 160         bracket L,R;
 161         .sort:brack5;
 162         keep brackets;
 163         id L(?a1)*L(?a2)*L(?a3)*L(?a4)*L(?a5)*R(?b1)*R(?b2)*R(?b3)*R(?b4)*R(?b5) =
 164         + phi(?a5,[p],?b5)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*R(?b1)*R(?b2)*R(?b3)*R(?b4)
 165         - phi(?a5,[p],?b1)*L(?a1)*L(?a2)*L(?a3)*L(?a4)*R(?b2)*R(?b3)*R(?b4)*R(?b5)
 166         ;
 167         #call MultiparticleConstraint(5)
 168         #endif
 169
 170         #if {'$Npts'-1-'n'} >= 4
 171         #include- BGphiN/BGphi4-ordered.h
 172         bracket L,R;
 173         .sort:brack4;
 174         keep brackets;
 175         id L(?a1)*L(?a2)*L(?a3)*L(?a4)*R(?b1)*R(?b2)*R(?b3)*R(?b4) =
 176         + phi(?a4,[p],?b4)*L(?a1)*L(?a2)*L(?a3)*R(?b1)*R(?b2)*R(?b3)
 177         - phi(?a4,[p],?b1)*L(?a1)*L(?a2)*L(?a3)*R(?b2)*R(?b3)*R(?b4)
 178         ;
 179         #call MultiparticleConstraint(4)
 180         #endif
 181
 182         #if {'$Npts'-1-'n'} >= 3
 183         #include- BGphiN/BGphi3-ordered.h
 184         bracket L,R;
 185         .sort:brack3;
 186         keep brackets;
 187         id L(?a1)*L(?a2)*L(?a3)*R(?b1)*R(?b2)*R(?b3) =
 188         + phi(?a3,[p],?b3)*L(?a1)*L(?a2)*R(?b1)*R(?b2)
 189         - phi(?a3,[p],?b1)*L(?a1)*L(?a2)*R(?b2)*R(?b3)
 190         ;
 191         #call MultiparticleConstraint(3)
 192         #endif
 193
 194         bracket L,R;
 195         .sort:brack2;
 196         keep brackets;
 197         id L(?a1)*L(?a2)*R(?b1)*R(?b2) =
 198         + phi(?a2,[p],?b2)*phi(?a1,[p],?b1)
 199         - phi(?a2,[p],?b1)*phi(?a1,[p],?b2)
 200         ;
 201         #call MultiparticleConstraint(2)
 202
 203         id ML(?m1)*ML(?m2)*MR(?n1)*MR(?n2) =
 204         + phi(?m1,[p],?n1)*phi(?m2,[p],?n2)
 205         - phi(?m2,[p],?n1)*phi(?m1,[p],?n2)
 206         ;
 207
 208         bracket phi;
 209         .sort:brackphi;
 210         keep brackets;
 211
 212         #call MultiparticleConstraint(1)
 213
 214 #endprocedure
 215
 216
 217 #procedure MultiparticleConstraint(n)
 218
 219         if (match(phi(?m,[c],?n)));
 220         print "bug: no subwords in phi %t";
 221         exit;
 222         endif;
 223
 224         id phi(i?,[p],j?) = delta_(i,j);
 225         id phi(i?,j?,[p],i?,j?) = Is(i,j);
 226         id phi(i?,j?,[p],j?,i?) = - Is(i,j);
 227
 228 * remove phi(A|B) where length(A) != length(B)
 229         id phi(?m,[p],?n) = phi(?m,[p],?n)*delta_(nargs_(?m),nargs_(?n));
 230
 231 * remove phi(A|B) when letters in A and B are not the same
 232         if (match(ctmpF1(?m))==0);
 233         id phi(?m,[p],?n) = phi(?m,[p],?n)*ctmpF1(?m,[p],?n);
 234         repeat id ctmpF1(?m1,i?,?m2,[p],?m3,i?,?m4) = ctmpF1(?m1,?m2,[p],?m3,?m4);
 235         id ctmpF1([p]) = 1;
 236         id ctmpF1(?m) = 0;
 237         endif;
 238
 239         .sort: MulCtrnt'n';
 240
 241 #endprocedure
 242
 243 #procedure BGexpand(w)
 244
 245         .sort:apExp'w';
 246         Symbol zwgt(:'w');
 247
 248         #if 'w' > 7
 249         print "Error: order higher than alpha^7 not implemented!";
 250         exit;
 251         #endif
 252
 253 * Sanity check: number of arguments on both sides must be equal
 254         id Zint(?m,[p],?n) = Zint(?m,[p],?n)*nargs(nargs_(?m) - nargs_(?n));
 255         id nargs(0) = 1;
 256         if (match(nargs(i?)));
 257         print "Error: unequal number of labels in domain and integrand: %t";
 258         exit;
 259         endif;
 260
 261         if (match(Zint(?m,[p],?n)) == 0);
 262         print "Error: specify the domain P and integrand Q as Zint(P,[p],Q)";
 263         exit;
 264         endif;
 265
 266 * Get the number of points
 267         id once Zint(?m,[p],?n) = Zint(?m,[p],?n)*nargs(nargs_(?n));
 268         id once nargs(x1?$Npts) = 1;
 269
 270         #define Npts "$Npts"
 271
 272 * Exploit cyclic symmetry of Zint(P|Q) in Q to match last label
 273         id once Zint(?m,i?,[p],?n,i?,?p) = s(?n,?p)*phi(?m,[p],?p,?n);
 274         .sort
 275
 276 * TODO: the number of calls to BiDeconcatenate(phi,i) should vary according to w
 277         #call BiDeconcatenate(phi,{'w'+2})
 278         #call CancelMandelstams('Npts')
 279
 280         #call tLalphaExpand(0,'w')
 281         #call BiDeconcatenate(phi,{'w'+2})
 282         #call tLalphaExpand(1,'w')
 283         #call BiDeconcatenate(phi,{'w'+2})
 284         #call tLalphaExpand(2,'w')
 285         #call BiDeconcatenate(phi,{'w'+2})
 286         #call tLalphaExpand(3,'w')
 287         #call BiDeconcatenate(phi,{'w'+2})
 288         #call tLalphaExpand(4,'w')
 289         #call BiDeconcatenate(phi,{'w'+2})
 290         #call tLalphaExpand(5,'w')
 291         #call BiDeconcatenate(phi,{'w'+2})
 292         #call tLalphaExpand(6,'w')
 293         #call BiDeconcatenate(phi,{'w'+2})
 294         #call tLalphaExpand(7,'w')
 295
 296         if (match(phi(?m)));
 297         print "Bug: cannot have phi() at the end of recursion: %t";
 298         print "Add more calls to tLalphaExpand()";
 299         exit;
 300         endif;
 301
 302         id zwgt = 1;
 303
 304         #call KKTos()
 305
 306 #endprocedure
 307
 308 #procedure KKTos()
 309
 310         id KK(i?,[c],j?) = s(i,j);
 311
 312         bracket KK;
 313         .sort:bracKK;
 314         keep brackets;
 315
 316         repeat id KK(?m,[c],?n,i?,j?) = KK(?m,[c],?n,i) + KK(?m,[c],j);
 317         repeat id KK(?n,i?,j?,[c],?m) = KK(?m,[c],?n,i) + KK(?m,[c],j);
 318         repeat id KK(?m,[c],?n,i?,j?) = KK(?m,[c],?n,i) + KK(?m,[c],j);
 319         repeat id KK(?n,i?,j?,[c],?m) = KK(?m,[c],?n,i) + KK(?m,[c],j);
 320         repeat id KK(?m,[c],?n,i?,j?) = KK(?m,[c],?n,i) + KK(?m,[c],j);
 321         repeat id KK(?n,i?,j?,[c],?m) = KK(?m,[c],?n,i) + KK(?m,[c],j);
 322
 323         id once KK(i?,[c],i?) = 0;
 324         id once s(i?,i?) = 0;
 325         id KK(i?,[c],j?) = s(i,j);
 326         id s(i?,j?)*Is(i?,j?) = 1;
 327
 328         .sort
 329 *        #call KNBasis('$Npts')
 330 *        #call MomentumConservationMandelstams('$Npts')
 331 *        #call NumericMandelstams('$Npts')
 332
 333 #endprocedure
 334
 335
 336 #ifndef 'pss2mandelstams'
 337
 338 #procedure CancelMandelstams(N)
 339
 340         #do i=2,'N'
 341         id s(i1?,...,i'i'?)^-1 = Is(i1,...,i'i');
 342         id s(i1?,...,i'i'?)*Is(i1?,...,i'i'?) = 1;
 343         #enddo
 344
 345 #endprocedure
 346
 347 #endif