| Trees | Indices | Help |
|---|
|
|
1 #
2 ################################################################################
3 # Copyright (c) 2009 The MadGraph5_aMC@NLO Development team and Contributors
4 #
5 # This file is a part of the MadGraph5_aMC@NLO project, an application which
6 # automatically generates Feynman diagrams and matrix elements for arbitrary
7 # high-energy processes in the Standard Model and beyond.
8 #
9 # It is subject to the MadGraph5_aMC@NLO license which should accompany this
10 # distribution.
11 #
12 # For more information, visit madgraph.phys.ucl.ac.be and amcatnlo.web.cern.ch
13 #
14 ################################################################################
15 """Classes for diagram generation with loop features.
16 """
17
18 import array
19 import copy
20 import itertools
21 import logging
22
23 import madgraph.loop.loop_base_objects as loop_base_objects
24 import madgraph.core.base_objects as base_objects
25 import madgraph.core.diagram_generation as diagram_generation
26 import madgraph.various.misc as misc
27
28 from madgraph import MadGraph5Error
29 from madgraph import InvalidCmd
30 logger = logging.getLogger('madgraph.loop_diagram_generation')
33 # This subroutine has typically quite large DEBUG info.
34 # So even in debug mode, they are turned off by default.
35 # Remove the line below for loop diagram generation diagnostic
36 if not force: return
37
38 flag = "LoopGenInfo: "
39 if len(msg)>40:
40 logger.debug(flag+msg[:35]+" [...] = %s"%str(val))
41 else:
42 logger.debug(flag+msg+''.join([' ']*(40-len(msg)))+' = %s'%str(val))
43
44 #===============================================================================
45 # LoopAmplitude
46 #===============================================================================
47 -class LoopAmplitude(diagram_generation.Amplitude):
48 """NLOAmplitude: process + list of diagrams (ordered)
49 Initialize with a process, then call generate_diagrams() to
50 generate the diagrams for the amplitude
51 """
52
54 """Default values for all properties"""
55
56 # The 'diagrams' entry from the mother class is inherited but will not
57 # be used in NLOAmplitude, because it is split into the four following
58 # different categories of diagrams.
59 super(LoopAmplitude, self).default_setup()
60 self['born_diagrams'] = None
61 self['loop_diagrams'] = None
62 self['loop_UVCT_diagrams'] = base_objects.DiagramList()
63 # This is in principle equal to self['born_diagram']==[] but it can be
64 # that for some reason the born diagram can be generated but do not
65 # contribute.
66 # This will decide wether the virtual is squared against the born or
67 # itself.
68 self['has_born'] = True
69 # This where the structures obtained for this amplitudes are stored
70 self['structure_repository'] = loop_base_objects.FDStructureList()
71
72 # A list that registers what Lcut particle have already been
73 # employed in order to forbid them as loop particles in the
74 # subsequent diagram generation runs.
75 self.lcutpartemployed=[]
76
78 """Allow initialization with Process.
79 If loop_filter is not None, then it will be applied to all subsequent
80 diagram generation from this LoopAmplitude."""
81
82 self.loop_filter = loop_filter
83
84 if isinstance(argument, base_objects.Process):
85 super(LoopAmplitude, self).__init__()
86 self.set('process', argument)
87 self.generate_diagrams()
88 elif argument != None:
89 # call the mother routine
90 super(LoopAmplitude, self).__init__(argument)
91 else:
92 # call the mother routine
93 super(LoopAmplitude, self).__init__()
94
96 """Return diagram property names as a nicely sorted list."""
97
98 return ['process', 'diagrams', 'has_mirror_process', 'born_diagrams',
99 'loop_diagrams','has_born',
100 'structure_repository']
101
103 """Filter for valid amplitude property values."""
104
105 if name == 'diagrams':
106 if not isinstance(value, base_objects.DiagramList):
107 raise self.PhysicsObjectError, \
108 "%s is not a valid DiagramList" % str(value)
109 for diag in value:
110 if not isinstance(diag,loop_base_objects.LoopDiagram) and \
111 not isinstance(diag,loop_base_objects.LoopUVCTDiagram):
112 raise self.PhysicsObjectError, \
113 "%s contains a diagram which is not an NLODiagrams." % str(value)
114 if name == 'born_diagrams':
115 if not isinstance(value, base_objects.DiagramList):
116 raise self.PhysicsObjectError, \
117 "%s is not a valid DiagramList" % str(value)
118 for diag in value:
119 if not isinstance(diag,loop_base_objects.LoopDiagram):
120 raise self.PhysicsObjectError, \
121 "%s contains a diagram which is not an NLODiagrams." % str(value)
122 if name == 'loop_diagrams':
123 if not isinstance(value, base_objects.DiagramList):
124 raise self.PhysicsObjectError, \
125 "%s is not a valid DiagramList" % str(value)
126 for diag in value:
127 if not isinstance(diag,loop_base_objects.LoopDiagram):
128 raise self.PhysicsObjectError, \
129 "%s contains a diagram which is not an NLODiagrams." % str(value)
130 if name == 'has_born':
131 if not isinstance(value, bool):
132 raise self.PhysicsObjectError, \
133 "%s is not a valid bool" % str(value)
134 if name == 'structure_repository':
135 if not isinstance(value, loop_base_objects.FDStructureList):
136 raise self.PhysicsObjectError, \
137 "%s is not a valid bool" % str(value)
138
139 else:
140 super(LoopAmplitude, self).filter(name, value)
141
142 return True
143
145 """Redefine set for the particular case of diagrams"""
146
147 if name == 'diagrams':
148 if self.filter(name, value):
149 self['born_diagrams']=base_objects.DiagramList([diag for diag in value if \
150 not isinstance(diag,loop_base_objects.LoopUVCTDiagram) and diag['type']==0])
151 self['loop_diagrams']=base_objects.DiagramList([diag for diag in value if \
152 not isinstance(diag,loop_base_objects.LoopUVCTDiagram) and diag['type']!=0])
153 self['loop_UVCT_diagrams']=base_objects.DiagramList([diag for diag in value if \
154 isinstance(diag,loop_base_objects.LoopUVCTDiagram)])
155
156 else:
157 return super(LoopAmplitude, self).set(name, value)
158
159 return True
160
162 """Redefine get for the particular case of '*_diagrams' property"""
163
164 if name == 'diagrams':
165 if self['process'] and self['loop_diagrams'] == None:
166 self.generate_diagrams()
167 return base_objects.DiagramList(self['born_diagrams']+\
168 self['loop_diagrams']+\
169 self['loop_UVCT_diagrams'])
170
171 if name == 'born_diagrams':
172 if self['born_diagrams'] == None:
173 # Have not yet generated born diagrams for this process
174 if self['process']['has_born']:
175 if self['process']:
176 self.generate_born_diagrams()
177 else:
178 self['born_diagrams']=base_objects.DiagramList()
179
180 return LoopAmplitude.__bases__[0].get(self, name) #return the mother routine
181
182 # Functions of the different tasks performed in generate_diagram
184 """ Choose the configuration of non-perturbed coupling orders to be
185 retained for all diagrams. This is used when the user did not specify
186 any order. """
187 chosen_order_config = {}
188 min_wgt = self['born_diagrams'].get_min_order('WEIGHTED')
189 # Scan the born diagrams of minimum weight to chose a configuration
190 # of non-perturbed orders.
191 min_non_pert_order_wgt = -1
192 for diag in [d for d in self['born_diagrams'] if \
193 d.get_order('WEIGHTED')==min_wgt]:
194 non_pert_order_wgt = min_wgt - sum([diag.get_order(order)*\
195 self['process']['model']['order_hierarchy'][order] for order in \
196 self['process']['perturbation_couplings']])
197 if min_non_pert_order_wgt == -1 or \
198 non_pert_order_wgt<min_non_pert_order_wgt:
199 chosen_order_config = self.get_non_pert_order_config(diag)
200 logger.info("Chosen coupling orders configuration: (%s)"\
201 %self.print_config(chosen_order_config))
202 return chosen_order_config
203
205 """If squared orders (other than WEIGHTED) are defined, then they can be
206 used for determining what is the expected upper bound for the order
207 restricting loop diagram generation."""
208 for order, value in self['process']['squared_orders'].items():
209 if order.upper()!='WEIGHTED' and order not in self['process']['orders']:
210 # If the bound is of type '>' we cannot say anything
211 if self['process'].get('sqorders_types')[order]=='>':
212 continue
213 # If there is no born, the min order will simply be 0 as it should.
214 bornminorder=self['born_diagrams'].get_min_order(order)
215 if value>=0:
216 self['process']['orders'][order]=value-bornminorder
217 elif self['process']['has_born']:
218 # This means the user want the leading if order=-1 or N^n
219 # Leading term if order=-n. If there is a born diag, we can
220 # infer the necessary maximum order in the loop:
221 # bornminorder+2*(n-1).
222 # If there is no born diag, then we cannot say anything.
223 self['process']['orders'][order]=bornminorder+2*(-value-1)
224
226 """Guess the upper bound for the orders for loop diagram generation
227 based on either no squared orders or simply 'Weighted'"""
228
229 hierarchy = self['process']['model']['order_hierarchy']
230
231 # Maximum of the hierarchy weigtht among all perturbed order
232 max_pert_wgt = max([hierarchy[order] for order in \
233 self['process']['perturbation_couplings']])
234
235 # In order to be sure to catch the corrections to all born diagrams that
236 # the user explicitly asked for with the amplitude orders, we take here
237 # the minimum weighted order as being the maximum between the min weighted
238 # order detected in the Born diagrams and the weight computed from the
239 # user input amplitude orders.
240 user_min_wgt = 0
241
242 # One can chose between the two behaviors below. It is debatable which
243 # one is best. The first one tries to only consider the loop which are
244 # dominant, even when the user selects the amplitude orders and the
245 # second chosen here makes sure that the user gets a correction of the
246 # desired type for all the born diagrams generated with its amplitude
247 # order specification.
248 # min_born_wgt=self['born_diagrams'].get_min_order('WEIGHTED')
249 min_born_wgt=max(self['born_diagrams'].get_min_order('WEIGHTED'),
250 sum([hierarchy[order]*val for order, val in user_orders.items() \
251 if order!='WEIGHTED']))
252
253 if 'WEIGHTED' not in [key.upper() for key in \
254 self['process']['squared_orders'].keys()]:
255 # Then we guess it from the born
256 self['process']['squared_orders']['WEIGHTED']= 2*(min_born_wgt+\
257 max_pert_wgt)
258
259 # Now we know that the remaining weighted orders which can fit in
260 # the loop diagram is (self['target_weighted_order']-
261 # min_born_weighted_order) so for each perturbed order we just have to
262 # take that number divided by its hierarchy weight to have the maximum
263 # allowed order for the loop diagram generation. Of course,
264 # we don't overwrite any order already defined by the user.
265 if self['process']['squared_orders']['WEIGHTED']>=0:
266 trgt_wgt=self['process']['squared_orders']['WEIGHTED']-min_born_wgt
267 else:
268 trgt_wgt=min_born_wgt+(-self['process']['squared_orders']['WEIGHTED']+1)*2
269 # We also need the minimum number of vertices in the born.
270 min_nvert=min([len([1 for vert in diag['vertices'] if vert['id']!=0]) \
271 for diag in self['born_diagrams']])
272 # And the minimum weight for the ordered declared as perturbed
273 min_pert=min([hierarchy[order] for order in \
274 self['process']['perturbation_couplings']])
275
276 for order, value in hierarchy.items():
277 if order not in self['process']['orders']:
278 # The four cases below come from a study of the maximal order
279 # needed in the loop for the weighted order needed and the
280 # number of vertices available.
281 if order in self['process']['perturbation_couplings']:
282 if value!=1:
283 self['process']['orders'][order]=\
284 int((trgt_wgt-min_nvert-2)/(value-1))
285 else:
286 self['process']['orders'][order]=int(trgt_wgt)
287 else:
288 if value!=1:
289 self['process']['orders'][order]=\
290 int((trgt_wgt-min_nvert-2*min_pert)/(value-1))
291 else:
292 self['process']['orders'][order]=\
293 int(trgt_wgt-2*min_pert)
294 # Now for the remaining orders for which the user has not set squared
295 # orders neither amplitude orders, we use the max order encountered in
296 # the born (and add 2 if this is a perturbed order).
297 # It might be that this upper bound is better than the one guessed
298 # from the hierarchy.
299 for order in self['process']['model']['coupling_orders']:
300 neworder=self['born_diagrams'].get_max_order(order)
301 if order in self['process']['perturbation_couplings']:
302 neworder+=2
303 if order not in self['process']['orders'].keys() or \
304 neworder<self['process']['orders'][order]:
305 self['process']['orders'][order]=neworder
306
308 """ Filter diags to select only the diagram with the non perturbed orders
309 configuration config and update discarded_configurations.Diags is the
310 name of the key attribute of this class containing the diagrams to
311 filter."""
312 newdiagselection = base_objects.DiagramList()
313 for diag in self[diags]:
314 diag_config = self.get_non_pert_order_config(diag)
315 if diag_config == config:
316 newdiagselection.append(diag)
317 elif diag_config not in discarded_configurations:
318 discarded_configurations.append(diag_config)
319 self[diags] = newdiagselection
320
322 """ Remove the loops which are zero because of Furry theorem. So as to
323 limit any possible mistake in case of BSM model, I limit myself here to
324 removing SM-quark loops with external legs with an odd number of photons,
325 possibly including exactly two gluons."""
326
327 new_diag_selection = base_objects.DiagramList()
328
329 n_discarded = 0
330 for diag in self['loop_diagrams']:
331 if diag.get('tag')==[]:
332 raise MadGraph5Error, "The loop diagrams should have been tagged"+\
333 " before going through the Furry filter."
334
335 loop_line_pdgs = diag.get_loop_lines_pdgs()
336 attached_pdgs = diag.get_pdgs_attached_to_loop(structs)
337 if (attached_pdgs.count(22)%2==1) and \
338 (attached_pdgs.count(21) in [0,2]) and \
339 (all(pdg in [22,21] for pdg in attached_pdgs)) and \
340 (abs(loop_line_pdgs[0]) in list(range(1,7))) and \
341 (all(abs(pdg)==abs(loop_line_pdgs[0]) for pdg in loop_line_pdgs)):
342 n_discarded += 1
343 else:
344 new_diag_selection.append(diag)
345
346 self['loop_diagrams'] = new_diag_selection
347
348 if n_discarded > 0:
349 logger.debug(("MadLoop discarded %i diagram%s because they appeared"+\
350 " to be zero because of Furry theorem.")%(n_discarded,'' if \
351 n_discarded<=1 else 's'))
352
353 @staticmethod
355 """ Returns a function which applies the filter corresponding to the
356 conditional expression encoded in filterdef."""
357
358 def filter(diag, structs, model, id):
359 """ The filter function generated '%s'."""%filterdef
360
361 loop_pdgs = diag.get_loop_lines_pdgs()
362 struct_pdgs = diag.get_pdgs_attached_to_loop(structs)
363 loop_masses = [model.get_particle(pdg).get('mass') for pdg in loop_pdgs]
364 struct_masses = [model.get_particle(pdg).get('mass') for pdg in struct_pdgs]
365 if not eval(filterdef.lower(),{'n':len(loop_pdgs),
366 'loop_pdgs':loop_pdgs,
367 'struct_pdgs':struct_pdgs,
368 'loop_masses':loop_masses,
369 'struct_masses':struct_masses,
370 'id':id}):
371 return False
372 else:
373 return True
374
375 return filter
376
378 """ User-defined user-filter. By default it is not called, but the expert
379 user can turn it on and code here is own filter. Some default examples
380 are provided here.
381 The tagging of the loop diagrams must be performed before using this
382 user loop filter"""
383
384 # By default the user filter does nothing if filter is not set,
385 # if you want to turn it on and edit it by hand, then set the
386 # variable edit_filter_manually to True
387 edit_filter_manually = False
388 if not edit_filter_manually and filter in [None,'None']:
389 return
390 if isinstance(filter,str) and filter.lower() == 'true':
391 edit_filter_manually = True
392 filter=None
393
394
395 if filter not in [None,'None']:
396 filter_func = LoopAmplitude.get_loop_filter(filter)
397 else:
398 filter_func = None
399
400 new_diag_selection = base_objects.DiagramList()
401 discarded_diags = base_objects.DiagramList()
402 i=0
403 for diag in self['loop_diagrams']:
404 if diag.get('tag')==[]:
405 raise MadGraph5Error, "Before using the user_filter, please "+\
406 "make sure that the loop diagrams have been tagged first."
407 valid_diag = True
408 i=i+1
409
410 # Apply the custom filter specified if any
411 if filter_func:
412 try:
413 valid_diag = filter_func(diag, structs, model, i)
414 except Exception as e:
415 raise InvalidCmd("The user-defined filter '%s' did not"%filter+
416 " returned the following error:\n > %s"%str(e))
417 # if any([abs(i)!=1000021 for i in diag.get_loop_lines_pdgs()]):
418 # valid_diag=False
419 # if len(diag.get_loop_lines_pdgs())<4:
420 # valid_diag = False
421
422 # Ex. 0: Chose a specific diagram number, here the 8th one for ex.
423 # if i not in [31]:
424 # valid_diag = False
425
426 # Ex. 0: Keeps only the top quark loops.
427 # if any([pdg not in [6,-6] for pdg in diag.get_loop_lines_pdgs()]):
428 # valid_diag = False
429
430 # Ex. 1: Chose the topology, i.e. number of loop line.
431 # Notice that here particles and antiparticles are not
432 # differentiated and always the particle PDG is returned.
433 # In this example, only boxes are selected.
434 # if len(diag.get_loop_lines_pdgs())>2 and \
435 # any([i in diag.get_loop_lines_pdgs() for i in[24,-24,23]]):
436 # valid_diag=False
437
438 # Ex. 2: Use the pdgs of the particles directly attached to the loop.
439 # In this example, we forbid the Z to branch off the loop.
440 # connected_id = diag.get_pdgs_attached_to_loop(structs)
441 # if 22 not in connected_id:
442 # valid_diag=False
443
444 # Ex. 3: Filter based on the mass of the particles running in the
445 # loop. It shows how to access the particles properties from
446 # the PDG.
447 # In this example, only massive parts. are allowed in the loop.
448 # if 'ZERO' in [model.get_particle(pdg).get('mass') for pdg in \
449 # diag.get_loop_lines_pdgs()]:
450 # valid_diag=False
451
452 # Ex. 4: Complicated filter which gets rid of all bubble diagrams made
453 # of two vertices being the four gluon vertex and the effective
454 # glu-glu-Higgs vertex.
455 # if len(diag.get_loop_lines_pdgs())==2:
456 # bubble_lines_pdgs=[abs(diag.get('canonical_tag')[0][0]),
457 # abs(diag.get('canonical_tag')[0][0])]
458 # first_vertex_pdgs=bubble_lines_pdgs+\
459 # [abs(structs.get_struct(struct_ID).get('binding_leg').get('id')) \
460 # for struct_ID in diag.get('canonical_tag')[0][1]]
461 # second_vertex_pdgs=bubble_lines_pdgs+\
462 # [abs(structs.get_struct(struct_ID).get('binding_leg').get('id')) \
463 # for struct_ID in diag.get('canonical_tag')[1][1]]
464 # first_vertex_pdgs.sort()
465 # second_vertex_pdgs.sort()
466 # bubble_vertices=[first_vertex_pdgs,second_vertex_pdgs]
467 # bubble_vertices.sort()
468 # if bubble_vertices==[[21,21,21,21],[21,21,25]]:
469 # valid_diag=False
470
471 # If you need any more advanced function for your filter and cannot
472 # figure out how to implement them. Just contact the authors.
473
474 if valid_diag:
475 new_diag_selection.append(diag)
476 else:
477 discarded_diags.append(diag)
478
479 self['loop_diagrams'] = new_diag_selection
480 if filter in [None,'None']:
481 warn_msg = """
482 The user-defined loop diagrams filter is turned on and discarded %d loops."""\
483 %len(discarded_diags)
484 else:
485 warn_msg = """
486 The loop diagrams filter '%s' is turned on and discarded %d loops."""\
487 %(filter,len(discarded_diags))
488 logger.warning(warn_msg)
489
491 """ Filter the loop diagrams to make sure they belong to the class
492 of coupling orders perturbed. """
493
494 # First define what are the set of particles allowed to run in the loop.
495 allowedpart=[]
496 for part in self['process']['model']['particles']:
497 for order in self['process']['perturbation_couplings']:
498 if part.is_perturbating(order,self['process']['model']):
499 allowedpart.append(part.get_pdg_code())
500 break
501
502 newloopselection=base_objects.DiagramList()
503 warned=False
504 warning_msg = ("Some loop diagrams contributing to this process"+\
505 " are discarded because they are not pure (%s)-perturbation.\nMake sure"+\
506 " you did not want to include them.")%\
507 ('+'.join(self['process']['perturbation_couplings']))
508 for i,diag in enumerate(self['loop_diagrams']):
509 # Now collect what are the coupling orders building the loop which
510 # are also perturbed order.
511 loop_orders=diag.get_loop_orders(self['process']['model'])
512 pert_loop_order=set(loop_orders.keys()).intersection(\
513 set(self['process']['perturbation_couplings']))
514 # Then make sure that the particle running in the loop for all
515 # diagrams belong to the set above. Also make sure that there is at
516 # least one coupling order building the loop which is in the list
517 # of the perturbed order.
518 valid_diag=True
519 if (diag.get_loop_line_types()-set(allowedpart))!=set() or \
520 pert_loop_order==set([]):
521 valid_diag=False
522 if not warned:
523 logger.warning(warning_msg)
524 warned=True
525 if len([col for col in [
526 self['process'].get('model').get_particle(pdg).get('color') \
527 for pdg in diag.get_pdgs_attached_to_loop(\
528 self['structure_repository'])] if col!=1])==1:
529 valid_diag=False
530
531 if valid_diag:
532 newloopselection.append(diag)
533 self['loop_diagrams']=newloopselection
534 # To monitor what are the diagrams filtered, simply comment the line
535 # directly above and uncomment the two directly below.
536 # self['loop_diagrams'] = base_objects.DiagramList(
537 # [diag for diag in self['loop_diagrams'] if diag not in newloopselection])
538
540 """ Makes sure that all non perturbed orders factorize the born diagrams
541 """
542 warning_msg = "All Born diagrams do not factorize the same sum of power(s) "+\
543 "of the the perturbed order(s) %s.\nThis is potentially dangerous"+\
544 " as the real-emission diagrams from aMC@NLO will not be consistent"+\
545 " with these virtual contributions."
546 if self['process']['has_born']:
547 trgt_summed_order = sum([self['born_diagrams'][0].get_order(order)
548 for order in self['process']['perturbation_couplings']])
549 for diag in self['born_diagrams'][1:]:
550 if sum([diag.get_order(order) for order in self['process']
551 ['perturbation_couplings']])!=trgt_summed_order:
552 logger.warning(warning_msg%' '.join(self['process']
553 ['perturbation_couplings']))
554 break
555
556 warning_msg = "All born diagrams do not factorize the same power of "+\
557 "the order %s which is not perturbed and for which you have not"+\
558 "specified any amplitude order. \nThis is potentially dangerous"+\
559 " as the real-emission diagrams from aMC@NLO will not be consistent"+\
560 " with these virtual contributions."
561 if self['process']['has_born']:
562 for order in self['process']['model']['coupling_orders']:
563 if order not in self['process']['perturbation_couplings'] and \
564 order not in user_orders.keys():
565 order_power=self['born_diagrams'][0].get_order(order)
566 for diag in self['born_diagrams'][1:]:
567 if diag.get_order(order)!=order_power:
568 logger.warning(warning_msg%order)
569 break
570
571 # Helper function
573 """ Return a dictionary of all the coupling orders of this diagram which
574 are not the perturbed ones."""
575 return dict([(order, diagram.get_order(order)) for \
576 order in self['process']['model']['coupling_orders'] if \
577 not order in self['process']['perturbation_couplings'] ])
578
580 """Return a string describing the coupling order configuration"""
581 res = []
582 for order in self['process']['model']['coupling_orders']:
583 try:
584 res.append('%s=%d'%(order,config[order]))
585 except KeyError:
586 res.append('%s=*'%order)
587 return ','.join(res)
588
590 """ Generates all diagrams relevant to this Loop Process """
591
592 # Description of the algorithm to guess the leading contribution.
593 # The summed weighted order of each diagram will be compared to
594 # 'target_weighted_order' which acts as a threshold to decide which
595 # diagram to keep. Here is an example on how MG5 sets the
596 # 'target_weighted_order'.
597 #
598 # In the sm process uu~ > dd~ [QCD, QED] with hierarchy QCD=1, QED=2 we
599 # would have at leading order contribution like
600 # (QED=4) , (QED=2, QCD=2) , (QCD=4)
601 # leading to a summed weighted order of respectively
602 # (4*2=8) , (2*2+2*1=6) , (4*1=4)
603 # at NLO in QCD and QED we would have the following possible contributions
604 # (QED=6), (QED=4,QCD=2), (QED=2,QCD=4) and (QCD=6)
605 # which translate into the following weighted orders, respectively
606 # 12, 10, 8 and 6
607 # So, now we take the largest weighted order at born level, 4, and add two
608 # times the largest weight in the hierarchy among the order for which we
609 # consider loop perturbation, in this case 2*2 wich gives us a
610 # target_weighted_order of 8. based on this we will now keep all born
611 # contributions and exclude the NLO contributions (QED=6) and (QED=4,QCD=2)
612
613 # Use the globally defined loop_filter if the locally defined one is empty
614 if (not self.loop_filter is None) and (loop_filter is None):
615 loop_filter = self.loop_filter
616
617 logger.debug("Generating %s "\
618 %self['process'].nice_string().replace('Process', 'process'))
619
620 # Hierarchy and model shorthands
621 model = self['process']['model']
622 hierarchy = model['order_hierarchy']
623
624 # Later, we will specify the orders for the loop amplitude.
625 # It is a temporary change that will be reverted after loop diagram
626 # generation. We then back up here its value prior modification.
627 user_orders=copy.copy(self['process']['orders'])
628 # First generate the born diagram if the user asked for it
629 if self['process']['has_born']:
630 bornsuccessful = self.generate_born_diagrams()
631 ldg_debug_info("# born diagrams after first generation",\
632 len(self['born_diagrams']))
633 else:
634 self['born_diagrams'] = base_objects.DiagramList()
635 bornsuccessful = True
636 logger.debug("Born diagrams generation skipped by user request.")
637
638 # Make sure that all orders specified belong to the model:
639 for order in self['process']['orders'].keys()+\
640 self['process']['squared_orders'].keys():
641 if not order in model.get('coupling_orders') and \
642 order != 'WEIGHTED':
643 raise InvalidCmd("Coupling order %s not found"%order +\
644 " in any interaction of the current model %s."%model['name'])
645
646 # The decision of whether the virtual must be squared against the born or the
647 # virtual is made based on whether there are Born or not unless the user
648 # already asked for the loop squared.
649 if self['process']['has_born']:
650 self['process']['has_born'] = self['born_diagrams']!=[]
651 self['has_born'] = self['process']['has_born']
652
653 ldg_debug_info("User input born orders",self['process']['orders'])
654 ldg_debug_info("User input squared orders",
655 self['process']['squared_orders'])
656 ldg_debug_info("User input perturbation",\
657 self['process']['perturbation_couplings'])
658
659 # Now, we can further specify the orders for the loop amplitude.
660 # Those specified by the user of course remain the same, increased by
661 # two if they are perturbed. It is a temporary change that will be
662 # reverted after loop diagram generation.
663 user_orders=copy.copy(self['process']['orders'])
664 user_squared_orders=copy.copy(self['process']['squared_orders'])
665
666 # If the user did not specify any order, we can expect him not to be an
667 # expert. So we must make sure the born all factorize the same powers of
668 # coupling orders which are not perturbed. If not we chose a configuration
669 # of non-perturbed order which has the smallest total weight and inform
670 # the user about this. It is then stored below for later filtering of
671 # the loop diagrams.
672 chosen_order_config={}
673 if self['process']['squared_orders']=={} and \
674 self['process']['orders']=={} and self['process']['has_born']:
675 chosen_order_config = self.choose_order_config()
676
677 discarded_configurations = []
678 # The born diagrams are now filtered according to the chose configuration
679 if chosen_order_config != {}:
680 self.filter_from_order_config('born_diagrams', \
681 chosen_order_config,discarded_configurations)
682
683 # Before proceeding with the loop contributions, we must make sure that
684 # the born diagram generated factorize the same sum of power of the
685 # perturbed couplings. If this is not true, then it is very
686 # cumbersome to get the real radiation contribution correct and consistent
687 # with the computations of the virtuals (for now).
688 # Also, when MadLoop5 guesses the a loop amplitude order on its own, it
689 # might decide not to include some subleading loop which might be not
690 # be consistently neglected for now in the MadFKS5 so that its best to
691 # warn the user that he should enforce that target born amplitude order
692 # to any value of his choice.
693 self.check_factorization(user_orders)
694
695 # Now find an upper bound for the loop diagram generation.
696 self.guess_loop_orders_from_squared()
697
698 # If the user had not specified any fixed squared order other than
699 # WEIGHTED, we will use the guessed weighted order to assign a bound to
700 # the loop diagram order. Later we will check if the order deduced from
701 # the max order appearing in the born diagrams is a better upper bound.
702 # It will set 'WEIGHTED' to the desired value if it was not already set
703 # by the user. This is why you see the process defined with 'WEIGHTED'
704 # in the squared orders no matter the user input. Leave it like this.
705 if [k.upper() for k in self['process']['squared_orders'].keys()] in \
706 [[],['WEIGHTED']] and self['process']['has_born']:
707 self.guess_loop_orders(user_orders)
708
709 # Finally we enforce the use of the orders specified for the born
710 # (augmented by two if perturbed) by the user, no matter what was
711 # the best guess performed above.
712 for order in user_orders.keys():
713 if order in self['process']['perturbation_couplings']:
714 self['process']['orders'][order]=user_orders[order]+2
715 else:
716 self['process']['orders'][order]=user_orders[order]
717 if 'WEIGHTED' in user_orders.keys():
718 self['process']['orders']['WEIGHTED']=user_orders['WEIGHTED']+\
719 2*min([hierarchy[order] for order in \
720 self['process']['perturbation_couplings']])
721
722 ldg_debug_info("Orders used for loop generation",\
723 self['process']['orders'])
724
725 # Make sure to warn the user if we already possibly excluded mixed order
726 # loops by smartly setting up the orders
727 warning_msg = ("Some loop diagrams contributing to this process might "+\
728 "be discarded because they are not pure (%s)-perturbation.\nMake sure"+\
729 " there are none or that you did not want to include them.")%(\
730 ','.join(self['process']['perturbation_couplings']))
731
732 if self['process']['has_born']:
733 for order in model['coupling_orders']:
734 if order not in self['process']['perturbation_couplings']:
735 try:
736 if self['process']['orders'][order]< \
737 self['born_diagrams'].get_max_order(order):
738 logger.warning(warning_msg)
739 break
740 except KeyError:
741 pass
742
743 # Now we can generate the loop diagrams.
744 totloopsuccessful=self.generate_loop_diagrams()
745
746 # If there is no born neither loop diagrams, return now.
747 if not self['process']['has_born'] and not self['loop_diagrams']:
748 self['process']['orders'].clear()
749 self['process']['orders'].update(user_orders)
750 return False
751
752 # We add here the UV renormalization contribution built in
753 # LoopUVCTDiagram. It is done before the squared order selection because
754 # it is possible that some UV-renorm. diagrams are removed as well.
755 if self['process']['has_born']:
756 self.set_Born_CT()
757
758 ldg_debug_info("#UVCTDiags generated",len(self['loop_UVCT_diagrams']))
759
760 # Reset the orders to their original specification by the user
761 self['process']['orders'].clear()
762 self['process']['orders'].update(user_orders)
763
764 # If there was no born, we will guess the WEIGHT squared order only now,
765 # based on the minimum weighted order of the loop contributions, if it
766 # was not specified by the user.
767 if not self['process']['has_born'] and not \
768 self['process']['squared_orders'] and not\
769 self['process']['orders'] and hierarchy:
770 pert_order_weights=[hierarchy[order] for order in \
771 self['process']['perturbation_couplings']]
772 self['process']['squared_orders']['WEIGHTED']=2*(\
773 self['loop_diagrams'].get_min_order('WEIGHTED')+\
774 max(pert_order_weights)-min(pert_order_weights))
775
776 ldg_debug_info("Squared orders after treatment",\
777 self['process']['squared_orders'])
778 ldg_debug_info("#Diags after diagram generation",\
779 len(self['loop_diagrams']))
780
781
782 # If a special non perturbed order configuration was chosen at the
783 # beginning because of the absence of order settings by the user,
784 # the corresponding filter is applied now to loop diagrams.
785 # List of discarded configurations
786 if chosen_order_config != {}:
787 self.filter_from_order_config('loop_diagrams', \
788 chosen_order_config,discarded_configurations)
789 # # Warn about discarded configurations.
790 if discarded_configurations!=[]:
791 msg = ("The contribution%s of th%s coupling orders "+\
792 "configuration%s %s discarded :%s")%(('s','ese','s','are','\n')\
793 if len(discarded_configurations)>1 else ('','is','','is',' '))
794 msg = msg + '\n'.join(['(%s)'%self.print_config(conf) for conf \
795 in discarded_configurations])
796 msg = msg + "\nManually set the coupling orders to "+\
797 "generate %sthe contribution%s above."%(('any of ','s') if \
798 len(discarded_configurations)>1 else ('',''))
799 logger.info(msg)
800
801 # The minimum of the different orders used for the selections can
802 # possibly increase, after some loop diagrams are selected out.
803 # So this check must be iterated until the number of diagrams
804 # remaining is stable.
805 # We first apply the selection rules without the negative constraint.
806 # (i.e. QCD=1 for LO contributions only)
807 regular_constraints = dict([(key,val) for (key,val) in
808 self['process']['squared_orders'].items() if val>=0])
809 negative_constraints = dict([(key,val) for (key,val) in
810 self['process']['squared_orders'].items() if val<0])
811 while True:
812 ndiag_remaining=len(self['loop_diagrams']+self['born_diagrams'])
813 self.check_squared_orders(regular_constraints)
814 if len(self['loop_diagrams']+self['born_diagrams'])==ndiag_remaining:
815 break
816 # And then only the negative ones
817 if negative_constraints!={}:
818 # It would be meaningless here to iterate because <order>=-X would
819 # have a different meaning every time.
820 # notice that this function will change the negative values of
821 # self['process']['squared_orders'] to their corresponding positive
822 # constraint for the present process.
823 # For example, u u~ > d d~ QCD^2=-2 becomes u u~ > d d~ QCD=2
824 # because the LO QCD contribution has QED=4, QCD=0 and the NLO one
825 # selected with -2 is QED=2, QCD=2.
826 self.check_squared_orders(negative_constraints,user_squared_orders)
827
828 ldg_debug_info("#Diags after constraints",len(self['loop_diagrams']))
829 ldg_debug_info("#Born diagrams after constraints",len(self['born_diagrams']))
830 ldg_debug_info("#UVCTDiags after constraints",len(self['loop_UVCT_diagrams']))
831
832 # Now the loop diagrams are tagged and filtered for redundancy.
833 tag_selected=[]
834 loop_basis=base_objects.DiagramList()
835 for diag in self['loop_diagrams']:
836 diag.tag(self['structure_repository'],model)
837 # Make sure not to consider wave-function renormalization, vanishing tadpoles,
838 # or redundant diagrams
839 if not diag.is_wf_correction(self['structure_repository'], \
840 model) and not diag.is_vanishing_tadpole(model) and \
841 diag['canonical_tag'] not in tag_selected:
842 loop_basis.append(diag)
843 tag_selected.append(diag['canonical_tag'])
844
845 self['loop_diagrams']=loop_basis
846
847 # Now select only the loops corresponding to the perturbative orders
848 # asked for.
849 self.filter_loop_for_perturbative_orders()
850
851 if len(self['loop_diagrams'])==0 and len(self['born_diagrams'])!=0:
852 raise InvalidCmd('All loop diagrams discarded by user selection.\n'+\
853 'Consider using a tree-level generation or relaxing the coupling'+\
854 ' order constraints.')
855 # If there is no born neither loop diagrams after filtering, return now.
856 if not self['process']['has_born'] and not self['loop_diagrams']:
857 self['process']['squared_orders'].clear()
858 self['process']['squared_orders'].update(user_squared_orders)
859 return False
860
861
862 # Discard diagrams which are zero because of Furry theorem
863 self.remove_Furry_loops(model,self['structure_repository'])
864
865 # Apply here some user-defined filter.
866 # For expert only, you can edit your own filter by modifying the
867 # user_filter() function which by default does nothing but in which you
868 # will find examples of common filters.
869 self.user_filter(model,self['structure_repository'], filter=loop_filter)
870
871 # Set the necessary UV/R2 CounterTerms for each loop diagram generated
872 self.set_LoopCT_vertices()
873
874 # Now revert the squared order. This function typically adds to the
875 # squared order list the target WEIGHTED order which has been detected.
876 # This is typically not desired because if the user types in directly
877 # what it sees on the screen, it does not get back the same process.
878 # for example, u u~ > d d~ [virt=QCD] becomes
879 # u u~ > d d~ [virt=QCD] WEIGHTED=6
880 # but of course the photon-gluon s-channel Born interference is not
881 # counted in.
882 # However, if you type it in generate again with WEIGHTED=6, you will
883 # get it.
884 self['process']['squared_orders'].clear()
885 self['process']['squared_orders'].update(user_squared_orders)
886
887 # The computation below is just to report what split order are computed
888 # and which one are considered (i.e. kept using the order specifications)
889 self.print_split_order_infos()
890
891 # Give some info about the run
892 nLoopDiag = 0
893 nCT={'UV':0,'R2':0}
894 for ldiag in self['loop_UVCT_diagrams']:
895 nCT[ldiag['type'][:2]]+=len(ldiag['UVCT_couplings'])
896 for ldiag in self['loop_diagrams']:
897 nLoopDiag+=1
898 nCT['UV']+=len(ldiag.get_CT(model,'UV'))
899 nCT['R2']+=len(ldiag.get_CT(model,'R2'))
900
901 # The identification of numerically equivalent diagrams is done here.
902 # Simply comment the line above to remove it for testing purposes
903 # (i.e. to make sure it does not alter the result).
904 nLoopsIdentified = self.identify_loop_diagrams()
905 if nLoopsIdentified > 0:
906 logger.debug("A total of %d loop diagrams "%nLoopsIdentified+\
907 "were identified with equivalent ones.")
908 logger.info("Contributing diagrams generated: "+\
909 "%d Born, %d%s loops, %d R2, %d UV"%(len(self['born_diagrams']),
910 len(self['loop_diagrams']),'(+%d)'%nLoopsIdentified \
911 if nLoopsIdentified>0 else '' ,nCT['R2'],nCT['UV']))
912
913 ldg_debug_info("#Diags after filtering",len(self['loop_diagrams']))
914 ldg_debug_info("# of different structures identified",\
915 len(self['structure_repository']))
916
917 return (bornsuccessful or totloopsuccessful)
918
920 """ Uses a loop_tag characterizing the loop with only physical
921 information about it (mass, coupling, width, color, etc...) so as to
922 recognize numerically equivalent diagrams and group them together,
923 such as massless quark loops in pure QCD gluon loop amplitudes."""
924
925 # This dictionary contains key-value pairs of the form
926 # (loop_tag, DiagramList) where the loop_tag key unambiguously
927 # characterizes a class of equivalent diagrams and the DiagramList value
928 # lists all the diagrams belonging to this class.
929 # In the end, the first diagram of this DiagramList will be used as
930 # the reference included in the numerical code for the loop matrix
931 # element computations and all the others will be omitted, being
932 # included via a simple multiplicative factor applied to the first one.
933 diagram_identification = {}
934
935 for i, loop_diag in enumerate(self['loop_diagrams']):
936 loop_tag = loop_diag.build_loop_tag_for_diagram_identification(
937 self['process']['model'], self.get('structure_repository'),
938 use_FDStructure_ID_for_tag = True)
939 # We store the loop diagrams in a 2-tuple that keeps track of 'i'
940 # so that we don't lose their original order. It is just for
941 # convenience, and not strictly necessary.
942 try:
943 diagram_identification[loop_tag].append((i+1,loop_diag))
944 except KeyError:
945 diagram_identification[loop_tag] = [(i+1,loop_diag)]
946
947 # Now sort the loop_tag keys according to their order of appearance
948 sorted_loop_tag_keys = sorted(diagram_identification.keys(),
949 key=lambda k:diagram_identification[k][0][0])
950
951 new_loop_diagram_base = base_objects.DiagramList([])
952 n_loops_identified = 0
953 for loop_tag in sorted_loop_tag_keys:
954 n_diag_in_class = len(diagram_identification[loop_tag])
955 n_loops_identified += n_diag_in_class-1
956 new_loop_diagram_base.append(diagram_identification[loop_tag][0][1])
957 # We must add the counterterms of all the identified loop diagrams
958 # to the reference one.
959 new_loop_diagram_base[-1]['multiplier'] = n_diag_in_class
960 for ldiag in diagram_identification[loop_tag][1:]:
961 new_loop_diagram_base[-1].get('CT_vertices').extend(
962 copy.copy(ldiag[1].get('CT_vertices')))
963 if n_diag_in_class > 1:
964 ldg_debug_info("# Diagram equivalence class detected","#(%s) -> #%d"\
965 %(','.join('%d'%diag[0] for diag in diagram_identification[loop_tag][1:])+
966 (',' if n_diag_in_class==2 else ''),diagram_identification[loop_tag][0][0]))
967
968
969 self.set('loop_diagrams',new_loop_diagram_base)
970 return n_loops_identified
971
973 """This function is solely for monitoring purposes. It reports what are
974 the coupling order combination which are obtained with the diagram
975 genarated and among those which ones correspond to those selected by
976 the process definition and which ones are the extra combinations which
977 comes as a byproduct of the computation of the desired one. The typical
978 example is that if you ask for d d~ > u u~ QCD^2==2 [virt=QCD, QED],
979 you will not only get (QCD,QED)=(2,2);(2,4) which are the desired ones
980 but the code output will in principle also be able to return
981 (QCD,QED)=(4,0);(4,2);(0,4);(0,6) because they involve the same amplitudes
982 """
983
984 hierarchy = self['process']['model']['order_hierarchy']
985
986 sqorders_types=copy.copy(self['process'].get('sqorders_types'))
987 # The WEIGHTED order might have been automatically assigned to the
988 # squared order constraints, so we must assign it a type if not specified
989 if 'WEIGHTED' not in sqorders_types:
990 sqorders_types['WEIGHTED']='<='
991
992 sorted_hierarchy = [order[0] for order in \
993 sorted(hierarchy.items(), key=lambda el: el[1])]
994
995 loop_SOs = set(tuple([d.get_order(order) for order in sorted_hierarchy])
996 for d in self['loop_diagrams']+self['loop_UVCT_diagrams'])
997
998 if self['process']['has_born']:
999 born_SOs = set(tuple([d.get_order(order) for order in \
1000 sorted_hierarchy]) for d in self['born_diagrams'])
1001 else:
1002 born_SOs = set([])
1003
1004 born_sqSOs = set(tuple([x + y for x, y in zip(b1_SO, b2_SO)]) for b1_SO
1005 in born_SOs for b2_SO in born_SOs)
1006 if self['process']['has_born']:
1007 ref_amps = born_SOs
1008 else:
1009 ref_amps = loop_SOs
1010 loop_sqSOs = set(tuple([x + y for x, y in zip(b_SO, l_SO)]) for b_SO in
1011 ref_amps for l_SO in loop_SOs)
1012
1013 # Append the corresponding WEIGHT of each contribution
1014 sorted_hierarchy.append('WEIGHTED')
1015 born_sqSOs = sorted([b_sqso+(sum([b*hierarchy[sorted_hierarchy[i]] for
1016 i, b in enumerate(b_sqso)]),) for b_sqso in born_sqSOs],
1017 key=lambda el: el[1])
1018 loop_sqSOs = sorted([l_sqso+(sum([l*hierarchy[sorted_hierarchy[i]] for
1019 i, l in enumerate(l_sqso)]),) for l_sqso in loop_sqSOs],
1020 key=lambda el: el[1])
1021
1022
1023 logger.debug("Coupling order combinations considered:"+\
1024 " (%s)"%','.join(sorted_hierarchy))
1025
1026 # Now check what is left
1027 born_considered = []
1028 loop_considered = []
1029 for i, sqSOList in enumerate([born_sqSOs,loop_sqSOs]):
1030 considered = []
1031 extra = []
1032 for sqSO in sqSOList:
1033 for sqo, constraint in self['process']['squared_orders'].items():
1034 sqo_index = sorted_hierarchy.index(sqo)
1035 # Notice that I assume here that the negative coupling order
1036 # constraint should have been replaced here (by its
1037 # corresponding positive value).
1038 if (sqorders_types[sqo]=='==' and
1039 sqSO[sqo_index]!=constraint ) or \
1040 (sqorders_types[sqo] in ['=','<='] and
1041 sqSO[sqo_index]>constraint) or \
1042 (sqorders_types[sqo] in ['>'] and
1043 sqSO[sqo_index]<=constraint):
1044 extra.append(sqSO)
1045 break;
1046
1047 # Set the ones considered to be the complement of the omitted ones
1048 considered = [sqSO for sqSO in sqSOList if sqSO not in extra]
1049
1050 if i==0:
1051 born_considered = considered
1052 name = "Born"
1053 if not self['process']['has_born']:
1054 logger.debug(" > No Born contributions for this process.")
1055 continue
1056 elif i==1:
1057 loop_considered = considered
1058 name = "loop"
1059
1060 if len(considered)==0:
1061 logger.debug(" > %s : None"%name)
1062 else:
1063 logger.debug(" > %s : %s"%(name,' '.join(['(%s,W%d)'%(
1064 ','.join(list('%d'%s for s in c[:-1])),c[-1])
1065 for c in considered])))
1066
1067 if len(extra)!=0:
1068 logger.debug(" > %s (not selected but available): %s"%(name,' '.
1069 join(['(%s,W%d)'%(','.join(list('%d'%s for s in e[:-1])),
1070 e[-1]) for e in extra])))
1071
1072 # In case it is needed, the considered orders are returned
1073 # (it is used by some of the unit tests)
1074 return (born_considered,
1075 [sqSO for sqSO in born_sqSOs if sqSO not in born_considered],
1076 loop_considered,
1077 [sqSO for sqSO in loop_sqSOs if sqSO not in loop_considered])
1078
1079
1081 """ Generates all born diagrams relevant to this NLO Process """
1082
1083 bornsuccessful, self['born_diagrams'] = \
1084 diagram_generation.Amplitude.generate_diagrams(self,True)
1085
1086 return bornsuccessful
1087
1089 """ Generates all loop diagrams relevant to this NLO Process """
1090
1091 # Reinitialize the loop diagram container
1092 self['loop_diagrams']=base_objects.DiagramList()
1093 totloopsuccessful=False
1094
1095 # Make sure to start with an empty l-cut particle list.
1096 self.lcutpartemployed=[]
1097
1098 for order in self['process']['perturbation_couplings']:
1099 ldg_debug_info("Perturbation coupling generated now ",order)
1100 lcutPart=[particle for particle in \
1101 self['process']['model']['particles'] if \
1102 (particle.is_perturbating(order, self['process']['model']) and \
1103 particle.get_pdg_code() not in \
1104 self['process']['forbidden_particles'])]
1105 # lcutPart = [lp for lp in lcutPart if abs(lp.get('pdg_code'))==6]
1106 # misc.sprint("lcutPart=",[part.get('name') for part in lcutPart])
1107 for part in lcutPart:
1108 if part.get_pdg_code() not in self.lcutpartemployed:
1109 # First create the two L-cut particles to add to the process.
1110 # Remember that in the model only the particles should be
1111 # tagged as contributing to the a perturbation. Never the
1112 # anti-particle. We chose here a specific orientation for
1113 # the loop momentum flow, say going IN lcutone and OUT
1114 # lcuttwo. We also define here the 'positive' loop fermion
1115 # flow by always setting lcutone to be a particle and
1116 # lcuttwo the corresponding anti-particle.
1117 ldg_debug_info("Generating loop diagram with L-cut type",\
1118 part.get_name())
1119 lcutone=base_objects.Leg({'id': part.get_pdg_code(),
1120 'state': True,
1121 'loop_line': True})
1122 lcuttwo=base_objects.Leg({'id': part.get_anti_pdg_code(),
1123 'state': True,
1124 'loop_line': True})
1125 self['process'].get('legs').extend([lcutone,lcuttwo])
1126 # WARNING, it is important for the tagging to notice here
1127 # that lcuttwo is the last leg in the process list of legs
1128 # and will therefore carry the highest 'number' attribute as
1129 # required to insure that it will never be 'propagated' to
1130 # any output leg.
1131
1132 # We generate the diagrams now
1133 loopsuccessful, lcutdiaglist = \
1134 super(LoopAmplitude, self).generate_diagrams(True)
1135
1136 # Now get rid of all the previously defined l-cut particles.
1137 leg_to_remove=[leg for leg in self['process']['legs'] \
1138 if leg['loop_line']]
1139 for leg in leg_to_remove:
1140 self['process']['legs'].remove(leg)
1141
1142 # The correct L-cut type is specified
1143 for diag in lcutdiaglist:
1144 diag.set('type',part.get_pdg_code())
1145 self['loop_diagrams']+=lcutdiaglist
1146
1147 # Update the list of already employed L-cut particles such
1148 # that we never use them again in loop particles
1149 self.lcutpartemployed.append(part.get_pdg_code())
1150 self.lcutpartemployed.append(part.get_anti_pdg_code())
1151
1152 ldg_debug_info("#Diags generated w/ this L-cut particle",\
1153 len(lcutdiaglist))
1154 # Accordingly update the totloopsuccessful tag
1155 if loopsuccessful:
1156 totloopsuccessful=True
1157
1158 # Reset the l-cut particle list
1159 self.lcutpartemployed=[]
1160
1161 return totloopsuccessful
1162
1163
1165 """ Scan all born diagrams and add for each all the corresponding UV
1166 counterterms. It creates one LoopUVCTDiagram per born diagram and set
1167 of possible coupling_order (so that QCD and QED wavefunction corrections
1168 are not in the same LoopUVCTDiagram for example). Notice that this takes
1169 care only of the UV counterterm which factorize with the born and the
1170 other contributions like the UV mass renormalization are added in the
1171 function setLoopCTVertices"""
1172
1173 # return True
1174 # ============================================
1175 # Including the UVtree contributions
1176 # ============================================
1177
1178 # The following lists the UV interactions potentially giving UV counterterms
1179 # (The UVmass interactions is accounted for like the R2s)
1180 UVCTvertex_interactions = base_objects.InteractionList()
1181 for inter in self['process']['model']['interactions'].get_UV():
1182 if inter.is_UVtree() and len(inter['particles'])>1 and \
1183 inter.is_perturbating(self['process']['perturbation_couplings']) \
1184 and (set(inter['orders'].keys()).intersection(\
1185 set(self['process']['perturbation_couplings'])))!=set([]) and \
1186 (any([set(loop_parts).intersection(set(self['process']\
1187 ['forbidden_particles']))==set([]) for loop_parts in \
1188 inter.get('loop_particles')]) or \
1189 inter.get('loop_particles')==[[]]):
1190 UVCTvertex_interactions.append(inter)
1191
1192 # Temporarly give the tagging order 'UVCT_SPECIAL' to those interactions
1193 self['process']['model'].get('order_hierarchy')['UVCT_SPECIAL']=0
1194 self['process']['model'].get('coupling_orders').add('UVCT_SPECIAL')
1195 for inter in UVCTvertex_interactions:
1196 neworders=copy.copy(inter.get('orders'))
1197 neworders['UVCT_SPECIAL']=1
1198 inter.set('orders',neworders)
1199 # Refresh the model interaction dictionary while including those special
1200 # interactions
1201 self['process']['model'].actualize_dictionaries(useUVCT=True)
1202
1203 # Generate the UVCTdiagrams (born diagrams with 'UVCT_SPECIAL'=0 order
1204 # will be generated along)
1205 self['process']['orders']['UVCT_SPECIAL']=1
1206
1207 UVCTsuccessful, UVCTdiagrams = \
1208 super(LoopAmplitude, self).generate_diagrams(True)
1209
1210 for UVCTdiag in UVCTdiagrams:
1211 if UVCTdiag.get_order('UVCT_SPECIAL')==1:
1212 newUVCTDiag = loop_base_objects.LoopUVCTDiagram({\
1213 'vertices':copy.deepcopy(UVCTdiag['vertices'])})
1214 UVCTinter = newUVCTDiag.get_UVCTinteraction(self['process']['model'])
1215 newUVCTDiag.set('type',UVCTinter.get('type'))
1216 # This interaction counter-term must be accounted for as many times
1217 # as they are list of loop_particles defined and allowed for by
1218 # the process.
1219 newUVCTDiag.get('UVCT_couplings').append((len([1 for loop_parts \
1220 in UVCTinter.get('loop_particles') if set(loop_parts).intersection(\
1221 set(self['process']['forbidden_particles']))==set([])])) if
1222 loop_parts!=[[]] else 1)
1223 self['loop_UVCT_diagrams'].append(newUVCTDiag)
1224
1225 # Remove the additional order requirement in the born orders for this
1226 # process
1227 del self['process']['orders']['UVCT_SPECIAL']
1228 # Remove the fake order added to the selected UVCT interactions
1229 del self['process']['model'].get('order_hierarchy')['UVCT_SPECIAL']
1230 self['process']['model'].get('coupling_orders').remove('UVCT_SPECIAL')
1231 for inter in UVCTvertex_interactions:
1232 del inter.get('orders')['UVCT_SPECIAL']
1233 # Revert the model interaction dictionaries to default
1234 self['process']['model'].actualize_dictionaries(useUVCT=False)
1235
1236 # Set the correct orders to the loop_UVCT_diagrams
1237 for UVCTdiag in self['loop_UVCT_diagrams']:
1238 UVCTdiag.calculate_orders(self['process']['model'])
1239
1240 # ============================================
1241 # Wavefunction renormalization
1242 # ============================================
1243
1244 if not self['process']['has_born']:
1245 return UVCTsuccessful
1246
1247 # We now scan each born diagram, adding the necessary wavefunction
1248 # renormalizations
1249 for bornDiag in self['born_diagrams']:
1250 # This dictionary takes for keys the tuple
1251 # (('OrderName1',power1),...,('OrderNameN',powerN) representing
1252 # the power brought by the counterterm and the value is the
1253 # corresponding LoopUVCTDiagram.
1254 # The last entry is of the form ('EpsilonOrder', value) to put the
1255 # contribution of each different EpsilonOrder to different
1256 # LoopUVCTDiagrams.
1257 LoopUVCTDiagramsAdded={}
1258 for leg in self['process']['legs']:
1259 counterterm=self['process']['model'].get_particle(abs(leg['id'])).\
1260 get('counterterm')
1261 for key, value in counterterm.items():
1262 if key[0] in self['process']['perturbation_couplings']:
1263 for laurentOrder, CTCoupling in value.items():
1264 # Create the order key of the UV counterterm
1265 orderKey=[(key[0],2),]
1266 orderKey.sort()
1267 orderKey.append(('EpsilonOrder',-laurentOrder))
1268 CTCouplings=[CTCoupling for loop_parts in key[1] if
1269 set(loop_parts).intersection(set(self['process']\
1270 ['forbidden_particles']))==set([])]
1271 if CTCouplings!=[]:
1272 try:
1273 LoopUVCTDiagramsAdded[tuple(orderKey)].get(\
1274 'UVCT_couplings').extend(CTCouplings)
1275 except KeyError:
1276 LoopUVCTDiagramsAdded[tuple(orderKey)]=\
1277 loop_base_objects.LoopUVCTDiagram({\
1278 'vertices':copy.deepcopy(bornDiag['vertices']),
1279 'type':'UV'+('' if laurentOrder==0 else
1280 str(-laurentOrder)+'eps'),
1281 'UVCT_orders':{key[0]:2},
1282 'UVCT_couplings':CTCouplings})
1283
1284 for LoopUVCTDiagram in LoopUVCTDiagramsAdded.values():
1285 LoopUVCTDiagram.calculate_orders(self['process']['model'])
1286 self['loop_UVCT_diagrams'].append(LoopUVCTDiagram)
1287
1288 return UVCTsuccessful
1289
1291 """ Scan each loop diagram and recognizes what are the R2/UVmass
1292 CounterTerms associated to them """
1293 #return # debug
1294 # We first create a base dictionary with as a key (tupleA,tupleB). For
1295 # each R2/UV interaction, tuple B is the ordered tuple of the loop
1296 # particles (not anti-particles, so that the PDG is always positive!)
1297 # listed in its loop_particles attribute. Tuple A is the ordered tuple
1298 # of external particles PDGs. making up this interaction. The values of
1299 # the dictionary are a list of the interaction ID having the same key
1300 # above.
1301 CT_interactions = {}
1302 for inter in self['process']['model']['interactions']:
1303 if inter.is_UVmass() or inter.is_UVloop() or inter.is_R2() and \
1304 len(inter['particles'])>1 and inter.is_perturbating(\
1305 self['process']['perturbation_couplings']):
1306 # This interaction might have several possible loop particles
1307 # yielding the same CT. So we add this interaction ID
1308 # for each entry in the list loop_particles.
1309 for i, lparts in enumerate(inter['loop_particles']):
1310 keya=copy.copy(lparts)
1311 keya.sort()
1312 if inter.is_UVloop():
1313 # If it is a CT of type UVloop, then do not specify the
1314 # keya (leave it empty) but make sure the particles
1315 # specified as loop particles are not forbidden before
1316 # adding this CT to CT_interactions
1317 if (set(self['process']['forbidden_particles']) & \
1318 set(lparts)) != set([]):
1319 continue
1320 else:
1321 keya=[]
1322 keyb=[part.get_pdg_code() for part in inter['particles']]
1323 keyb.sort()
1324 key=(tuple(keyb),tuple(keya))
1325 # We keep track of 'i' (i.e. the position of the
1326 # loop_particle list in the inter['loop_particles']) so
1327 # that each coupling in a vertex of type 'UVloop' is
1328 # correctly accounted for since the keya is always replaced
1329 # by an empty list since the constraint on the loop particles
1330 # is simply that there is not corresponding forbidden
1331 # particles in the process definition and not that the
1332 # actual particle content of the loop generate matches.
1333 #
1334 # This can also happen with the type 'UVmass' or 'R2'
1335 # CTvertex ex1(
1336 # type='UVmass'
1337 # [...]
1338 # loop_particles=[[[d,g],[d,g]]])
1339 # Which is a bit silly but can happen and would mean that
1340 # we must account twice for the coupling associated to each
1341 # of these loop_particles.
1342 # One might imagine someone doing it with
1343 # loop_particles=[[[],[]]], for example, because he wanted
1344 # to get rid of the loop particle constraint for some reason.
1345 try:
1346 CT_interactions[key].append((inter['id'],i))
1347 except KeyError:
1348 CT_interactions[key]=[(inter['id'],i),]
1349
1350 # The dictionary CTmass_added keeps track of what are the CounterTerms of
1351 # type UVmass or R2 already added and prevents us from adding them again.
1352 # For instance, the fermion boxes with four external gluons exists in 6 copies
1353 # (with different crossings of the external legs each time) and the
1354 # corresponding R2 must be added only once. The key of this dictionary
1355 # characterizing the loop is (tupleA,tupleB). Tuple A is made from the
1356 # list of the ID of the external structures attached to this loop and
1357 # tuple B from list of the pdg of the particles building this loop.
1358
1359 # Notice that when a CT of type UVmass is specified with an empty
1360 # loop_particles attribute, then it means it must be added once for each
1361 # particle with a matching topology, irrespectively of the loop content.
1362 # Whenever added, such a CT is put in the dictionary CT_added with a key
1363 # having an empty tupleB.
1364 # Finally, because CT interactions of type UVloop do specify a
1365 # loop_particles attribute, but which serves only to be filtered against
1366 # particles forbidden in the process definition, they will also be added
1367 # with an empty tupleB.
1368 CT_added = {}
1369
1370 for diag in self['loop_diagrams']:
1371 # First build the key from this loop for the CT_interaction dictionary
1372 # (Searching Key) and the key for the CT_added dictionary (tracking Key)
1373 searchingKeyA=[]
1374 # Notice that searchingKeyB below also serves as trackingKeyB
1375 searchingKeyB=[]
1376 trackingKeyA=[]
1377 for tagElement in diag['canonical_tag']:
1378 for structID in tagElement[1]:
1379 trackingKeyA.append(structID)
1380 searchingKeyA.append(self['process']['model'].get_particle(\
1381 self['structure_repository'][structID]['binding_leg']['id']).\
1382 get_pdg_code())
1383 searchingKeyB.append(self['process']['model'].get_particle(\
1384 tagElement[0]).get('pdg_code'))
1385 searchingKeyA.sort()
1386 # We do not repeat particles present many times in the loop
1387 searchingKeyB=list(set(searchingKeyB))
1388 searchingKeyB.sort()
1389 trackingKeyA.sort()
1390 # I repeat, they are two kinds of keys:
1391 # searchingKey:
1392 # This serves to scan the CT interactions defined and then find
1393 # which ones match a given loop topology and particle.
1394 # trackingKey:
1395 # Once some CT vertices are identified to be a match for a loop,
1396 # the trackingKey is used in conjunction with the dictionary
1397 # CT_added to make sure that this CT has not already been included.
1398
1399 # Each of these two keys above, has the format
1400 # (tupleA, tupleB)
1401 # with tupleB being the loop_content and either contains the set of
1402 # loop particles PDGs of the interaction (for the searchingKey)
1403 # or of the loops already scanned (trackingKey). It can also be
1404 # empty when considering interactions of type UVmass or R2 which
1405 # have an empty loop_particle attribute or those of type UVloop.
1406 # TupleA is the set of external particle PDG (for the searchingKey)
1407 # and the unordered list of structID attached to the loop (for the
1408 # trackingKey)
1409 searchingKeySimple=(tuple(searchingKeyA),())
1410 searchingKeyLoopPart=(tuple(searchingKeyA),tuple(searchingKeyB))
1411 trackingKeySimple=(tuple(trackingKeyA),())
1412 trackingKeyLoopPart=(tuple(trackingKeyA),tuple(searchingKeyB))
1413 # Now we look for a CT which might correspond to this loop by looking
1414 # for its searchingKey in CT_interactions
1415
1416 # misc.sprint("I have the following CT_interactions=",CT_interactions)
1417 try:
1418 CTIDs=copy.copy(CT_interactions[searchingKeySimple])
1419 except KeyError:
1420 CTIDs=[]
1421 try:
1422 CTIDs.extend(copy.copy(CT_interactions[searchingKeyLoopPart]))
1423 except KeyError:
1424 pass
1425 if not CTIDs:
1426 continue
1427 # We have found some CT interactions corresponding to this loop
1428 # so we must make sure we have not included them already
1429 try:
1430 usedIDs=copy.copy(CT_added[trackingKeySimple])
1431 except KeyError:
1432 usedIDs=[]
1433 try:
1434 usedIDs.extend(copy.copy(CT_added[trackingKeyLoopPart]))
1435 except KeyError:
1436 pass
1437
1438 for CTID in CTIDs:
1439 # Make sure it has not been considered yet and that the loop
1440 # orders match
1441 if CTID not in usedIDs and diag.get_loop_orders(\
1442 self['process']['model'])==\
1443 self['process']['model']['interaction_dict'][CTID[0]]['orders']:
1444 # Create the amplitude vertex corresponding to this CT
1445 # and add it to the LoopDiagram treated.
1446 CTleglist = base_objects.LegList()
1447 for tagElement in diag['canonical_tag']:
1448 for structID in tagElement[1]:
1449 CTleglist.append(\
1450 self['structure_repository'][structID]['binding_leg'])
1451 CTVertex = base_objects.Vertex({'id':CTID[0], \
1452 'legs':CTleglist})
1453 diag['CT_vertices'].append(CTVertex)
1454 # Now add this CT vertex to the CT_added dictionary so that
1455 # we are sure it will not be double counted
1456 if self['process']['model']['interaction_dict'][CTID[0]]\
1457 ['loop_particles'][CTID[1]]==[] or \
1458 self['process']['model']['interaction_dict'][CTID[0]].\
1459 is_UVloop():
1460 try:
1461 CT_added[trackingKeySimple].append(CTID)
1462 except KeyError:
1463 CT_added[trackingKeySimple] = [CTID, ]
1464 else:
1465 try:
1466 CT_added[trackingKeyLoopPart].append(CTID)
1467 except KeyError:
1468 CT_added[trackingKeyLoopPart] = [CTID, ]
1469
1471 """ Return a LoopDiagram created."""
1472 return loop_base_objects.LoopDiagram({'vertices':vertexlist})
1473
1475 """ Returns a DGLoopLeg list instead of the default copy_leglist
1476 defined in base_objects.Amplitude """
1477
1478 dgloopleglist=base_objects.LegList()
1479 for leg in leglist:
1480 dgloopleglist.append(loop_base_objects.DGLoopLeg(leg))
1481
1482 return dgloopleglist
1483
1485 """ Overloaded here to convert back all DGLoopLegs into Legs. """
1486 for vertexlist in vertexdoublelist:
1487 for vertex in vertexlist:
1488 if not isinstance(vertex['legs'][0],loop_base_objects.DGLoopLeg):
1489 continue
1490 vertex['legs'][:]=[leg.convert_to_leg() for leg in \
1491 vertex['legs']]
1492 return True
1493
1495 """Create a set of new legs from the info given."""
1496
1497 looplegs=[leg for leg in legs if leg['loop_line']]
1498
1499 # Get rid of all vanishing tadpoles
1500 #Ease the access to the model
1501 model=self['process']['model']
1502 exlegs=[leg for leg in looplegs if leg['depth']==0]
1503 if(len(exlegs)==2):
1504 if(any([part['mass'].lower()=='zero' for pdg,part in model.get('particle_dict').items() if pdg==abs(exlegs[0]['id'])])):
1505 return []
1506
1507 # Correctly propagate the loopflow
1508 loopline=(len(looplegs)==1)
1509 mylegs = []
1510 for i, (leg_id, vert_id) in enumerate(leg_vert_ids):
1511 # We can now create the set of possible merged legs.
1512 # However, we make sure that its PDG is not in the list of
1513 # L-cut particles we already explored. If it is, we simply reject
1514 # the diagram.
1515 if not loopline or not (leg_id in self.lcutpartemployed):
1516 # Reminder: The only purpose of the "depth" flag is to get rid
1517 # of (some, not all) of the wave-function renormalization
1518 # already during diagram generation. We reckognize a wf
1519 # renormalization diagram as follows:
1520 if len(legs)==2 and len(looplegs)==2:
1521 # We have candidate
1522 depths=(looplegs[0]['depth'],looplegs[1]['depth'])
1523 if (0 in depths) and (-1 not in depths) and depths!=(0,0):
1524 # Check that the PDG of the outter particle in the
1525 # wavefunction renormalization bubble is equal to the
1526 # one of the inner particle.
1527 continue
1528
1529 # If depth is not 0 because of being an external leg and not
1530 # the propagated PDG, then we set it to -1 so that from that
1531 # point we are sure the diagram will not be reckognized as a
1532 # wave-function renormalization.
1533 depth=-1
1534 # When creating a loop leg from exactly two external legs, we
1535 # set the depth to the PDG of the external non-loop line.
1536 if len(legs)==2 and loopline and (legs[0]['depth'],\
1537 legs[1]['depth'])==(0,0):
1538 if not legs[0]['loop_line']:
1539 depth=legs[0]['id']
1540 else:
1541 depth=legs[1]['id']
1542 # In case of two point interactions among two same particle
1543 # we propagate the existing depth
1544 if len(legs)==1 and legs[0]['id']==leg_id:
1545 depth=legs[0]['depth']
1546 # In all other cases we set the depth to -1 since no
1547 # wave-function renormalization diagram can arise from this
1548 # side of the diagram construction.
1549
1550 mylegs.append((loop_base_objects.DGLoopLeg({'id':leg_id,
1551 'number':number,
1552 'state':state,
1553 'from_group':True,
1554 'depth': depth,
1555 'loop_line': loopline}),
1556 vert_id))
1557 return mylegs
1558
1560 """Allow for selection of vertex ids."""
1561
1562 looplegs=[leg for leg in legs if leg['loop_line']]
1563 nonlooplegs=[leg for leg in legs if not leg['loop_line']]
1564
1565 # Get rid of all vanishing tadpoles
1566 model=self['process']['model']
1567 exlegs=[leg for leg in looplegs if leg['depth']==0]
1568 if(len(exlegs)==2):
1569 if(any([part['mass'].lower()=='zero' for pdg,part in \
1570 model.get('particle_dict').items() if pdg==abs(exlegs[0]['id'])])):
1571 return []
1572
1573
1574 # Get rid of some wave-function renormalization diagrams already during
1575 # diagram generation already.In a similar manner as in get_combined_legs.
1576 if(len(legs)==3 and len(looplegs)==2):
1577 depths=(looplegs[0]['depth'],looplegs[1]['depth'])
1578 if (0 in depths) and (-1 not in depths) and depths!=(0,0):
1579 return []
1580
1581 return vert_ids
1582
1583 # Helper function
1584
1586 """ Filters the diagrams according to the constraints on the squared
1587 orders in argument and wether the process has a born or not. """
1588
1589 diagRef=base_objects.DiagramList()
1590 AllLoopDiagrams=base_objects.DiagramList(self['loop_diagrams']+\
1591 self['loop_UVCT_diagrams'])
1592
1593 AllBornDiagrams=base_objects.DiagramList(self['born_diagrams'])
1594 if self['process']['has_born']:
1595 diagRef=AllBornDiagrams
1596 else:
1597 diagRef=AllLoopDiagrams
1598
1599 sqorders_types=copy.copy(self['process'].get('sqorders_types'))
1600
1601 # The WEIGHTED order might have been automatically assigned to the
1602 # squared order constraints, so we must assign it a type if not specified
1603 if 'WEIGHTED' not in sqorders_types:
1604 sqorders_types['WEIGHTED']='<='
1605
1606 if len(diagRef)==0:
1607 # If no born contributes but they were supposed to ( in the
1608 # case of self['process']['has_born']=True) then it means that
1609 # the loop cannot be squared against anything and none should
1610 # contribute either. The squared order constraints are just too
1611 # tight for anything to contribute.
1612 AllLoopDiagrams = base_objects.DiagramList()
1613
1614
1615 # Start by filtering the loop diagrams
1616 AllLoopDiagrams = AllLoopDiagrams.apply_positive_sq_orders(diagRef,
1617 sq_order_constrains, sqorders_types)
1618 # And now the Born ones if there are any
1619 if self['process']['has_born']:
1620 # We consider both the Born*Born and Born*Loop squared terms here
1621 AllBornDiagrams = AllBornDiagrams.apply_positive_sq_orders(
1622 AllLoopDiagrams+AllBornDiagrams, sq_order_constrains, sqorders_types)
1623
1624 # Now treat the negative squared order constraint (at most one)
1625 neg_orders = [(order, value) for order, value in \
1626 sq_order_constrains.items() if value<0]
1627 if len(neg_orders)==1:
1628 neg_order, neg_value = neg_orders[0]
1629 # If there is a Born contribution, then the target order will
1630 # be computed over all Born*Born and Born*loop contributions
1631 if self['process']['has_born']:
1632 AllBornDiagrams, target_order =\
1633 AllBornDiagrams.apply_negative_sq_order(
1634 base_objects.DiagramList(AllLoopDiagrams+AllBornDiagrams),
1635 neg_order,neg_value,sqorders_types[neg_order])
1636 # Now we must filter the loop diagrams using to the target_order
1637 # computed above from the LO and NLO contributions
1638 AllLoopDiagrams = AllLoopDiagrams.apply_positive_sq_orders(
1639 diagRef,{neg_order:target_order},
1640 {neg_order:sqorders_types[neg_order]})
1641
1642 # If there is no Born, then the situation is completely analoguous
1643 # to the tree level case since it is simply Loop*Loop
1644 else:
1645 AllLoopDiagrams, target_order = \
1646 AllLoopDiagrams.apply_negative_sq_order(
1647 diagRef,neg_order,neg_value,sqorders_types[neg_order])
1648
1649 # Substitute the negative value to this positive one
1650 # (also in the backed up values in user_squared_orders so that
1651 # this change is permanent and we will still have access to
1652 # it at the output stage)
1653 self['process']['squared_orders'][neg_order]=target_order
1654 user_squared_orders[neg_order]=target_order
1655
1656 elif len(neg_orders)>1:
1657 raise MadGraph5Error('At most one negative squared order constraint'+\
1658 ' can be specified, not %s.'%str(neg_orders))
1659
1660 if self['process']['has_born']:
1661 self['born_diagrams'] = AllBornDiagrams
1662 self['loop_diagrams']=[diag for diag in AllLoopDiagrams if not \
1663 isinstance(diag,loop_base_objects.LoopUVCTDiagram)]
1664 self['loop_UVCT_diagrams']=[diag for diag in AllLoopDiagrams if \
1665 isinstance(diag,loop_base_objects.LoopUVCTDiagram)]
1666
1668 """ This is a helper function for order_diagrams_according_to_split_orders
1669 and intended to be used from LoopHelasAmplitude only"""
1670
1671 # The dictionary below has keys being the tuple (split_order<i>_values)
1672 # and values being diagram lists sharing the same split orders.
1673 diag_by_so = {}
1674
1675 for diag in diag_set:
1676 so_key = tuple([diag.get_order(order) for order in split_orders])
1677 try:
1678 diag_by_so[so_key].append(diag)
1679 except KeyError:
1680 diag_by_so[so_key]=base_objects.DiagramList([diag,])
1681
1682 so_keys = diag_by_so.keys()
1683 # Complete the order hierarchy by possibly missing defined order for
1684 # which we set the weight to zero by default (so that they are ignored).
1685 order_hierarchy = self.get('process').get('model').get('order_hierarchy')
1686 order_weights = copy.copy(order_hierarchy)
1687 for so in split_orders:
1688 if so not in order_hierarchy.keys():
1689 order_weights[so]=0
1690
1691 # Now order the keys of diag_by_so by the WEIGHT of the split_orders
1692 # (and only those, the orders not included in the split_orders do not
1693 # count for this ordering as they could be mixed in any given group).
1694 so_keys = sorted(so_keys, key = lambda elem: (sum([power*order_weights[\
1695 split_orders[i]] for i,power in enumerate(elem)])))
1696
1697 # Now put the diagram back, ordered this time, in diag_set
1698 diag_set[:] = []
1699 for so_key in so_keys:
1700 diag_set.extend(diag_by_so[so_key])
1701
1702
1704 """ Reorder the loop and Born diagrams (if any) in group of diagrams
1705 sharing the same coupling orders are put together and these groups are
1706 order in decreasing WEIGHTED orders.
1707 Notice that this function is only called for now by the
1708 LoopHelasMatrixElement instances at the output stage.
1709 """
1710
1711 # If no split order is present (unlikely since the 'corrected order'
1712 # normally is a split_order by default, then do nothing
1713 if len(split_orders)==0:
1714 return
1715
1716 self.order_diagram_set(self['born_diagrams'], split_orders)
1717 self.order_diagram_set(self['loop_diagrams'], split_orders)
1718 self.order_diagram_set(self['loop_UVCT_diagrams'], split_orders)
1719
1720 #===============================================================================
1721 # LoopMultiProcess
1722 #===============================================================================
1723 -class LoopMultiProcess(diagram_generation.MultiProcess):
1724 """LoopMultiProcess: MultiProcess with loop features.
1725 """
1726
1727 @classmethod
1729 """ Return the correct amplitude type according to the characteristics
1730 of the process proc """
1731 return LoopAmplitude({"process": proc},**opts)
1732
1733 #===============================================================================
1734 # LoopInducedMultiProcess
1735 #===============================================================================
1736 -class LoopInducedMultiProcess(diagram_generation.MultiProcess):
1737 """Special mode for the LoopInduced."""
1738
1739 @classmethod
1741 """ Return the correct amplitude type according to the characteristics of
1742 the process proc """
1743 return LoopAmplitude({"process": proc, 'has_born':False},**opts)
1744
| Trees | Indices | Help |
|---|
| Generated by Epydoc 3.0.1 on Fri May 26 09:50:18 2017 | http://epydoc.sourceforge.net |