Merge lp:~maddevelopers/mg5amcnlo/PY8_parallelization into lp:~maddevelopers/mg5amcnlo/2.5.1

[Valentin Hirschi]

This branch implements PY8 LO parallelization.

It is fully functional and has been tested for a few cases (not exhaustive yet).

The remaining issues are:

1) [Not related to this branch] Systematics parallelization crash when using a cluster_temp_directory.

2) [Not related to this branch] The PY8 HTML is screwed up past the first run/tag.

3) The .lhe splitting is slow. It would be nice to have an advanced function for this in lhe_parser.py that bypasses the (full) parsing of the event files.
Also, obtaining the number of events in the event file is slow. Again an optimized static method for this bypassing the full parsing would be nice.

4) The merging of the split HEPMC is done in a very efficient way in this branch (very important given their size). However it use two system calls which are not secure. They need to be made secure.

5) More testing must be made, especially a comparison of the results between parallel and sequential runs for the merged_x_secs, HwU plots et hepmc event files must be performed so as to guarantee the correctness of the implementation.

6) The new bits of code in do_pythia8() could be a bit more refactored. In particular the two parts of the code related to parallel submission and merging of split results could be factored out in dedicated functions.

Olivier, could you review this and fix what you can already.
I you manage to clean it up all, then don't hesitate to merge this already to 2.5.1 (or even 2.5.0 since it is a nice feature and introduces some important bug fixing).
If there is still something that needs to be discussed with me, then we will only be able to do this on tuesday, since I will be mostly unavailable 'til then now.

Thanks,

[Valentin Hirschi (valentin-hirschi)]

To test the above on a condor cluster, one must re-install MG5aMC_PY8_interface, because I modified its installation so that it links *statically* against HEPMC2 so that it doesn't have to be found on the worker nodes at run time.

[Olivier Mattelaer (olivier-mattelaer)]

Hi Valentin,

1) I do not think that I will fix that.

2) this is now fixed in 2.5.0

3) I have improved this point. However, you should avoid to use the len() function.
By bypassing the parsing of the event, you are indeed much faster from 50s to s

4) we can discuss that.

5) I only did one simple test:
but it did not work:

Pythia8 shower jobs: 0 Idle, 1 Running, 7 Done [49m33s]
Pythia8 shower jobs: 0 Idle, 0 Running, 8 Done [51m10s]
Merging results from the split PY8 runs...
Fail to produce a pythia8 output. More info in
     /Users/omatt/Documents/eclipse/PY8_parallelization/PROC_sm_1/Events/run_01/tag_1_pythia8.log

The timing seems to indicate that the splitting did not occur at all...

Cheers,

Olivier

review: Needs Fixing

[lp:~maddevelopers/mg5amcnlo/PY8_parallelization]

301. By Olivier Mattelaer on 2016-09-05

faster parsing for splitting event/get number of events

302. By Olivier Mattelaer on 2016-09-06

also apply the bypass of parsing for systematics

303. By Valentin Hirschi on 2016-09-08

1. Fixed the sanity check of PY8 log file which was not ok with the parallelization.

304. By Valentin Hirschi on 2016-09-08

1. fixed an issue with the warning about failing PY8 log. (needed to close the log stream).
2. Sandboxed the HEPMC merging syscalls.

305. By Valentin Hirschi on 2016-09-08

1. Merged with latest version of 2.5.1

[Valentin Hirschi (valentin-hirschi)]

Olivier,

I addressed your remaining points. Hepmc concatenation is now secure, although the event number are no longer unique nor sequential; but this appears to be irrelevant and can be easily fixed if necessary in the future (I'd rather not if not necessary as I want to avoid parsing the HEPMC).

Concerning timing, I tested 10k events for p p > t t~ with the 'simplepy8' shower and '/dev/null' for the HEPMC output. The outcome is (for the shower alone):

2.5.0, single core:

INFO: Pythia8 shower finished after 2m15s.

PY8_parallelization branch on 6 cores:

INFO: Pythia8 shower finished after 57 seconds.

So this is not linear because of the remaining initialization overhead, but it clearly shows that parallelization occurs.

This branch should therefore be ready to be merged in 2.5.1.

=== modified file 'madgraph/interface/common_run_interface.py'
--- madgraph/interface/common_run_interface.py	2016-09-08 09:15:22 +0000
+++ madgraph/interface/common_run_interface.py	2016-09-08 21:48:50 +0000
@@ -3896,7 +3896,7 @@
     'lhc'   : 'syntax: set lhc VALUE:\n      Set for a proton-proton collision with that given center of mass energy (in TeV)',
     'lep'   : 'syntax: set lep VALUE:\n      Set for a electron-positron collision with that given center of mass energy (in GeV)',
     'fixed_scale' : 'syntax: set fixed_scale VALUE:\n      Set all scales to the give value (in GeV)',
-    'simplepy8' : 'Turn off non-perturbative slow features of Pythia8.'
+    'simplepy8' : 'Turn off non-perturbative slow features of Pythia8.',
     'mpi' : 'syntax: set mpi value: allow to turn mpi in Pythia8 on/off'         
     }
     

