diff --git a/check/diagrams/compare.frm b/check/diagrams/compare.frm new file mode 100644 index 00000000..766c43cf --- /dev/null +++ b/check/diagrams/compare.frm @@ -0,0 +1,1060 @@ +* Compares the diagram generator with others. +#ifndef `TEST' + #message Use -D TEST=XXX + #terminate +#else + #include `NAME_' # `TEST' +#endif +.end + +*--#[ diagram_compare_include : + +*#define Verbose +*#define CacheRead +*#define CacheWrite +#ifndef `CacheDir' + #define CacheDir "." +#endif + +#ifdef `TESTFILEDIR' + #redefine CacheRead + #redefine CacheDir "`TESTFILEDIR'/cache" +#endif + +#define MaxLegs "9" +#define MaxPropagators "19" + +#define CurrentModel "" +#define Bosons "" +#define AntiParticles "" + +#define DiagramGenerator "" +#define DiagramGeneratorVersion "" + +#define QgrafVersion "" +#define FeynGraphVersion "" + +#define QgrafStyleFile "qgraf.sty" +#define QgrafModelFile "qgraf.mdl" +#define QgrafInputFile "qgraf.dat" +#define QgrafOutputFile "qgraf.out" +#define FeynGraphTemplateFile "feyngraph.jinja" +#define FeynGraphInputFile "feyngraph.toml" +#define FeynGraphOutputFile "diagrams.out" + +Format nospaces; +Format 255; + +* Momentum convention: all incoming +* incoming particles: p1, p2, ... +* outgoing particles: -q1, -q2, ... +* internal particles: k1, k2, ... + +S n1,n2,x1,x2; +CF topo,node,counter,vx,vxs(symmetric),acc,replace; +V p1,...,p`MaxLegs',q1,...,q`MaxLegs',k1,...,k`MaxPropagators',l1,...,l`MaxPropagators'; +Set kk:k1,...,k`MaxPropagators'; +Set ll:l1,...,l`MaxPropagators'; + +#OpenDictionary wildmom + #do i = 1, `MaxPropagators' + #add k`i': "k`i'?$k`i'" + #enddo +#CloseDictionary + +#do i = 1, `MaxPropagators' + #$k`i' = 0; +#enddo + +Model TEMP; + Particle f1; + Particle f2; +EndModel; + +*--#[ SetupPhi3Model : + +** +* Prepares the PHI3 model. +* +#procedure SetupPhi3Model() + #define filename "`QgrafModelFile'" + #create <`filename'> + #write <`filename'> "[phi,phi,+]" + #write <`filename'> "[phi,phi,phi]" + #close <`filename'> + + Model PHI3; + Particle phi,1; + Vertex phi,phi,phi:g; + EndModel; + + #redefine CurrentModel "PHI3" + #redefine Bosons "phi" + #redefine AntiParticles "" +#endprocedure + +*--#] SetupPhi3Model : +*--#[ SetupPhi4Model : + +** +* Prepares the PHI4 model. +* +#procedure SetupPhi4Model() + #define filename "`QgrafModelFile'" + #create <`filename'> + #write <`filename'> "[phi,phi,+]" + #write <`filename'> "[phi,phi,phi,phi]" + #close <`filename'> + + Model PHI4; + Particle phi,1; + Vertex phi,phi,phi,phi:g^2; + EndModel; + + #redefine CurrentModel "PHI4" + #redefine Bosons "phi" + #redefine AntiParticles "" +#endprocedure + +*--#] SetupPhi4Model : +*--#[ SetupQcdModel : + +** +* Prepares the QCD model. +* +#procedure SetupQcdModel() + #define filename "`QgrafModelFile'" + #create <`filename'> + #write <`filename'> "[qua,QUA,-]" + #write <`filename'> "[glu,glu,+]" + #write <`filename'> "[gho,GHO,-]" + #write <`filename'> "[QUA,qua,glu]" + #write <`filename'> "[glu,glu,glu]" + #write <`filename'> "[glu,glu,glu,glu]" + #write <`filename'> "[GHO,gho,glu]" + #close <`filename'> + + Model QCD; + Particle qua,QUA,-2; + Particle gho,GHO,-1; + Particle glu,+3; + Vertex QUA,qua,glu:g; + Vertex GHO,gho,glu:g; + Vertex glu,glu,glu:g; + Vertex glu,glu,glu,glu:g^2; + EndModel; + + #redefine CurrentModel "QCD" + #redefine Bosons "glu" + #redefine AntiParticles "QUA,GHO" +#endprocedure + +*--#] SetupQcdModel : +*--#[ MakeDiagrams : + +** +* Manipulates diagrams for further manipulation. +* +* Parameters +* ---------- +* numIn: positive integer +* The number of incoming particles. +* numOut: positive integer +* The number of outgoing particles. +* +* Example Input +* ------------- +* topo_(1) +* *node_(1,1,QUA(-p1)) +* *node_(2,1,qua(-q1)) +* *node_(3,g,QUA(k2),qua(p1),glu(k1)) +* *node_(4,g,QUA(q1),qua(-k2),glu(-k1)) +* +* vx(QUA(-2),qua(1),glu(2)) +* *vx(QUA(1),qua(-1),glu(2)) +* +* Example Output +* -------------- +* node(QUA(q1),qua(-k2),glu(-k1)) +* *node(QUA(k2),qua(p1),glu(k1)) +* +#procedure MakeDiagrams(numIn,numOut) +* for diagrams_ + multiply replace_(topo_,dummy_); + id node_(n1?,n2?,?a) = node_(?a); + id node_(f1?(p1?)) = 1; + multiply replace_(node_,node); + +* for QGRAF +* Incoming fields have indices -1, -3, -5, ..., - 2j+1, ... +* Outgoing fields have indices -2, -4, -6, ..., - 2j, ... +* See the manual of QGRAF section 5.1 "The indices". + #do i = 1, `numIn' + repeat id vx(?a,f1?(1-2*`i'),?b) = vx(?a,f1(p`i'),?b); + #enddo + #do i = 1, `numOut' + repeat id vx(?a,f1?(-2*`i'),?b) = vx(?a,f1(q`i'),?b); + #enddo + repeat id vx(?a,f1?(-1),?b) = vx(?a,f1(p1),?b); + repeat id vx(?a,f1?(-2),?b) = vx(?a,f1(q1),?b); + multiply counter(1); + repeat id vx(?a,f1?(n1?),?b) * vx(?c,f2?(n1?),?d) * counter(n2?) + = vx(?a,f1(kk[n2]),?b) * vx(?c,f2(-kk[n2]),?d) * counter(n2+1); + repeat id vx(?a,f1?(n1?),?b,f2?(n1?),?c) * counter(n2?) + = vx(?a,f1(kk[n2]),?b,f2(-kk[n2]),?c) * counter(n2+1); + multiply replace_(counter,dummy_); + multiply replace_(vx,node); + +#ifdef `Verbose' + P +sss; +#endif + .sort:mkdiag; +#endprocedure + +*--#] MakeDiagrams : +*--#[ CanonicalizeDiagrams : + +** +* Canonicalizes diagrams. +* +* This procedure attaches a unique topology ID and reassigns +* the loop momenta k1, k2, ... in a canonical way. +* +* Example Input +* ------------- +* node(QUA(q1),qua(-k2),glu(-k1)) +* *node(QUA(k2),qua(p1),glu(k1)) +* +* Example Output +* -------------- +* topo(1) +* *node(QUA(q1),qua(-k1),glu(k2)) +* *node(QUA(k1),qua(p1),glu(-k2)) +* +#procedure CanonicalizeDiagrams() + #define numactiveexprs "`numactiveexprs_'" + #define activeexprnames "`activeexprnames_'" + +* Label each term uniquely. + #$counter = 0; + #define i "0" + #do e = {`activeexprnames',} + #ifdef `e' + #redefine i "{`i'+1}" + #$counter`i' = 0; + inexpression `e'; + $counter`i' = $counter`i' + 1; + $counter = $counter`i'; + multiply counter($counter`i'); + goto end; + endinexpression; + #endif + #enddo + #undefine i +label end; + B node; + ModuleOption maximum, $counter; + #do i = 1, `numactiveexprs' + ModuleOption local, $counter`i'; + #enddo + Moduleoption inparallel; + .sort:canondiag-1; + +* Extract vertex structure into acc(vxs(v1,...)*...). + Keep Brackets; + id node(?a) = node(?a) * vxs(?a); + argument vxs; + id f1?(p1?) = p1; + endargument; + multiply acc(1); + multiply replace_(vxs,vx); + id vx(?a) = acc(vxs(?a)); + repeat id acc(x1?) * acc(x2?) = acc(x1 * x2); + B acc; + .sort:canondiag-2; + + #do e = {`activeexprnames',} + #ifdef `e' + CreateSpectator [`e'.spec], "`e'.`PID_'.spec"; + #endif + #enddo + + #$numtopo = 0; + + #do loop = 1, 1 +* Find a new topology. + Keep Brackets; + #$found = 0; + if ($found == 0); + if (match(acc(x1?$x1))); + $found = 1; + $numtopo = $numtopo + 1; + $topo = $x1; + inside $x1; + $numv = count_(vxs,1); + endinside; + endif; + endif; + B acc; + ModuleOption noparallel; + .sort:canondiag-3-{`$numtopo'+1}; + #if `$found' +* Try to match the topology with the current term in all possible ways. + Keep Brackets; + #$found = 0; + $replaced = 0; + argument acc; +#UseDictionary wildmom($) + if ((count(vxs,1) == `$numv') && match(`$topo')); +#CloseDictionary + $found = 1; + $replaced = 1; + $t = term_; + inside $t; + id all, $topo * replace_(,...,) + = $topo * *...*; + multiply replace_(vxs,dummy_); + id replace(-k1?vector_,k2?) = replace(k1,-k2); + chainin replace; + multiply replace_(,...,); + endinside; + multiply replace_(vxs,dummy_); + multiply topo(`$numtopo'); + multiply $t; + endif; + endargument; + if ($replaced); + id acc(x1?) = x1; + endif; + ModuleOption maximum, $found; + ModuleOption local, $t, $replaced, <$k1>,...,<$k`MaxPropagators'>; + .sort:canondiag-4-`$numtopo'; + #if `$found' +* Perform all the generated replacements. + #redefine loop "0" + id ifnomatch-> notreplaced, replace(?a) = replace_(?a); + #do e = {`activeexprnames',} + #ifdef `e' + inexpression `e'; + tospectator [`e'.spec]; + endinexpression; + #endif + #enddo +label notreplaced; + B acc; + .sort:canondiag-5-`$numtopo'; + #endif + #endif + #enddo + +* Restore expressions from the spectators. + #do e = {`activeexprnames',} + #ifdef `e' + CopySpectator [`e'.copy] = [`e'.spec]; + #exchange `e',[`e'.copy] + #endif + #enddo + .sort:canondiag-6; + #do e = {`activeexprnames',} + #ifdef `e' + Drop [`e'.copy]; + EmptySpectator [`e'.spec]; + #endif + #enddo + +* Drop duplicate terms (keep the first one). + #define maxunrolling "100" + #define unrolling + #define i "0" + #do k = 1, `$counter' + #if {`k'+`maxunrolling'-1} <= `$counter' + #redefine unrolling "`maxunrolling'" + #else + #redefine unrolling "{`$counter'-`k'+1}" + #endif + #do j = 0, {`unrolling'-1} + if (match(counter({`k'+`j'}))); + #redefine i "0" + #do e = {`activeexprnames',} + #ifdef `e' + #redefine i "{`i'+1}" + #$found`i'x`j' = 0; + inexpression `e'; + if ($found`i'x`j'); + discard; + else; + $found`i'x`j' = 1; + multiply replace_(counter,dummy_); + tospectator [`e'.spec]; + endif; + endinexpression; + #endif + #enddo + endif; + #enddo + #do i = 1, `numactiveexprs' + #do j = 0, {`unrolling'-1} + ModuleOption local, $found`i'x`j'; + #enddo + #enddo + ModuleOption inparallel; + .sort:canondiag-7-`k'; + #redefine k "{`k'+`unrolling'-1}" + #enddo + #undefine i + #if `$counter' == 0 + .sort:canondiag-7; + #endif + +* Restore expressions from the spectators. + #do e = {`activeexprnames',} + #ifdef `e' + CopySpectator [`e'.copy] = [`e'.spec]; + #exchange `e',[`e'.copy] + #endif + #enddo + .sort:canondiag-8; + #do e = {`activeexprnames',} + #ifdef `e' + Drop [`e'.copy]; + RemoveSpectator [`e'.spec]; + #endif + #enddo + +* Reorder identical boson fields at each vertex to align the momenta. + #do b = {`Bosons',} + #ifdef `b' +* Assume that they appear in adjacent positions. + repeat id node(?a,`b'(?b),`b'(?c),?d) = node(?a,`b'(?b,?c),?d); + multiply replace_(`b',vxs); + multiply replace_(vxs,`b'); + repeat id node(?a,`b'(p1?,p2?,?b),?c) = node(?a,`b'(p1),`b'(p2,?b),?c); + #endif + #enddo + +* Align the loop momenta with the flow of particles. + #ifdef `AntiParticles' + #do i = 1, `MaxPropagators' + if (match(node(?a,f1?{`AntiParticles'}(-k`i'),?b))); + multiply replace_(k`i',-k`i'); + endif; + #enddo + #endif + + .sort:canondiag-9; +#endprocedure + +*--#] CanonicalizeDiagrams : +*--#[ CompareDiagrams : + +** +* Compares diagrams in two expressions. +* +* The common diagrams are extracted from F1 and F2, and stored in F0. +* Diagrams must consist of "topo" and "node" functions. +* +* Parameters +* ---------- +* F1, F2: expression [in/out] +* The expressions to be compared. Common diagrams will be removed from +* these expressions. +* F0: expression [out] +* The extracted common diagrams. +* +#procedure CompareDiagrams(F1,F2,F0) + #call CanonicalizeDiagrams() + + B+ topo,node; + .sort:compdiag-1; + Keep Brackets; + L `F0' = 0 + #do t = `F1'; + #$t = `t'; + #inside $t; + dropcoefficient; + #endinside; + #$t1 = `F1'[`$t']; + #$t2 = `F2'[`$t']; + #$tt = $t1 - $t2; + #if termsin($tt) == 0 + + (`t') + #endif + #enddo + ; + .sort + L `F1' = `F1' - `F0'; + L `F2' = `F2' - `F0'; +#ifdef `Verbose' + P +sss; +#endif + .sort:compdiag-3; +#endprocedure + +*--#] CompareDiagrams : +*--#[ DoComparison : + +** +* Performs comparison. +* +* Example +* ------- +* #call DoComparison(qcd,in=qua,out=qua,loops=1,form_options=`NOSNAIL_',qgraf_options=nosnail) +* +#procedure DoComparison(model,?a) + #switch `DiagramGenerator' + #case qgraf + #define DiagramGeneratorName "QGRAF" + #break + #case feyngraph + #define DiagramGeneratorName "FeynGraph" + #break + #default + #ifndef `DiagramGenerator' + #message DiagramGenerator not defined + #else + #message Unknown DiagramGenerator: `DiagramGenerator' + #endif + #terminate + #endswitch + + #define in "" + #define out "" + #define loops "" + #define formoptions "" + #define qgrafoptions "" + #define feyngraphoptions "" + + #define next "" + #do a = {`?a',} + #ifdef `a' + #if "`keepleft_(`a',3)'" == "in=" + #redefine next "in" + #redefine a "`takeleft_(`a',3)'" + #elseif "`keepleft_(`a',4)'" == "out=" + #redefine next "out" + #redefine a "`takeleft_(`a',4)'" + #elseif "`keepleft_(`a',6)'" == "loops=" + #redefine next "loops" + #redefine a "`takeleft_(`a',6)'" + #elseif "`keepleft_(`a',13)'" == "form_options=" + #redefine next "formoptions" + #redefine a "`takeleft_(`a',13)'" + #elseif "`keepleft_(`a',14)'" == "qgraf_options=" + #redefine next "qgrafoptions" + #redefine a "`takeleft_(`a',14)'" + #elseif "`keepleft_(`a',18)'" == "feyngraph_options=" + #redefine next "feyngraphoptions" + #redefine a "`takeleft_(`a',18)'" + #endif + #ifdef `a' + #ifdef `next' + #ifdef ``next'' + #redefine `next' "``next'',`a'" + #else + #redefine `next' "`a'" + #endif + #endif + #endif + #endif + #enddo + + #ifdef `formoptions' + #redefine formoptions "{`formoptions'}" + #else + #redefine formoptions "0" + #endif + + #define numIn "0" + #define numOut "0" + + #do a = {`in',} + #ifdef `a' + #redefine numIn "{`numIn'+1}" + #endif + #enddo + + #do a = {`out',} + #ifdef `a' + #redefine numOut "{`numOut'+1}" + #endif + #enddo + + #define momenta "" + + #do i = 1, `numIn' + #redefine momenta "`momenta',p`i'" + #enddo + + #do i = 1, `numOut' + #redefine momenta "`momenta',q`i'" + #enddo + + #define cachekey "`DiagramGenerator'-`model'" + #do b = {in,out,loops,formoptions,`DiagramGenerator'options} + #redefine cachekey "`cachekey'-" + #ifdef ``b'' + #do a = {``b'',} + #ifdef `a' + #redefine cachekey "`cachekey'`a'" + #endif + #enddo + #else + #redefine cachekey "`cachekey'none" + #endif + #enddo + + #define cachefile "`CacheDir'/`cachekey'.dat" + #define cachetmpfile "`cachekey'.dat.`PID_'.tmp" + #define cachetmpfile2 "`CacheDir'/`cachekey'.dat2.`PID_'.tmp" + +* Set up the model. + + #call Setup`toupper_(`keepleft_(`model',1)')'`takeleft_(`model',1)'Model() + +* Generate diagrams. + + L F1 = diagrams_(`CurrentModel',{`in'},{`out'},{`momenta'},kk,`loops',`formoptions'); +#ifdef `Verbose' + P +sss; +#endif + .sort:gen-diag; + #$nF1 = termsin_(F1); + +* Check if the cache exists. + + #ifdef `CacheRead' + #$cacheexists = + #pipe if [ -f "`cachefile'" ]; then echo 1; else echo 0; fi + ; + #else + #$cacheexists = 0; + #endif + + #if `$cacheexists' == 1 +* If the cache exists, then compare the generated diagrams with the cache, +* which has already been confirmed to be consistent with the other generator. + + #message Read cache file: `cachefile' + + Skip; + + L F2 = + #include `cachefile' + ; + .sort:load-cache; + Drop; + L oldF1 = F1; + L oldF2 = F2; + .sort; + Drop; + L F0 = 0; + L F1 = oldF1 - oldF2; + L F2 = - F1; + #else +* Otherwise, we must perform a full comparison. + + #ifdef `CacheWrite' +* Output the generated diagrams into a temporary file. + .sort:save-tmp; + Format 255; + #write <`cachetmpfile'> " %+E", F1 + Format 72; + #close <`cachetmpfile'> + #endif + + Skip; + + L F2 = + #switch `DiagramGenerator' + #case qgraf + #call RunQgraf(in=`in',out=`out',loops=`loops',options=`qgrafoptions') + #break + #case feyngraph + #call RunFeynGraph(in=`in',out=`out',loops=`loops',options=`feyngraphoptions') + #break + #endswitch + ; +#ifdef `Verbose' + P +sss; +#endif + .sort:run-`DiagramGenerator'; + #$nF2 = termsin_(F2); + + #call MakeDiagrams(`numIn',`numOut') + #call CompareDiagrams(F1,F2,F0) + + #ifdef `CacheWrite' +* If the generated diagrams are consistent with those by the other generator +* then write the cache file. + #if (`ZERO_F1' == 1) && (`ZERO_F2' == 1) + #message Write cache file: `cachefile' + #system mkdir -p "`CacheDir'" + #system echo "* Compared with `DiagramGeneratorVersion'" >"`cachetmpfile2'" + #system echo "* Number of diagrams: FORM: `$nF1', `DiagramGeneratorName': `$nF2'" >>"`cachetmpfile2'" + #system sed "s/^ */ /" "`cachetmpfile'" >>"`cachetmpfile2'" + #system mv "`cachetmpfile2'" "`cachefile'" + #endif + #remove <`cachetmpfile'> + #endif + #endif +#endprocedure + +*--#] DoComparison : +*--#] diagram_compare_include : +*--#[ qgraf_compare_include : +#- +#include compare.frm # diagram_compare_include + +*--#[ RunQgraf : + +** +* Runs QGRAF and reads its output. +* +* Example +* ------- +* #call SetupQcdModel() +* L F = +* #call RunQgraf(in=qua,out=qua,loops=1,options=nosnail) +* ; +* +* Example Output +* -------------- +* F = +* vx(QUA(-2),qua(1),glu(2))*vx(QUA(1),qua(-1),glu(2)); +* +#procedure RunQgraf(?a) + #define in "" + #define out "" + #define loops "" + #define loopmomentum "" + #define options "" + #define next "options" + + #do a = {`?a',} + #ifdef `a' + #if "`keepleft_(`a',3)'" == "in=" + #redefine next "in" + #redefine a "`takeleft_(`a',3)'" + #elseif "`keepleft_(`a',4)'" == "out=" + #redefine next "out" + #redefine a "`takeleft_(`a',4)'" + #elseif "`keepleft_(`a',6)'" == "loops=" + #redefine next "loops" + #redefine a "`takeleft_(`a',6)'" + #elseif "`keepleft_(`a',14)'" == "loop_momentum=" + #redefine next "loopmomentum" + #redefine a "`takeleft_(`a',14)'" + #elseif "`keepleft_(`a',8)'" == "options=" + #redefine next "options" + #redefine a "`takeleft_(`a',8)'" + #endif + #ifdef `a' + #ifdef ``next'' + #redefine `next' "``next'',`a'" + #else + #redefine `next' "`a'" + #endif + #endif + #endif + #enddo + + #ifndef `loopmomentum' + #redefine loopmomentum "k" + #endif + + #remove <`QgrafOutputFile'> + + #define filename "`QgrafInputFile'" + #create <`filename'> + #write <`filename'> "output = \'`QgrafOutputFile'\';" + #write <`filename'> "style = \'`QgrafStyleFile'\';" + #write <`filename'> "model = \'`QgrafModelFile'\';" + #write <`filename'> "in = `in';" + #write <`filename'> "out = `out';" + #write <`filename'> "loops = `loops';" + #write <`filename'> "loop_momentum = `loopmomentum';" + #write <`filename'> "options = `options';" + #close <`filename'> + + #system qgraf + +#ifdef `Verbose' + #system cat `QgrafOutputFile' +#endif + #include `QgrafOutputFile' +#endprocedure + +*--#] RunQgraf : +*--#[ SetupQgrafFormStyle : + +** +* Prepares the QGRAF style file for FORM. +* +#procedure SetupQgrafFormStyle() + #define filename "`QgrafStyleFile'" + #create <`filename'> + #write <`filename'> "" + #write <`filename'> "*" + #write <`filename'> "* file generated by " + #write <`filename'> "*" + #write <`filename'> "* *" + #write <`filename'> "*" + #write <`filename'> "#redefine QgrafVersion \"\"" + #write <`filename'> "*" + #write <`filename'> "" + #write <`filename'> "*--#[[ d:" + #write <`filename'> " " + #write <`filename'> " *vx((),)" + #write <`filename'> "" + #write <`filename'> "*--#]] d: " + #write <`filename'> "" + #write <`filename'> "*" + #write <`filename'> "* end" + #write <`filename'> "*" + #write <`filename'> "" + #close <`filename'> +#endprocedure + +*--#] SetupQgrafFormStyle : + +#call SetupQgrafFormStyle() +#redefine DiagramGenerator "qgraf" +#redefine DiagramGeneratorVersion "`~QgrafVersion'" + +*--#] qgraf_compare_include : +*--#[ feyngraph_compare_include : +#- +#include compare.frm # diagram_compare_include + +*--#[ RunFeynGraph : + +** +* Runs FeynGraph and reads its output. +* +* Example +* ------- +* #call SetupQcdModel() +* L F = +* #call RunFeynGraph(in=qua,out=qua,loops=1,options=self_loops=0) +* ; +* +* Example Output +* -------------- +* F = +* vx(QUA(-2),qua(3),glu(2))*vx(QUA(3),qua(-1),glu(2)); +* +#procedure RunFeynGraph(?a) + #define in "" + #define out "" + #define loops "" + #define loopmomentum "" + #define options "" + #define next "options" + + #do a = {`?a',} + #ifdef `a' + #if "`keepleft_(`a',3)'" == "in=" + #redefine next "in" + #redefine a "`takeleft_(`a',3)'" + #elseif "`keepleft_(`a',4)'" == "out=" + #redefine next "out" + #redefine a "`takeleft_(`a',4)'" + #elseif "`keepleft_(`a',6)'" == "loops=" + #redefine next "loops" + #redefine a "`takeleft_(`a',6)'" + #elseif "`keepleft_(`a',14)'" == "loop_momentum=" + #redefine next "loopmomentum" + #redefine a "`takeleft_(`a',14)'" + #elseif "`keepleft_(`a',8)'" == "options=" + #redefine next "options" + #redefine a "`takeleft_(`a',8)'" + #endif + #ifdef `a' + #ifdef ``next'' + #redefine `next' "``next'',`a'" + #else + #redefine `next' "`a'" + #endif + #endif + #endif + #enddo + + #ifndef `loopmomentum' + #redefine loopmomentum "k" + #endif + + #remove <`FeynGraphOutputFile'> + + #define filename "`FeynGraphInputFile'" + #create <`filename'> + #write <`filename'> "template = \"`FeynGraphTemplateFile'\"" + #write <`filename'> "[process]" + #write <`filename'> "in = [%" + #do a = {`in',} + #ifdef `a' + #write <`filename'> "\"`a'\", %" + #endif + #enddo + #write <`filename'> "]" + #write <`filename'> "out = [%" + #do a = {`out',} + #ifdef `a' + #write <`filename'> "\"`a'\", %" + #endif + #enddo + #write <`filename'> "]" + #write <`filename'> "loops = `loops'" + #write <`filename'> "model = \"`QgrafModelFile'\"" + #write <`filename'> "momenta = [%" + #define i "1" + #do a = {`in',} + #ifdef `a' + #write <`filename'> "\"p`i++'\", %" + #endif + #enddo + #redefine i "1" + #do a = {`out',} + #ifdef `a' + #write <`filename'> "\"q`i++'\", %" + #endif + #enddo + #do i = 1, `loops' + #write <`filename'> "\"k`i'\", %" + #enddo + #write <`filename'> "]" + #write <`filename'> "[filter]" + #do a = {`options',} + #ifdef `a' + #switch `a' + #case onshell + #write <`filename'> "`a' = true" + #break + #default + #write <`filename'> "`a'" + #endswitch + #endif + #enddo + #close <`filename'> + + #system feyngraph `FeynGraphInputFile' + #pipe echo \#redefine FeynGraphVersion \"$(feyngraph --version)\" + +#ifdef `Verbose' + #system cat `FeynGraphOutputFile' +#endif + #include `FeynGraphOutputFile' +#endprocedure + +*--#] RunFeynGraph : +*--#[ SetupFeynGraphTemplate : + +** +* Prepares the FeynGraph template file for FORM. +* +#procedure SetupFeynGraphTemplate() + #define filename "`FeynGraphTemplateFile'" + #create <`filename'> + #write <`filename'> "{\%- for diag in diagrams -\%}" + #write <`filename'> " {\%- set n_in = diag.n_in() -\%}" + #write <`filename'> " {\%- set n_ext = diag.n_ext() -\%}" + #write <`filename'> "*--#[[ d{{ loop.index }}:" + #write <`filename'> " {{ \"+\" ~ diag.sign() if diag.sign() >= 0 else diag.sign() -}}" + #write <`filename'> " {{ \"/\" ~ diag.symmetry_factor() if diag.symmetry_factor() != 1 else \"\" }}" + #write <`filename'> "{\%- for vertex in diag.vertices() \%}" + #write <`filename'> " *vx(" + #write <`filename'> "{\%- for par in vertex.particles_ordered() -\%}" + #write <`filename'> " {\%- set j = vertex.propagators_ordered()[loop.index0].id() -\%}" + #write <`filename'> " {\%- if j < n_in -\%}" + #write <`filename'> " {\%- set j = - 2 * j - 1 -\%}" + #write <`filename'> " {\%- elif j < n_ext -\%}" + #write <`filename'> " {\%- set j = - 2 * (j - n_in + 1) -\%}" + #write <`filename'> " {\%- endif -\%}" + #write <`filename'> " {{ par.name() }}({{ j }})" + #write <`filename'> " {\%- if not loop.last -\%},{\%- endif -\%}" + #write <`filename'> "{\%- endfor -\%}" + #write <`filename'> ")" + #write <`filename'> "{\%- endfor \%}" + #write <`filename'> "*--#]] d{{ loop.index }}: " + #write <`filename'> "{\% endfor \%}" + #close <`filename'> +#endprocedure + +*--#] SetupFeynGraphTemplate : + +#call SetupFeynGraphTemplate() +#redefine DiagramGenerator "feyngraph" +#redefine DiagramGeneratorVersion "feyngraph `~FeynGraphVersion'" + +*--#] feyngraph_compare_include : +*--#[ qgraf_phi4_phiphi_phiphi_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(phi4,in=phi,phi,out=phi,phi,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 7; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_phi4_phiphi_phiphi_1 : +*--#[ qgraf_phi4_phiphi_phiphi_2 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(phi4,in=phi,phi,out=phi,phi,loops=2,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 42; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_phi4_phiphi_phiphi_2 : +*--#[ qgraf_qcd_qua_qua_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(qcd,in=qua,out=qua,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 4; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_qcd_qua_qua_1 : +*--#[ qgraf_qcd_glu_glu_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(qcd,in=glu,out=glu,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 7; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_qcd_glu_glu_1 : +*--#[ qgraf_qcd_gho_gho_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(qcd,in=gho,out=gho,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 4; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_qcd_gho_gho_1 : +*--#[ qgraf_qcd_qua_quaglu_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(qcd,in=qua,out=qua,glu,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 17; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_qcd_qua_quaglu_1 : +*--#[ qgraf_qcd_glu_gluglu_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(qcd,in=glu,out=glu,glu,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 32; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_qcd_glu_gluglu_1 : +*--#[ qgraf_qcd_gluglu_gluglu_1 : +#include- compare.frm # qgraf_compare_include +#call DoComparison(qcd,in=glu,glu,out=glu,glu,loops=1,form_options=,qgraf_options=) +.end +assert succeeded? +assert nterms("F1", 0) == 223; +assert nterms("F1") == 0; +assert nterms("F2") == 0; +*--#] qgraf_qcd_gluglu_gluglu_1 : diff --git a/check/examples.frm b/check/examples.frm index 4162b757..b1ac873a 100644 --- a/check/examples.frm +++ b/check/examples.frm @@ -1830,129 +1830,4 @@ Print; assert result("aPLUSbTO2") =~ expr("b^2 + 2*a*b + a^2") assert result("aPLUSbTO3") =~ expr("b^3 + 3*a*b^2 + 3*a^2*b + a^3") *--#] ExtComm_1 : -*--#[ Diagrams_1 : - Vectors Q1,Q2,p1,...,p8; - Set QQ:Q1,Q2; - Set PP:p1,...,p8; - #define LOOPS "2" - Local F = topologies_(`LOOPS',2,{3,},QQ,PP); - Print +f +s; - .end -# TODO: enable it -#pend_if true - assert succeeded? - assert result("F") =~ expr(" - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,p1,-p3)* - node_(4,p2,-p4,-p5)*node_(5,p3,p4,p5) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,p1,-p3,-p4)* - node_(4,p2,p3,-p5)*node_(5,Q2,p4,p5) - ") -*--#] Diagrams_1 : -*--#[ Diagrams_2 : - Vectors Q1,Q2,p1,...,p8; - Set QQ:Q1,Q2; - Set PP:p1,...,p8; - #define LOOPS "2" - Local F = topologies_(`LOOPS',2,{3,4},QQ,PP); - Print +f +s; - .end -# TODO: enable it -#pend_if true - assert succeeded? - assert result("F") =~ expr(" - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,Q2,-p1,-p2)*node_(3,p1,p2,-p3,p3 - ) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,p1,-p3)* - node_(4,p2,p3,-p4,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,p1,-p3)* - node_(4,p2,-p4,-p5)*node_(5,p3,p4,p5) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,-p3,-p4)* - node_(4,p1,p2,p3,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,p1,-p3,-p4)* - node_(4,Q2,p2,p3,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,p1,-p3,-p4)* - node_(4,p2,p3,-p5)*node_(5,Q2,p4,p5) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2,-p3)*node_(3,Q2,p1,p2,p3 - ) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,-p1,-p2,-p3)*node_(3,Q2,p1,-p4)* - node_(4,Q1,p2,p3,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,-p1,-p2,-p3)*node_(3,p1,p2,-p4)* - node_(4,Q1,Q2,p3,p4) - ") -*--#] Diagrams_2 : -*--#[ Diagrams_3 : - Vectors Q1,Q2,p1,...,p8; - Set QQ:Q1,Q2; - Set PP:p1,...,p8; - #define LOOPS "2" - Local F = topologies_(`LOOPS',-2,{3,4},QQ,PP); - Print +f +s; - .end -# TODO: enable it -#pend_if true - assert succeeded? - assert result("F") =~ expr(" - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,Q2,-p1,-p2)*node_(3,p1,p2,-p3,p3 - ) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,p1,-p3)* - node_(4,p2,p3,-p4,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,p1,-p3)* - node_(4,p2,-p4,-p5)*node_(5,p3,p4,p5) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,-p3,-p4)* - node_(4,p1,p2,p3,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,Q2,-p3,-p4)* - node_(4,p1,p3,-p5)*node_(5,p2,p4,p5) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2)*node_(3,p1,-p3,-p4)* - node_(4,Q2,p2,p3,p4) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,Q1,-p1,-p2,-p3)*node_(3,Q2,p1,p2,p3 - ) - + node_(0,-Q1)*node_(1,-Q2)*node_(2,-p1,-p2,-p3)*node_(3,p1,p2,-p4)* - node_(4,Q1,Q2,p3,p4) - ") -*--#] Diagrams_3 : -*--#[ Diagrams_4 : - Vectors Q1,Q2,p1,...,p17; - Set QQ:Q1,Q2; - Set PP:p1,...,p17; - #define LOOPS "6" - Local F = topologies_(`LOOPS',-2,{3,},QQ,PP); - .end -# TODO: enable it -#pend_if true - #time_dilation 2.0 - assert succeeded? - assert nterms("F") == 2793 -*--#] Diagrams_4 : -*--#[ Diagrams_5 : - #define LOOPS "6" - Model PHI3; - Particle phi,+1; - Vertex phi,phi,phi:g; - EndModel; - Vector Q1,Q2,p1,...,p{3*`LOOPS'}; - Set QQ:Q1,Q2; - Set pp:p1,...,p{3*`LOOPS'}; - Set empty:; - Local F = diagrams_(PHI3,{phi,phi},empty,QQ,pp,`LOOPS', - `OnePI_'+`NoTadpoles_'+`Symmetrize_'+`TopologiesOnly_'); - .end - assert succeeded? - assert nterms("F") == 2793 -*--#] Diagrams_5 : -*--#[ Diagrams_6 : -* #define LOOPS "6" -* Model PHI3; -* Particle phi,+1; -* Vertex phi,phi,phi:g; -* EndModel; -* Vector Q1,Q2,p1,...,p{3*`LOOPS'}; -* Set QQ:Q1,Q2; -* Set pp:p1,...,p{3*`LOOPS'}; -* Set empty:; -* Local F = diagrams_(PHI3,{phi,phi},empty,QQ,pp,`LOOPS', -* `OnePI_'+`NoTadpoles_'+`TopologiesOnly_'); -* .end -* assert succeeded? -* assert nterms("F") == 4999 -*--#] Diagrams_6 : diff --git a/sources/diawrap.cc b/sources/diawrap.cc index 723acc9c..e9b0800f 100644 --- a/sources/diawrap.cc +++ b/sources/diawrap.cc @@ -14,7 +14,7 @@ typedef struct ToPoTyPe { WORD *vertmax; Options *opt; int cldeg[MAXPOINTS], clnum[MAXPOINTS], clext[MAXPOINTS]; - int cmind[MAXPOINTS], cmaxd[MAXPOINTS]; + int cmind[MAXLEGS+1],cmaxd[MAXLEGS+1]; int ncl, nvert; } TOPOTYPE; @@ -199,6 +199,21 @@ int ReConvertParticle(Model *model,int grccnum) } // #] ReConvertParticle : +// #[ numParticle : + +int numParticle(MODEL *m,WORD n) +{ + int i; + for ( i = 0; i < m->nparticles; i++ ) { + if ( m->vertices[i]->particles[0].number == n ) return(i); + if ( m->vertices[i]->particles[1].number == n ) return(i); + } + MesPrint("numParticle: particle %d not found in model",n); + Terminate(-1); + return(-1); +} + +// #] numParticle : // #[ ProcessDiagram : void ProcessDiagram(EGraph *eg, void *ti) @@ -230,7 +245,8 @@ void ProcessDiagram(EGraph *eg, void *ti) // // Now get the nodes // - for ( i = 0; i < eg->nNodes; i++ ) { + if ( ( info->flags & WITHOUTNODES ) == 0 ) { + for ( i = 0; i < eg->nNodes; i++ ) { // // node_(number,coupling,particle_1(momentum_1),...,particle_n(momentum_n)) // @@ -245,7 +261,7 @@ void ProcessDiagram(EGraph *eg, void *ti) // function for when we want to work with counterterms. // if ( !eg->isExternal(i) ) { - afill = fill; *fill++ = 0; *fill++ = 0; FILLARG(fill) + afill = fill; *fill++ = 0; *fill++ = 1; FILLARG(fill) cfill = fill; *fill++ = 0; intr = eg->nodes[i]->intrct; for ( j = 0; j < model->interacts[intr]->nclist; j++ ) { @@ -273,8 +289,10 @@ void ProcessDiagram(EGraph *eg, void *ti) edge = eg->nodes[i]->edges[j]; vect = ABS(edge)-1; *fill++ = 6+FUNHEAD; - - *fill++ = ReConvertParticle(model,eg->edges[vect]->ptcl); // code of particle_j + int a; + if ( edge < 0 ) { a = ReConvertParticle(model,-eg->edges[vect]->ptcl); } + else { a = ReConvertParticle(model,eg->edges[vect]->ptcl); } + *fill++ = a; *fill++ = FUNHEAD+2; FILLFUN(fill) *fill++ = edge < 0 ? -MINVECTOR: -VECTOR; @@ -287,26 +305,131 @@ void ProcessDiagram(EGraph *eg, void *ti) *fill++ = 1; *fill++ = 1; *fill++ = 3; } startfill[1] = fill-startfill; + } + } + if ( ( info->flags & WITHEDGES ) == WITHEDGES ) { + for ( i = 0; i < eg->nEdges; i++ ) { + int n1 = eg->edges[i]->nodes[0]; + int n2 = eg->edges[i]->nodes[1]; +// int l1 = eg->edges[i]->nlegs[0]; +// int l2 = eg->edges[i]->nlegs[1]; + + startfill = fill; + *fill++ = EDGE; + *fill++ = 0; + FILLFUN(fill) + *fill++ = -SNUMBER; *fill++ = i+1; // number of the edge +// + *fill++ = ARGHEAD+FUNHEAD+6; + *fill++ = 0; + FILLARG(fill) + *fill++ = 6+FUNHEAD; + int a = ReConvertParticle(model,eg->edges[i]->ptcl); + *fill++ = a; + *fill++ = FUNHEAD+2; + FILLFUN(fill) + *fill++ = -VECTOR; + // Look up in set of internal momenta set + if ( i < info->numextern ) { // Look up in set of external momenta + *fill++ = SetElements[Sets[info->externalset].first+i]; + } + else { // Look up in set of internal momenta set + *fill++ = SetElements[Sets[info->internalset].first+(i-eg->nExtern)]; + } + *fill++ = 1; *fill++ = 1; *fill++ = 3; +// + *fill++ = -SNUMBER; *fill++ = n1+1; // number of the node from + *fill++ = -SNUMBER; *fill++ = n2+1; // number of the node to + startfill[1] = fill - startfill; + } + } + if ( ( info->flags & WITHBLOCKS ) == WITHBLOCKS ) { + for ( i = 0; i < eg->econn->nblocks; i++ ) { + startfill = fill; + *fill++ = BLOCK; + *fill++ = 0; + FILLFUN(fill); + *fill++ = -SNUMBER; + *fill++ = i+1; + *fill++ = -SNUMBER; + *fill++ = eg->econn->blocks[i].loop; +// +// Now we have to make a list of all nodes inside this block +// + int bnodes[GRCC_MAXNODES], k; + WORD *argfill = fill, *funfill; + *fill++ = 0; *fill++ = 0; FILLARG(fill) + *fill++ = 0; + for ( k = 0; k < GRCC_MAXNODES; k++ ) bnodes[k] = 0; + for ( k = 0; k < eg->econn->blocks[i].nmedges; k++ ) { + bnodes[eg->econn->blocks[i].edges[k][0]] = 1; + bnodes[eg->econn->blocks[i].edges[k][1]] = 1; + } + for ( k = 0; k < GRCC_MAXNODES; k++ ) { + if ( bnodes[k] == 0 ) continue; +// +// Now we put the node inside this argument. +// + funfill = fill; + *fill++ = NODEFUNCTION; + *fill++ = 0; + FILLFUN(fill) + *fill++ = -SNUMBER; *fill++ = k+1; + numlegs = eg->nodes[k]->deg; + for ( j = 0; j < numlegs; j++ ) { + edge = eg->nodes[k]->edges[j]; + vect = ABS(edge)-1; + *fill++ = -VECTOR; + if ( numlegs == 1 || vect < info->numextern ) { // Look up in set of external momenta + *fill++ = SetElements[Sets[info->externalset].first+vect]; + } + else { // Look up in set of internal momenta set + *fill++ = SetElements[Sets[info->internalset].first+(vect-info->numextern)]; + } + } + funfill[1] = fill-funfill; + } + *fill++ = 1; *fill++ = 1; *fill++ = 3; + *argfill = fill - argfill; + argfill[ARGHEAD] = argfill[0] - ARGHEAD; + startfill[1] = fill-startfill; + } + } + if ( ( info->flags & WITHONEPISETS ) == WITHONEPISETS ) { + for ( i = 0; i < eg->econn->nopic; i++ ) { + startfill = fill; + *fill++ = ONEPI; + *fill++ = 0; + FILLFUN(fill); + *fill++ = -SNUMBER; + *fill++ = i+1; + *fill++ = -ONEPI; + for ( j = 0; j < eg->econn->opics[i].nnodes; j++ ) { + *fill++ = -SNUMBER; + *fill++ = eg->econn->opics[i].nodes[j]+1; + } + startfill[1] = fill-startfill; + } } // // Topology counter. We have exaggerated a bit with the eye on the far future. // - if ( info->numtopo < MAXPOSITIVE ) { + if ( info->numtopo-1 < MAXPOSITIVE ) { *fill++ = TOPO; *fill++ = FUNHEAD+2; FILLFUN(fill) - *fill++ = -SNUMBER; *fill++ = (WORD)(info->numtopo); + *fill++ = -SNUMBER; *fill++ = (WORD)(info->numtopo-1); } - else if ( info->numtopo < FULLMAX-1 ) { + else if ( info->numtopo-1 < FULLMAX-1 ) { *fill++ = TOPO; *fill++ = FUNHEAD+ARGHEAD+4; FILLFUN(fill) *fill++ = ARGHEAD+4; *fill++ = 0; FILLARG(fill) *fill++ = 4; - *fill++ = (WORD)(info->numtopo & WORDMASK); + *fill++ = (WORD)((info->numtopo-1) & WORDMASK); *fill++ = 1; *fill++ = 3; } else { // for now: science fiction *fill++ = TOPO; *fill++ = FUNHEAD+ARGHEAD+6; FILLFUN(fill) *fill++ = ARGHEAD+6; *fill++ = 0; FILLARG(fill) - *fill++ = 6; *fill++ = (WORD)(info->numtopo >> BITSINWORD); - *fill++ = (WORD)(info->numtopo & WORDMASK); + *fill++ = 6; *fill++ = (WORD)((info->numtopo-1) >> BITSINWORD); + *fill++ = (WORD)((info->numtopo-1) & WORDMASK); *fill++ = 0; *fill++ = 1; *fill++ = 5; } // @@ -334,14 +457,51 @@ void ProcessDiagram(EGraph *eg, void *ti) } // #] ProcessDiagram : +// #[ fendMG : + +Bool fendMG(EGraph *eg, void *ti) +{ + DUMMYUSE(eg); + DUMMYUSE(ti); + return True; +} + +// #] fendMG : // #[ ProcessTopology : Bool ProcessTopology(EGraph *eg, void *ti) { // // This routine is called for each new topology. +// New convention: return True; generate the corresponding diagrams if needed +// return False; skip diagram generation (when asked for). // TERMINFO *info = (TERMINFO *)ti; + +// This seems to work properly. It was disabled before. +#define WITHEARLYVETO +#ifdef WITHEARLYVETO + if ( ( ( info->flags & CHECKEXTERN ) == CHECKEXTERN ) && info->currentMODEL != NULL ) { + int i, j; + int numlegs, vect, edge; + for ( i = 0; i < eg->nNodes; i++ ) { + if ( eg->isExternal(i) ) continue; + numlegs = eg->nodes[i]->deg; + for ( j = 0; j < numlegs; j++ ) { + edge = eg->nodes[i]->edges[j]; + vect = ABS(edge)-1; + if ( vect < info->numextern && info->legcouple[vect][numlegs] == 0 ) { + +// This cannot be. + + info->numtopo++; + return False; + } + } + } + } +#endif + if ( ( info->flags & TOPOLOGIESONLY ) == 0 ) { info->numtopo++; return True; @@ -355,6 +515,8 @@ Bool ProcessTopology(EGraph *eg, void *ti) WORD *tail = tdia + tdia[1]; WORD *tend = term + *term; WORD *fill, *startfill; + Model *model = (Model *)info->currentModel; + MODEL *m = (MODEL *)info->currentMODEL; int i, j; int numlegs, vect, edge; @@ -374,6 +536,26 @@ Bool ProcessTopology(EGraph *eg, void *ti) *fill++ = 0; FILLFUN(fill) *fill++ = -SNUMBER; *fill++ = i+1; + if ( model != NULL && m != NULL ) { + if ( !eg->isExternal(i) ) { + WORD *afill = fill; *fill++ = 0; *fill++ = 1; FILLARG(fill) + WORD *cfill = fill; *fill++ = 0; + int cpl = 2*eg->nodes[i]->extloop+eg->nodes[i]->deg-2; + *fill++ = SYMBOL; *fill++ = 4; + *fill++ = m->couplings[0]; + *fill++ = cpl; + *fill++ = 1; *fill++ = 1; *fill++ = 3; + *cfill = fill - cfill; + *afill = fill - afill; + if ( *afill == ARGHEAD+8 && afill[ARGHEAD+4] == 1 ) { + fill = afill; *fill++ = -SYMBOL; *fill++ = afill[ARGHEAD+3]; + } + } + else { + *fill++ = -SNUMBER; + *fill++ = 1; + } + } // // Now the momenta. // @@ -390,6 +572,116 @@ Bool ProcessTopology(EGraph *eg, void *ti) } startfill[1] = fill-startfill; } + if ( ( info->flags & WITHEDGES ) == WITHEDGES ) { + for ( i = 0; i < eg->nEdges; i++ ) { + int n1 = eg->edges[i]->nodes[0]; + int n2 = eg->edges[i]->nodes[1]; +// int l1 = eg->edges[i]->nlegs[0]; +// int l2 = eg->edges[i]->nlegs[1]; + startfill = fill; + *fill++ = EDGE; + *fill++ = 0; + FILLFUN(fill) + *fill++ = -SNUMBER; *fill++ = i+1; // number of the edge +// + *fill++ = -VECTOR; + if ( i < info->numextern ) { // Look up in set of external momenta + *fill++ = SetElements[Sets[info->externalset].first+i]; + } + else { // Look up in set of internal momenta set + *fill++ = SetElements[Sets[info->internalset].first+(i-eg->nExtern)]; + } +// + *fill++ = -SNUMBER; *fill++ = n1+1; // number of the node from + *fill++ = -SNUMBER; *fill++ = n2+1; // number of the node to + startfill[1] = fill - startfill; + } + } +// +// Block information +// + if ( ( info->flags & WITHBLOCKS ) == WITHBLOCKS ) { + for ( i = 0; i < eg->econn->nblocks; i++ ) { + startfill = fill; + *fill++ = BLOCK; + *fill++ = 0; + FILLFUN(fill); + *fill++ = -SNUMBER; + *fill++ = i+1; + *fill++ = -SNUMBER; + *fill++ = eg->econn->blocks[i].loop; +// +// Now we have to make a list of all nodes inside this block +// + int bnodes[GRCC_MAXNODES], k; + WORD *argfill = fill, *funfill; + *fill++ = 0; *fill++ = 0; FILLARG(fill) + *fill++ = 0; + for ( k = 0; k < GRCC_MAXNODES; k++ ) bnodes[k] = 0; + for ( k = 0; k < eg->econn->blocks[i].nmedges; k++ ) { + bnodes[eg->econn->blocks[i].edges[k][0]] = 1; + bnodes[eg->econn->blocks[i].edges[k][1]] = 1; + } + for ( k = 0; k < GRCC_MAXNODES; k++ ) { + if ( bnodes[k] == 0 ) continue; +// +// Now we put the node inside this argument. +// + funfill = fill; + *fill++ = NODEFUNCTION; + *fill++ = 0; + FILLFUN(fill) + *fill++ = -SNUMBER; *fill++ = k+1; + numlegs = eg->nodes[k]->deg; + for ( j = 0; j < numlegs; j++ ) { + edge = eg->nodes[k]->edges[j]; + vect = ABS(edge)-1; + *fill++ = -VECTOR; + if ( numlegs == 1 || vect < info->numextern ) { // Look up in set of external momenta + *fill++ = SetElements[Sets[info->externalset].first+vect]; + } + else { // Look up in set of internal momenta set + *fill++ = SetElements[Sets[info->internalset].first+(vect-info->numextern)]; + } + } + funfill[1] = fill-funfill; + } + *fill++ = 1; *fill++ = 1; *fill++ = 3; + *argfill = fill - argfill; + argfill[ARGHEAD] = argfill[0] - ARGHEAD; + startfill[1] = fill-startfill; + } + +// if ( eg->econn->narticuls > 0 ) { +// startfill = fill; +// *fill++ = BLOCK; +// *fill++ = 0; +// FILLFUN(fill); +// for ( i = 0; i < eg->econn->snodes; i++ ) { +// if ( eg->econn->articuls[i] != 0 ) { +// *fill++ = -SNUMBER; +// *fill++ = i+1; +// } +// } +// startfill[1] = fill-startfill; +// } + } + if ( ( info->flags & WITHONEPISETS ) == WITHONEPISETS ) { + for ( i = 0; i < eg->econn->nopic; i++ ) { + startfill = fill; + *fill++ = ONEPI; + *fill++ = 0; + FILLFUN(fill); + *fill++ = -SNUMBER; + *fill++ = i+1; + *fill++ = -ONEPI; + for ( j = 0; j < eg->econn->opics[i].nnodes; j++ ) { + *fill++ = -SNUMBER; + *fill++ = eg->econn->opics[i].nodes[j]+1; + } + startfill[1] = fill-startfill; + } + } // // Topology counter. We have exaggerated a bit with the eye on the far future. // @@ -431,10 +723,36 @@ Bool ProcessTopology(EGraph *eg, void *ti) Generator(BHEAD newterm,info->level); AT.WorkPointer = oldworkpointer; info->numtopo++; - return True; + return False; } -// #] ProcessTopology : +// #] ProcessTopology : +// #[ SetDualOpts : +void SetDualOpts(int *opt, const WORD num, const int key, const char* key_name, + const int dual, const char* dual_name, const int val, const int dval) { + + if ( ( num & key ) == key ) { + if ( ( num & dual ) == dual ) { + MLOCK(ErrorMessageLock); + MesPrint("&Conflicting diagram filters: %s and %s.", key_name, dual_name); + MUNLOCK(ErrorMessageLock); + Terminate(-1); + } + else { + *opt = val; + } + } + else { + if ( ( num & dual ) == dual ) { + *opt = dval; + } + else { + // The default value is always 0. + *opt = 0; + } + } +} +// #] SetDualOpts : // #[ GenDiagrams : int GenDiagrams(PHEAD WORD *term, WORD level) @@ -447,7 +765,7 @@ int GenDiagrams(PHEAD WORD *term, WORD level) int babble = 0; // Later we may set this at the FORM code level TERMINFO info; WORD inset,outset,*coupl,setnum,optionnumber = 0; - int i, cpl[GRCC_MAXNCPLG]; + int i, j, cpl[GRCC_MAXNCPLG]; int ninitl, initlPart[GRCC_MAXLEGS], nfinal, finalPart[GRCC_MAXLEGS]; for ( i = 0; i < GRCC_MAXNCPLG; i++ ) cpl[i] = 0; // @@ -459,6 +777,7 @@ int GenDiagrams(PHEAD WORD *term, WORD level) info.diaoffset = AR.funoffset; info.externalset = term[info.diaoffset+FUNHEAD+7]; info.internalset = term[info.diaoffset+FUNHEAD+9]; + info.flags = 0; inset = term[info.diaoffset+FUNHEAD+3]; outset = term[info.diaoffset+FUNHEAD+5]; coupl = term + info.diaoffset + FUNHEAD + 10; @@ -471,7 +790,6 @@ int GenDiagrams(PHEAD WORD *term, WORD level) if ( term[info.diaoffset+1] > *coupl+FUNHEAD+10 ) optionnumber = term[info.diaoffset+*coupl+FUNHEAD+11]; } - setnum = term[info.diaoffset+FUNHEAD+1]; m = AC.models[SetElements[Sets[setnum].first]]; @@ -486,33 +804,57 @@ int GenDiagrams(PHEAD WORD *term, WORD level) opt = new Options(); - opt->setOutAG(ProcessDiagram, (void*) &info); - opt->setOutMG(ProcessTopology, (void*) &info); + opt->setOutAG(ProcessDiagram, &info); + opt->setOutMG(ProcessTopology, &info); + + opt->values[GRCC_OPT_SymmInitial] = ( optionnumber & WITHSYMMETRIZEI ) == WITHSYMMETRIZEI; + opt->values[GRCC_OPT_SymmFinal] = ( optionnumber & WITHSYMMETRIZEF ) == WITHSYMMETRIZEF; - opt->values[GRCC_OPT_1PI] = ( optionnumber & ONEPARTICLEIRREDUCIBLE ) == ONEPARTICLEIRREDUCIBLE; - opt->values[GRCC_OPT_NoTadpole] = ( optionnumber & NOTADPOLES ) == NOTADPOLES; - opt->values[GRCC_OPT_No1PtBlock] = ( optionnumber & NOTADPOLES ) == NOTADPOLES; - opt->values[GRCC_OPT_No2PtL1PI] = ( optionnumber & WITHINSERTIONS ) == WITHINSERTIONS; - opt->values[GRCC_OPT_SymmInitial] = ( optionnumber & WITHSYMMETRIZE ) == WITHSYMMETRIZE; - if ( ( optionnumber & TOPOLOGIESONLY ) == TOPOLOGIESONLY ) - info.flags |= TOPOLOGIESONLY; + // WITHBLOCKS controls output formatting. We could introduce an extra filtering option + // corresponding to GRCC_OPT_Block, which is somewhat like Qgraf "onevi" but not quite + // the same currently. + //opt->values[GRCC_OPT_Block] = ; + + // Now the "qgraf-compatible filtering options": + int qgopt[GRCC_QGRAF_OPT_Size]; + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_ONEPI], optionnumber,ONEPARTI, "ONEPI_", ONEPARTR, "ONEPR_", 1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_ONSHELL], optionnumber,ONSHELL, "ONSHELL_", OFFSHELL, "OFFSHELL_", 1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_NOSIGMA], optionnumber,NOSIGMA, "NOSIGMA_", SIGMA, "SIGMA_", 1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_NOSNAIL], optionnumber,NOSNAIL, "NOSNAIL_", SNAIL, "SNAIL_", 1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_NOTADPOLE],optionnumber,NOTADPOLE,"NOTADPOLE_",TADPOLE , "TADPOLE_", 1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_SIMPLE], optionnumber,SIMPLE, "SIMPLE_", NOTSIMPLE,"NOTSIMPLE_",1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_BIPART], optionnumber,BIPART, "BIPART_", NONBIPART,"NONBIPART_",1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_CYCLI], optionnumber,CYCLI, "CYCLI_", CYCLR, "CYCLR_", 1,-1); + SetDualOpts(&qgopt[GRCC_QGRAF_OPT_FLOOP], optionnumber,FLOOP, "FLOOP_", NOTFLOOP, "NOTFLOOP_", 1,-1); + // Now set the options internally: + opt->setQGrafOpt(qgopt); opt->setOutputF(False,""); + opt->setOutputP(False,""); opt->printLevel(babble); // Load the various arrays. - - ninitl = Sets[inset].last - Sets[inset].first; - for ( i = 0; i < ninitl; i++ ) { - x = SetElements[Sets[inset].first+i]; - initlPart[i] = ConvertParticle(model,x); - } - nfinal = Sets[outset].last - Sets[outset].first; - for ( i = 0; i < nfinal; i++ ) { - x = SetElements[Sets[outset].first+i]; - finalPart[i] = ConvertParticle(model,x); - } + ninitl = Sets[inset].last - Sets[inset].first; + for ( i = 0; i < ninitl; i++ ) { + x = SetElements[Sets[inset].first+i]; + initlPart[i] = ConvertParticle(model,x); + info.legcouple[i] = m->vertices[numParticle(m,x)]->couplings; + } + nfinal = Sets[outset].last - Sets[outset].first; + for ( i = 0; i < nfinal; i++ ) { + x = SetElements[Sets[outset].first+i]; + finalPart[i] = ConvertParticle(model,x); + info.legcouple[i+ninitl] = m->vertices[numParticle(m,x)]->couplings; + } info.numextern = ninitl + nfinal; + for ( i = 2; i <= MAXLEGS; i++ ) { + if ( m->legcouple[i] == 1 ) { + for ( j = 0; j < info.numextern; j++ ) { + if ( info.legcouple[j][i] == 0 ) { info.flags |= CHECKEXTERN; goto Go_on; } + } + } + } +Go_on:; // // Now we have to sort out the coupling constants. // The argument at coupl can be of type -SNUMBER, -SYMBOL or generic @@ -529,11 +871,18 @@ int GenDiagrams(PHEAD WORD *term, WORD level) int nc = coupl[1]*2 + ninitl + nfinal - 2; int *scratch = (int *)Malloc1(nc*sizeof(int),"DistrN"); scratch[0] = -2; // indicating startup cq first call. - while ( DistrN(nc,cpl,m->ncouplings,scratch) ) { - proc = new Process(pid, model, opt, - ninitl, initlPart, nfinal, finalPart, cpl); + + if ( ( info.flags & TOPOLOGIESONLY ) == 0 ) { + while ( DistrN(nc,cpl,m->ncouplings,scratch) ) { + proc = new Process(pid, model, opt, ninitl, initlPart, nfinal, finalPart, cpl); + delete proc; + info.numtopo = 1; + } + } + else { + cpl[0] = nc; + proc = new Process(pid, model, opt, ninitl, initlPart, nfinal, finalPart, cpl); delete proc; - info.numtopo = 1; } M_free(scratch,"DistrN"); opt->end(); @@ -562,8 +911,10 @@ int GenDiagrams(PHEAD WORD *term, WORD level) t += 2; } } - proc = new Process(pid, model, opt, - ninitl, initlPart, nfinal, finalPart, cpl); +/* + And now the generation: +*/ + proc = new Process(pid, model, opt, ninitl, initlPart, nfinal, finalPart, cpl); opt->end(); delete proc; delete opt; @@ -596,15 +947,11 @@ int processVertex(TOPOTYPE *TopoInf, int pointsremaining, int level) TopoInf->cldeg[TopoInf->ncl] = TopoInf->vert[level]; TopoInf->clnum[TopoInf->ncl] = j; TopoInf->clext[TopoInf->ncl] = 0; - TopoInf->cmind[TopoInf->ncl] = 0; - TopoInf->cmaxd[TopoInf->ncl] = 0; TopoInf->ncl++; - MGraph *mgraph = new MGraph(1, - TopoInf->ncl, TopoInf->cldeg, - TopoInf->clnum, TopoInf->clext, - TopoInf->cmind, TopoInf->cmaxd, - TopoInf->opt); + MGraph *mgraph = new MGraph(1, TopoInf->ncl, TopoInf->cldeg, + TopoInf->clnum, TopoInf->clext, + TopoInf->cmind, TopoInf->cmaxd, TopoInf->opt); mgraph->generate(); @@ -663,6 +1010,8 @@ int GenTopologies(PHEAD WORD *term, WORD level) info.numextern = nlegs; + for ( i = 0; i <= MAXLEGS; i++ ) { TopoInf.cmind[i] = TopoInf.cmaxd[i] = 0; } + t1 += 10; if ( t1 < tstop && t1[0] == -SETSET ) { TopoInf.vertmax = &(SetElements[Sets[t1[1]].first]); @@ -672,21 +1021,30 @@ int GenTopologies(PHEAD WORD *term, WORD level) info.flags |= TOPOLOGIESONLY; // this is the topologies_ function after all. if ( t1 < tstop && t1[0] == -SNUMBER ) { - if ( ( t1[1] & NONODES ) == NONODES ) info.flags |= NONODES; + if ( ( t1[1] & WITHOUTNODES ) == WITHOUTNODES ) info.flags |= WITHOUTNODES; if ( ( t1[1] & WITHEDGES ) == WITHEDGES ) info.flags |= WITHEDGES; - opt->values[GRCC_OPT_1PI] = ( t1[1] & ONEPARTICLEIRREDUCIBLE ) == ONEPARTICLEIRREDUCIBLE; - opt->values[GRCC_OPT_NoTadpole] = ( t1[1] & NOTADPOLES ) == NOTADPOLES; - opt->values[GRCC_OPT_No1PtBlock] = ( t1[1] & NOTADPOLES ) == NOTADPOLES; - opt->values[GRCC_OPT_No2PtL1PI] = ( t1[1] & WITHINSERTIONS ) == WITHINSERTIONS; - opt->values[GRCC_OPT_SymmInitial] = ( t1[1] & WITHSYMMETRIZE ) == WITHSYMMETRIZE; + if ( ( t1[1] & WITHBLOCKS ) == WITHBLOCKS ) info.flags |= WITHBLOCKS; + if ( ( t1[1] & WITHONEPISETS ) == WITHONEPISETS ) info.flags |= WITHONEPISETS; + opt->values[GRCC_OPT_1PI] = ( t1[1] & ONEPARTI ) == ONEPARTI; +// opt->values[GRCC_OPT_NoTadpole] = ( t1[1] & NOTADPOLE ) == NOTADPOLE; + opt->values[GRCC_OPT_NoTadpole] = ( t1[1] & NOSNAIL ) == NOSNAIL; + opt->values[GRCC_OPT_No1PtBlock] = ( t1[1] & NOTADPOLE ) == NOTADPOLE; +// opt->values[GRCC_OPT_NoExtSelf] = ( t1[1] & NOEXTSELF ) == NOEXTSELF; + +// if ( ( t1[1] & WITHINSERTIONS ) == WITHINSERTIONS ) { +// opt->values[GRCC_OPT_No2PtL1PI] = True; +// opt->values[GRCC_OPT_NoAdj2PtV] = True; +// opt->values[GRCC_OPT_No2PtL1PI] = True; +// } + opt->values[GRCC_OPT_SymmInitial] = ( t1[1] & WITHSYMMETRIZEI ) == WITHSYMMETRIZEI; + opt->values[GRCC_OPT_SymmFinal] = ( t1[1] & WITHSYMMETRIZEF ) == WITHSYMMETRIZEF; } info.numdia = 0; info.numtopo = 1; - opt->setOutAG(ProcessDiagram, (void*) &info); - opt->setOutMG(ProcessTopology, (void*) &info); - + opt->setOutAG(ProcessDiagram, &info); + opt->setOutMG(ProcessTopology, &info); // // Now we should sum over all possible vertices and run MGraph for // each combination. This is done by recursion in the processVertex routine @@ -701,7 +1059,6 @@ int GenTopologies(PHEAD WORD *term, WORD level) } for ( i = 0; i < nlegs; i++ ) { TopoInf.cldeg[i] = 1; TopoInf.clnum[i] = 1; TopoInf.clext[i] = -1; - TopoInf.cmind[i] = 0; TopoInf.cmaxd[i] = 0; } int points = 2*nloops-2+nlegs; diff --git a/sources/ftypes.h b/sources/ftypes.h index 80246c44..be1a8439 100644 --- a/sources/ftypes.h +++ b/sources/ftypes.h @@ -1102,16 +1102,36 @@ typedef int (*TFUN1)(); #define DENSETABLE 1 #define SPARSETABLE 0 -#define ONEPARTICLEIRREDUCIBLE 1 -#define WITHINSERTIONS 2 -#define NOTADPOLES 4 -#define WITHSYMMETRIZE 8 -#define TOPOLOGIESONLY 16 -#define NONODES 32 -#define WITHEDGES 64 -#define CHECKEXTERN 128 -#define WITHBLOCKS 256 -#define WITHONEPI 512 -#define NOSNAILS 1024 -#define NOEXTSELF 2048 +// Diagram generator flags. They should be powers of two, since they are added +// to pass to diagrams and masked in diawrap.cc to check the flags. +// We use a "stringified" version of these when defining the FORM preprocessor +// variables, so can't neatly use "1<<10" etc here. +#define TOPOLOGIESONLY 1 +#define WITHOUTNODES 2 +#define WITHEDGES 4 +#define WITHBLOCKS 8 +#define WITHONEPISETS 16 +#define WITHSYMMETRIZEI 32 +#define WITHSYMMETRIZEF 64 +// This is not an "option" preprocessor var but is set by "external" particle definitions +#define CHECKEXTERN 128 +// The "qgraf compatible" filtering flags: +#define ONEPARTI 256 +#define ONEPARTR 512 +#define ONSHELL 1024 +#define OFFSHELL 2048 +#define NOSIGMA 4096 +#define SIGMA 8192 +#define NOSNAIL 16384 +#define SNAIL 32768 +#define NOTADPOLE 65536 +#define TADPOLE 131072 +#define SIMPLE 262144 +#define NOTSIMPLE 524288 +#define BIPART 1048576 +#define NONBIPART 2097152 +#define CYCLI 4194304 +#define CYCLR 8388608 +#define FLOOP 16777216 +#define NOTFLOOP 33554432 diff --git a/sources/grcc.cc b/sources/grcc.cc index 18089d89..8272f86e 100644 --- a/sources/grcc.cc +++ b/sources/grcc.cc @@ -7242,16 +7242,22 @@ Bool EGraph::optQGrafA(Options *opt) printf("optQGrafA: %8ld\n", mId); econn->print(); #endif + Bool retval = True; if (opt->qgopt[GRCC_QGRAF_OPT_FLOOP] != 0) { for (int fl=0; fl < nFlines; fl++) { if (flines[fl]->ftype == FL_Closed) { if (flines[fl]->nlist % 2 != 0) { - return False; + retval = False; + break; } } } + // `notfloop_' is the dual of `floop_', so flip the decision: + if (opt->qgopt[GRCC_QGRAF_OPT_FLOOP] == -1) { + retval = (retval == True ? False : True); + } } - return True; + return retval; } //-------------------------------------------------------------- diff --git a/sources/grccparam.h b/sources/grccparam.h index e769b0b5..88d03e5d 100644 --- a/sources/grccparam.h +++ b/sources/grccparam.h @@ -129,7 +129,8 @@ typedef struct { #define GRCC_QGRAF_OPT_BIPART 6 #define GRCC_QGRAF_OPT_CYCLI 7 #define GRCC_QGRAF_OPT_FLOOP 8 -#define GRCC_QGRAF_OPT_TOPOL 9 +//TODO what does this do? Does it work? +//#define GRCC_QGRAF_OPT_TOPOL 9 #ifdef GRCC_QGRAF_OPT_TOPOL #define GRCC_QGRAF_OPT_Size 10 diff --git a/sources/startup.c b/sources/startup.c index 7de4474b..1b39678f 100644 --- a/sources/startup.c +++ b/sources/startup.c @@ -1205,20 +1205,37 @@ void StartVariables(void) PutPreVar((UBYTE *)"keepright_",(UBYTE *)("0"),0,0); PutPreVar((UBYTE *)"SYSTEMERROR_",(UBYTE *)("0"),0,0); /* - Next are a few 'constants' for diagram generation + Next are the flags to control diagram generation filters */ - PutPreVar((UBYTE *)"ONEPI_",(UBYTE *)("1"),0,0); - PutPreVar((UBYTE *)"WITHOUTINSERTIONS_",(UBYTE *)("2"),0,0); - PutPreVar((UBYTE *)"NOTADPOLES_",(UBYTE *)("4"),0,0); - PutPreVar((UBYTE *)"SYMMETRIZE_",(UBYTE *)("8"),0,0); - PutPreVar((UBYTE *)"TOPOLOGIESONLY_",(UBYTE *)("16"),0,0); - PutPreVar((UBYTE *)"NONODES_",(UBYTE *)("32"),0,0); - PutPreVar((UBYTE *)"WITHEDGES_",(UBYTE *)("64"),0,0); -/* Note that CHECKEXTERN is 128 */ - PutPreVar((UBYTE *)"WITHBLOCKS_",(UBYTE *)("256"),0,0); - PutPreVar((UBYTE *)"WITHONEPISETS_",(UBYTE *)("512"),0,0); - PutPreVar((UBYTE *)"NOSNAILS_",(UBYTE *)("1024"),0,0); - PutPreVar((UBYTE *)"NOEXTSELF_",(UBYTE *)("2048"),0,0); + #define STR2(x) #x + #define STR(x) STR2(x) + PutPreVar((UBYTE *)"TOPOLOGIESONLY_" ,(UBYTE*)(STR(TOPOLOGIESONLY)) ,0,0); + PutPreVar((UBYTE *)"WITHOUTNODES_" ,(UBYTE*)(STR(WITHOUTNODES)) ,0,0); + PutPreVar((UBYTE *)"WITHEDGES_" ,(UBYTE*)(STR(WITHEDGES)) ,0,0); + PutPreVar((UBYTE *)"WITHBLOCKS_" ,(UBYTE*)(STR(WITHBLOCKS)) ,0,0); + PutPreVar((UBYTE *)"WITHONEPISETS_" ,(UBYTE*)(STR(WITHONEPISETS)) ,0,0); + PutPreVar((UBYTE *)"WITHSYMMETRIZEI_",(UBYTE*)(STR(WITHSYMMETRIZEI)),0,0); + PutPreVar((UBYTE *)"WITHSYMMETRIZEF_",(UBYTE*)(STR(WITHSYMMETRIZEF)),0,0); + // This is not an "option" preprocessor var but is set by "external" particle definitions +// PutPreVar((UBYTE *)"CHECKEXTERN_" ,(UBYTE*)(STR(CHECKEXTERN)) ,0,0); + PutPreVar((UBYTE *)"ONEPI_" ,(UBYTE*)(STR(ONEPARTI)) ,0,0); + PutPreVar((UBYTE *)"ONEPR_" ,(UBYTE*)(STR(ONEPARTR)) ,0,0); + PutPreVar((UBYTE *)"ONSHELL_" ,(UBYTE*)(STR(ONSHELL)) ,0,0); + PutPreVar((UBYTE *)"OFFSHELL_" ,(UBYTE*)(STR(OFFSHELL)) ,0,0); + PutPreVar((UBYTE *)"NOSIGMA_" ,(UBYTE*)(STR(NOSIGMA)) ,0,0); + PutPreVar((UBYTE *)"SIGMA_" ,(UBYTE*)(STR(SIGMA)) ,0,0); + PutPreVar((UBYTE *)"NOSNAIL_" ,(UBYTE*)(STR(NOSNAIL)) ,0,0); + PutPreVar((UBYTE *)"SNAIL_" ,(UBYTE*)(STR(SNAIL)) ,0,0); + PutPreVar((UBYTE *)"NOTADPOLE_" ,(UBYTE*)(STR(NOTADPOLE)) ,0,0); + PutPreVar((UBYTE *)"TADPOLE_" ,(UBYTE*)(STR(TADPOLE)) ,0,0); + PutPreVar((UBYTE *)"SIMPLE_" ,(UBYTE*)(STR(SIMPLE)) ,0,0); + PutPreVar((UBYTE *)"NOTSIMPLE_" ,(UBYTE*)(STR(NOTSIMPLE)) ,0,0); + PutPreVar((UBYTE *)"BIPART_" ,(UBYTE*)(STR(BIPART)) ,0,0); + PutPreVar((UBYTE *)"NONBIPART_" ,(UBYTE*)(STR(NONBIPART)) ,0,0); + PutPreVar((UBYTE *)"CYCLI_" ,(UBYTE*)(STR(CYCLI)) ,0,0); + PutPreVar((UBYTE *)"CYCLR_" ,(UBYTE*)(STR(CYCLR)) ,0,0); + PutPreVar((UBYTE *)"FLOOP_" ,(UBYTE*)(STR(FLOOP)) ,0,0); + PutPreVar((UBYTE *)"NOTFLOOP_" ,(UBYTE*)(STR(NOTFLOOP)) ,0,0); { char buf[41]; /* up to 128-bit */