=== modified file 'madgraph/interface/madevent_interface.py'
--- madgraph/interface/madevent_interface.py	2016-09-08 09:15:22 +0000
+++ madgraph/interface/madevent_interface.py	2016-09-08 21:48:50 +0000
@@ -4024,7 +4024,6 @@
                 
                 # Now merge logs
                 pythia_log_file = open(pythia_log,'w')
-                
                 n_added = 0
                 for split_dir in split_dirs:
                     log_file = pjoin(split_dir,'PY8_log.txt')
@@ -4049,6 +4048,8 @@
                            PY8_extracted_information['Ntry'] = Ntry
                         else:
                            PY8_extracted_information['Ntry'] += Ntry
+                pythia_log_file.close()
+                
                 # Normalize the values added
                 if n_added>0:
                     PY8_extracted_information['sigma_m'] /= float(n_added)               
@@ -4129,7 +4130,7 @@
                                 break
                         header.close()
                         tail = open(pjoin(tmp_dir,'tail.hepmc'),'w')
-                        n_tail = 0               
+                        n_tail = 0
                         for line in misc.BackRead(all_hepmc_files[-1]):
                             if line.startswith('HepMC::'):
                                 n_tail += 1
@@ -4139,16 +4140,16 @@
                         tail.close()
                         if n_tail>1:
                             raise MadGraph5Error,'HEPMC files should only have one trailing command.'
-                        ###################################################################### 
-                        # This is the most efficient way of putting together HEPMC's, *BUT*  #
-                        #    WARNING: NEED TO RENDER THE CODE BELOW SAFE TOWARDS INJECTION   #
-                        ######################################################################
+
                         for hepmc_file in all_hepmc_files:
                             # Remove in an efficient way the starting and trailing HEPMC tags
-                            os.system(' '.join(['sed','-i',"''","'%s;$d'"%
-                                        (';'.join('%id'%(i+1) for i in range(n_head))),hepmc_file]))                            
-                        os.system(' '.join(['cat',pjoin(tmp_dir,'header.hepmc')]+all_hepmc_files+
-                                                    [pjoin(tmp_dir,'tail.hepmc'),'>',hepmc_output]))
+                            misc.call(['sed','-i',"''",
+                                       "%s;$d"%(';'.join('%id'%(i+1) for i in range(n_head))),
+                                       hepmc_file],cwd=parallelization_dir)
+                        misc.call(['cat',pjoin(tmp_dir,'header.hepmc')]+all_hepmc_files+
+                                                [pjoin(tmp_dir,'tail.hepmc')],
+                                                stdout=open(hepmc_output,'w'),
+                                                cwd=os.path.dirname(hepmc_output))
 
                 # We are done with the parallelization directory. Clean it.
                 if os.path.isdir(parallelization_dir):
@@ -4163,9 +4164,9 @@
         if os.path.isfile(pt_output):
             shutil.move(pt_output, pjoin(self.me_dir,'Events',
                                             self.run_name, '%s_pts.dat' % tag))
-            
+
         if not os.path.isfile(pythia_log) or \
-             'PYTHIA Abort' in '\n'.join(open(pythia_log,'r').readlines()[:-20]):
+             'cross section' not in '\n'.join(open(pythia_log,'r').readlines()[-20:]):
             logger.warning('Fail to produce a pythia8 output. More info in \n     %s'%pythia_log)
             return
         

=== modified file 'madgraph/various/lhe_parser.py'
--- madgraph/various/lhe_parser.py	2016-09-08 09:15:22 +0000
+++ madgraph/various/lhe_parser.py	2016-09-08 21:48:50 +0000
@@ -153,6 +153,8 @@
 class EventFile(object):
     """A class to allow to read both gzip and not gzip file"""
 
+    parse_event= True #specify if the parsing of the events has to be done
+
     def __new__(self, path, mode='r', *args, **opt):
         
         if not path.endswith(".gz"):
@@ -183,6 +185,7 @@
                 raise
                 
         self.banner = ''
+
         if mode == 'r':
             line = ''
             while '</init>' not in line.lower():
@@ -227,14 +230,15 @@
         if hasattr(self,"len"):
             return self.len
 
-        init_pos = self.tell()
-        self.seek(0)
-        nb_event=0
-        for _ in self:
-            nb_event +=1
-        self.len = nb_event
-        self.seek(init_pos)
-        return self.len
+        with misc.TMP_variable(self, 'parse_event', False):
+            init_pos = self.tell()
+            self.seek(0)
+            nb_event=0
+            for _ in self:
+                nb_event +=1
+            self.len = nb_event
+            self.seek(init_pos)
+            return self.len
 
     def next(self):
         """get next event"""
@@ -248,8 +252,11 @@
                 text = ''
             if mode:
                 text += line
-        return Event(text)
-    
+        if self.parse_event:
+            return Event(text)
+        else:
+            return text
+        
     def initialize_unweighting(self, get_wgt, trunc_error):
         """ scan once the file to return 
             - the list of the hightest weight (of size trunc_error*NB_EVENT
@@ -516,29 +523,30 @@
     def split(self, nb_event=0, partition=None, cwd=os.path.curdir, zip=False):
         """split the file in multiple file. Do not change the weight!"""
 
-        nb_file = -1
-        for i, event in enumerate(self):
-            if (not (partition is None) and i==sum(partition[:nb_file+1])) or \
-                                   (partition is None and i % nb_event == 0):
-                if i:
-                    #close previous file
-                    current.write('</LesHouchesEvent>\n')
-                    current.close()
-                # create the new file
-                nb_file +=1
-                # If end of partition then finish writing events here.
-                if not partition is None and (nb_file+1>len(partition)):
-                    return nb_file+1
-                if zip:
-                    current = EventFile(pjoin(cwd,'%s_%s.lhe.gz' % (self.name, nb_file)),'w')
-                else:
-                    current = open(pjoin(cwd,'%s_%s.lhe' % (self.name, nb_file)),'w')                    
-                current.write(self.banner)
-            current.write(str(event))
-        if i!=0:
-            current.write('</LesHouchesEvent>\n')
-            current.close()   
-        return nb_file +1
+        with misc.TMP_variable(self, 'parse_event', False):
+            nb_file = -1
+            for i, event in enumerate(self):
+                if (not (partition is None) and i==sum(partition[:nb_file+1])) or \
+                                       (partition is None and i % nb_event == 0):
+                    if i:
+                        #close previous file
+                        current.write('</LesHouchesEvent>\n')
+                        current.close()
+                    # create the new file
+                    nb_file +=1
+                    # If end of partition then finish writing events here.
+                    if not partition is None and (nb_file+1>len(partition)):
+                        return nb_file+1
+                    if zip:
+                        current = EventFile(pjoin(cwd,'%s_%s.lhe.gz' % (self.name, nb_file)),'w')
+                    else:
+                        current = open(pjoin(cwd,'%s_%s.lhe' % (self.name, nb_file)),'w')                    
+                    current.write(self.banner)
+                current.write(str(event))
+            if i!=0:
+                current.write('</LesHouchesEvent>\n')
+                current.close()   
+            return nb_file +1
 
     def update_HwU(self, hwu, fct, name='lhe', keep_wgt=False, maxevents=sys.maxint):
         """take a HwU and add this event file for the function fct"""

=== modified file 'madgraph/various/systematics.py'
--- madgraph/various/systematics.py	2016-09-07 23:55:03 +0000
+++ madgraph/various/systematics.py	2016-09-08 21:48:50 +0000
@@ -218,8 +218,13 @@
         ids = [lowest_id+i for i in range(len(self.args)-1)]
         all_cross = [0 for i in range(len(self.args))]
         
+        if self.start_event !=0:
+            self.input.parse_event = False
+
         for nb_event,event in enumerate(self.input):
             if nb_event < self.start_event:
+                if nb_event == self.start_event-1:
+                    self.input.parse_event = True
                 continue
             elif nb_event >= self.stop_event:
                 if self.force_write_banner: