Merge lp:~maddevelopers/mg5amcnlo/HwU_allocated into lp:~mg5core2/mg5amcnlo/2.5.5

Proposed by Rikkert Frederix
Status: Merged
Merge reported by: Olivier Mattelaer
Merged at revision: not available
Proposed branch: lp:~maddevelopers/mg5amcnlo/HwU_allocated
Merge into: lp:~mg5core2/mg5amcnlo/2.5.5
Diff against target: 449 lines (+202/-83)
3 files modified
Template/NLO/FixedOrderAnalysis/HwU.f (+200/-59)
Template/NLO/FixedOrderAnalysis/HwU.inc (+0/-22)
Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f (+2/-2)
To merge this branch: bzr merge lp:~maddevelopers/mg5amcnlo/HwU_allocated
Reviewer Review Type Date Requested Status
Valentin Hirschi (community) Approve
Review via email: mp+323759@code.launchpad.net

Description of the change

Converted a large HwU common block to a module. This allows for dynamic (de-)allocation of the memory required by all variables in the common block, reducing the footprint and size of the executable in case large analyses are used. Moreover, this leads to more flexible plotting options (more histograms can be 'booked' with as many bins as one requires).

Changes are very much contained to (almost) a single file, hence, review should be straight-forward.

To post a comment you must log in.
295. By Rikkert Frederix

small fix for linux

Revision history for this message
Valentin Hirschi (valentin-hirschi) wrote :

Hi Rik,

Sorry for the late review.
The code seems OK and definitely a welcome addition.
It seems that all arrays are correctly deallocated in HwU_deallocate_all, but I see that it is only called in the "output" subroutine, which is in principle not ideal.
The usual structure is to have a master "clean_up" function in MadFKS that is called at the end of any driver and takes care of calling the deallocation of all possible modules used.
Now in this case, it of course works since we know that the user will output once only and at the end of the run basically.

I think move_alloc is supported in gfortran 4.6 (our min. requirement), but it is always worth checking that it compiles with v4.6 when coding beyond f77.

Also, ideally, one could turn the whole HwU.f into a module (not just its variables) and use it as library, putting it in the Source directory, but it's not crucial (and would also mean adding 'use HwU' in all the existing FO analysis which is not optimal).

Also note that in fortran2003+, one can use "derived types" which are very useful to mimic C++ 'struct' in fortran.
This could be especially useful in cases like this one where one would really want to allocate and de-allocate arrays of histograms directly, which could themselves contain allocatable arrays for their various properties (such as bins).
I had done something similar for the the dimensions of my DiscreteSampler (where bins are themselves derived types):

      type sampledDimension
        integer :: grid_mode
        real*8 :: small_contrib_threshold
        real*8 :: damping_power
        integer :: min_bin_probing_points
        real*8 :: norm
        real*8 :: abs_norm
        real*8 :: variance_norm
        real*8 :: norm_sqr
        integer :: n_tot_entries
        character, dimension(:), allocatable :: dimension_name
        type(bin) , dimension(:), allocatable :: bins
      endtype sampledDimension

Anyway, none of this is essential here, so thanks for this improvement!

review: Approve

Preview Diff

[H/L] Next/Prev Comment, [J/K] Next/Prev File, [N/P] Next/Prev Hunk
1=== modified file 'Template/NLO/FixedOrderAnalysis/HwU.f'
2--- Template/NLO/FixedOrderAnalysis/HwU.f 2017-02-03 13:37:36 +0000
3+++ Template/NLO/FixedOrderAnalysis/HwU.f 2017-05-09 14:02:00 +0000
4@@ -1,7 +1,7 @@
5 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
6 C C
7 C HwU: Histograms with Uncertainties C
8-C By Rikkert Frederix, Dec. 2014 C
9+C By Rikkert Frederix, 12-2014--05-2017 C
10 C C
11 C Book, fill and write out histograms. Particularly suited for NLO C
12 C computations with correlations between points (ie. event and C
13@@ -9,26 +9,41 @@
14 C and PDF uncertainties through reweighting). C
15 C C
16 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
17+
18+c The module contains effectively the common block with allocatable
19+c variables (something not possible in old fortran version)
20+ module HwU_variables
21+ implicit none
22+ integer :: max_plots,max_points,max_bins,nwgts,np
23+ integer :: error_estimation=3
24+ logical, allocatable :: booked(:)
25+ integer, allocatable :: nbin(:),histi(:,:),p_bin(:),p_label(:)
26+ character(len=50), allocatable :: title(:)
27+ character(len=50), allocatable :: wgts_info(:)
28+ double precision, allocatable :: histy(:,:,:),histy_acc(:,:,:)
29+ $ ,histy2(:,:),histy_err(:,:),histxl(:,:),histxm(:,:)
30+ $ ,step(:),p_wgts(:,:)
31+ save
32+ end module HwU_variables
33+
34
35 c To be called once at the start of each run. Initialises the packages
36 c and sets the number of weights that need to be included for each point.
37 subroutine HwU_inithist(nweights,wgt_info)
38+ use HwU_variables
39 implicit none
40- include "HwU.inc"
41 integer i,nweights
42 character*(*) wgt_info(*)
43- do i=1,max_plots
44- booked(i)=.false.
45- enddo
46+ call HwU_deallocate_all
47+ max_plots=0
48+ max_points=0
49+ max_bins=0
50+ np=0
51 c Number of weights associated to each point. Note that the first
52 c weight should always be the 'central value' and it should not be
53 c zero if any of the other weights are non-zero.
54 nwgts=nweights
55- if (nwgts.gt.max_wgts) then
56- write (*,*) 'ERROR: increase max_wgts in HwU histogramming'
57- $ ,max_wgts,nwgts
58- stop 1
59- endif
60+ allocate(wgts_info(nwgts))
61 do i=1,nwgts
62 wgts_info(i)=wgt_info(i)
63 enddo
64@@ -48,10 +63,9 @@
65 c square root of the sum of the squares. Perform a weighted average
66 c iteration-by-iteration
67 c input 3: Same as input 2, but weighted average is same as from MINT
68+ use HwU_variables
69 implicit none
70 integer input
71- integer error_estimation
72- common /HwU_common2/ error_estimation
73 if (input.ge.0 .and. input.le.3) then
74 error_estimation=input
75 else
76@@ -60,39 +74,20 @@
77 endif
78 return
79 end
80-
81- block data HwU
82-c set the default for the error estimation method. To reduce the size of
83-c the executable, put the error_estimation variable in a separate common
84-c block. If we would have included the 'HwU.inc' file here, that
85-c complete common block seems to be included in the size executable
86-c (approx. 110 MB).
87- integer error_estimation
88- common /HwU_common2/ error_estimation
89- data error_estimation /3/
90- end
91-
92+
93 c Book the histograms at the start of the run. Give a 'label' (an
94 c integer) that identifies the plot when filling it and a title
95 c ('title_l') for each plot. Also the number of bins ('nbin_l') and the
96 c plot range (from 'xmin' to 'xmax') should be given.
97 subroutine HwU_book(label,title_l,nbin_l,xmin,xmax)
98+ use HwU_variables
99 implicit none
100- include "HwU.inc"
101 integer label,nbin_l,i,j
102 character*(*) title_l
103 double precision xmin,xmax
104-c Check that label and number of bins are reasonable
105- if (label.gt.max_plots) then
106- write (*,*) 'ERROR: increase max_plots in HwU histogramming'
107- $ ,max_plots, label
108- stop 1
109- endif
110- if (nbin_l.gt.max_bins) then
111- write (*,*) 'ERROR: increase max_bins in HwU histogramming'
112- $ ,max_bins,nbin_l
113- stop 1
114- endif
115+c Allocate space for new histograms if needed
116+ call HwU_allocate_histo(label,nbin_l)
117+c Setup the histogram
118 booked(label)=.true.
119 title(label)=title_l
120 nbin(label)=nbin_l
121@@ -111,7 +106,6 @@
122 histy2(label,i)=0d0
123 histy_err(label,i)=0d0
124 enddo
125- np=0
126 return
127 end
128
129@@ -121,8 +115,8 @@
130 c the 'HwU_inithist' subroutine. That means that each point should have
131 c the same number of weights.
132 subroutine HwU_fill(label,x,wgts)
133+ use HwU_variables
134 implicit none
135- include "HwU.inc"
136 integer label,i,j,bin
137 double precision x, wgts(*)
138 c If central weight is zero do not add this point.
139@@ -147,11 +141,7 @@
140 enddo
141 c If a new bin, add it to the list of points
142 np=np+1
143- if (np.gt.max_points) then
144- write (*,*) 'ERROR: increase max_points in HwU histogramming'
145- $ ,max_points
146- stop 1
147- endif
148+ call HwU_allocate_p
149 p_label(np)=label
150 p_bin(np)=bin
151 do j=1,nwgts
152@@ -159,7 +149,7 @@
153 enddo
154 return
155 end
156-
157+
158 c Call after all correlated contributions for a give phase-space
159 c point. I.e., every time you get a new set of random numbers from
160 c MINT/VEGAS. It adds the current list of points to the histograms. Add
161@@ -168,8 +158,8 @@
162 c this way, correlations between events and counter-events can be
163 c correctly taken into account.
164 subroutine HwU_add_points
165+ use HwU_variables
166 implicit none
167- include "HwU.inc"
168 integer i,j
169 do i=1,np
170 do j=1,nwgts
171@@ -192,12 +182,11 @@
172 c the current iteration so that they can be filled with the next
173 c iteration.
174 subroutine HwU_accum_iter(inclde,nPSpoints,values)
175+ use HwU_variables
176 implicit none
177- include "HwU.inc"
178 logical inclde
179 integer nPSpoints,label,i,j
180- double precision nPSinv,etot,vtot(max_wgts),niter,y_squared
181- $ ,values(2)
182+ double precision nPSinv,etot,niter,y_squared,values(2)
183 data niter /0d0/
184 nPSinv = 1d0/dble(nPSpoints)
185 if (inclde) niter = niter+1d0
186@@ -226,8 +215,8 @@
187 c intermediate stages this function can be called (together with
188 c HwU_output) to write intermediate plots to disk.
189 subroutine finalize_histograms(nPSpoints)
190+ use HwU_variables
191 implicit none
192- include "HwU.inc"
193 integer label,nPSpoints,i,j
194 double precision nPSinv,niter,dummy(2)
195 nPSinv=1d0/dble(nPSpoints)
196@@ -257,13 +246,13 @@
197 c histograms the central value should not be zero if any of the other
198 c weights are non-zero.
199 subroutine accumulate_results(label,nPSinv,niter,values)
200+ use HwU_variables
201 implicit none
202- include "HwU.inc"
203 integer label,i,j
204- double precision nPSinv,etot,vtot(max_wgts),niter,y_squared
205+ double precision nPSinv,etot,niter,y_squared
206 $ ,values(2),a1,a2
207- integer error_estimation
208- common /HwU_common2/ error_estimation
209+ double precision,allocatable :: vtot(:)
210+ if (.not. allocated(vtot)) allocate(vtot(nwgts))
211 if (error_estimation.eq.2) then
212 c Use the weighted average bin-by-bin. This is not really justified
213 c for fNLO computations, because for bins with low statistics, the
214@@ -395,14 +384,15 @@
215 c Write the histograms to disk at the end of the run, multiplying the
216 c output by 'xnorm'
217 subroutine HwU_output(unit,xnorm)
218+ use HwU_variables
219 implicit none
220- include "HwU.inc"
221 integer unit,i,j,label
222 integer max_length
223- parameter (max_length=(max_wgts+3)*17)
224- character*(max_length) buffer
225+ character(len=:), allocatable :: buffer
226 character*4 str_nbin
227 double precision xnorm
228+ if (.not. allocated(buffer))
229+ & allocate(character(len=(nwgts+3)*17) :: buffer)
230 c column info: x_min, x_max, y (central value), dy, {extra
231 c weights}.
232 write (unit,'(a$)') '##& xmin'
233@@ -440,9 +430,157 @@
234 write (unit,'(a)') ''
235 write (unit,'(a)') ''
236 enddo
237- return
238- end
239-
240+ call HwU_deallocate_all
241+ deallocate(buffer)
242+ return
243+ end
244+
245+c Clean all the allocatable variables:
246+ subroutine HwU_deallocate_all
247+ use HwU_variables
248+ implicit none
249+ if (allocated(wgts_info)) deallocate(wgts_info)
250+ if (allocated(booked)) deallocate(booked)
251+ if (allocated(title)) deallocate(title)
252+ if (allocated(nbin)) deallocate(nbin)
253+ if (allocated(step)) deallocate(step)
254+ if (allocated(histxl)) deallocate(histxl)
255+ if (allocated(histxm)) deallocate(histxm)
256+ if (allocated(histy)) deallocate(histy)
257+ if (allocated(histy_acc)) deallocate(histy_acc)
258+ if (allocated(histi)) deallocate(histi)
259+ if (allocated(histy2)) deallocate(histy2)
260+ if (allocated(histy_err)) deallocate(histy_err)
261+ if (allocated(p_bin)) deallocate(p_bin)
262+ if (allocated(p_label)) deallocate(p_label)
263+ if (allocated(p_wgts)) deallocate(p_wgts)
264+ return
265+ end
266+
267+
268+ subroutine HwU_allocate_p
269+ use HwU_variables
270+ implicit none
271+ integer,allocatable :: itemp1(:)
272+ double precision, allocatable :: temp2(:,:)
273+ if (.not. allocated(p_bin)) then
274+ allocate(p_bin(max_plots))
275+ allocate(p_label(max_plots))
276+ allocate(p_wgts(nwgts,max_plots))
277+ max_points=max_plots
278+ else
279+ if (np.gt.max_points) then
280+c p_bin
281+ allocate(itemp1(np+max_plots))
282+ itemp1(1:max_points)=p_bin
283+ call move_alloc(itemp1,p_bin)
284+
285+c p_label
286+ allocate(itemp1(np+max_plots))
287+ itemp1(1:max_points)=p_label
288+ call move_alloc(itemp1,p_label)
289+c p_wgts
290+ allocate(temp2(nwgts,np+max_plots))
291+ temp2(1:nwgts,1:max_points)=p_wgts
292+ call move_alloc(temp2,p_wgts)
293+ max_points=np+max_plots
294+ endif
295+ endif
296+ return
297+ end
298+
299+ subroutine HwU_allocate_histo(label,nbin_l)
300+ use HwU_variables
301+ implicit none
302+ logical,allocatable :: ltemp(:)
303+ integer,allocatable :: itemp1(:),itemp2(:,:)
304+ character(len=50),allocatable :: ctemp(:)
305+ double precision, allocatable :: temp1(:),temp2(:,:),temp3(:,:,:)
306+ integer label,i,nbin_l,label_max,nbin_max
307+ logical debug
308+ parameter (debug=.true.)
309+c Check if variables are already allocated. If not, simply allocate a
310+c single histogram
311+ if (.not. allocated(booked)) then
312+ allocate(booked(1))
313+ booked(1)=.false.
314+ allocate(title(1))
315+ allocate(nbin(1))
316+ allocate(step(1))
317+ allocate(histxl(1,nbin_l))
318+ allocate(histxm(1,nbin_l))
319+ allocate(histy(nwgts,1,nbin_l))
320+ allocate(histy_acc(nwgts,1,nbin_l))
321+ allocate(histi(1,nbin_l))
322+ allocate(histy2(1,nbin_l))
323+ allocate(histy_err(1,nbin_l))
324+ max_plots=1
325+ max_bins=nbin_l
326+ endif
327+c If current label is greater than the plots already allocated, increase
328+c the size of the allocated arrays. This is kind of slow, but shouldn't
329+c really matter since it's only done at the start of a run.
330+ if (label.gt.max_plots .or. nbin_l.gt.max_bins) then
331+ label_max=max(label,max_plots)
332+ nbin_max=max(nbin_l,max_bins)
333+c booked
334+ allocate(ltemp(label_max))
335+ ltemp(1:max_plots)=booked
336+ call move_alloc(ltemp,booked)
337+ do i=max_plots+1,label_max
338+ booked(i)=.false. ! histos have not yet been setup
339+ enddo
340+c title
341+ allocate(ctemp(label_max))
342+ ctemp(1:max_plots)=title
343+ call move_alloc(ctemp,title)
344+c nbin
345+ allocate(itemp1(label_max))
346+ itemp1(1:max_plots)=nbin
347+ call move_alloc(itemp1,nbin)
348+c step
349+ allocate(temp1(label_max))
350+ temp1(1:max_plots)=step
351+ call move_alloc(temp1,step)
352+c histxl
353+ allocate(temp2(label_max,nbin_max))
354+ temp2(1:max_plots,1:max_bins)=histxl
355+ call move_alloc(temp2,histxl)
356+c histxm
357+ allocate(temp2(label_max,nbin_max))
358+ temp2(1:max_plots,1:max_bins)=histxm
359+ call move_alloc(temp2,histxm)
360+c histy
361+ allocate(temp3(nwgts,label_max,nbin_max))
362+ temp3(1:nwgts,1:max_plots,1:max_bins)=histy
363+ call move_alloc(temp3,histy)
364+c histy_acc
365+ allocate(temp3(nwgts,label_max,nbin_max))
366+ temp3(1:nwgts,1:max_plots,1:max_bins)=histy_acc
367+ call move_alloc(temp3,histy_acc)
368+c histi
369+ allocate(itemp2(label_max,nbin_max))
370+ itemp2(1:max_plots,1:max_bins)=histi
371+ call move_alloc(itemp2,histi)
372+c histy2
373+ allocate(temp2(label_max,nbin_max))
374+ temp2(1:max_plots,1:max_bins)=histy2
375+ call move_alloc(temp2,histy2)
376+c histy_err
377+ allocate(temp2(label_max,nbin_max))
378+ temp2(1:max_plots,1:max_bins)=histy_err
379+ call move_alloc(temp2,histy_err)
380+c Update maximums
381+ max_plots=label_max
382+ max_bins=nbin_max
383+ elseif (booked(label)) then
384+ write (*,*) 'ERROR in HwU.f: histogram already booked',label
385+ stop
386+ endif
387+ return
388+ end
389+
390+
391 c dummy subroutine
392 subroutine accum(idummy)
393 integer idummy
394@@ -451,3 +589,6 @@
395 subroutine addfil(string)
396 character*(*) string
397 end
398+
399+
400+
401
402=== removed file 'Template/NLO/FixedOrderAnalysis/HwU.inc'
403--- Template/NLO/FixedOrderAnalysis/HwU.inc 2016-02-18 14:05:45 +0000
404+++ Template/NLO/FixedOrderAnalysis/HwU.inc 1970-01-01 00:00:00 +0000
405@@ -1,22 +0,0 @@
406-* -*-fortran-*-
407-
408- integer max_plots,max_bins,max_wgts,max_points
409- parameter (max_plots=200)
410- parameter (max_bins=100)
411- parameter (max_wgts=1024)
412- parameter (max_points=max_plots*40)
413-
414- logical booked(max_plots)
415- integer nbin(max_plots),nwgts,np,p_bin(max_points)
416- & ,p_label(max_points),histi(max_plots,max_bins)
417- character*50 title(max_plots)
418- character*50 wgts_info(max_wgts)
419- double precision histy(max_wgts,max_plots,max_bins)
420- $ ,histy_acc(max_wgts,max_plots,max_bins),histy2(max_plots
421- $ ,max_bins),histy_err(max_plots,max_bins),histxl(max_plots
422- $ ,max_bins),histxm(max_plots,max_bins),step(max_plots)
423- $ ,p_wgts(max_wgts,max_points)
424-
425- common/HwU_common/histy,histy_acc,histy2,histy_err,histxl,histxm
426- & ,p_wgts,step,histi,nbin,p_bin,p_label,np,nwgts
427- & ,booked,title,wgts_info
428
429=== modified file 'Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f'
430--- Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f 2014-12-03 11:51:14 +0000
431+++ Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f 2017-05-09 14:02:00 +0000
432@@ -17,12 +17,12 @@
433 c o) The first argument is an integer that labels the histogram. In
434 c the analysis_end and analysis_fill subroutines this label is used
435 c to keep track of the histogram. The label should be a number
436-c between 1 and max_plots=200 (can be increased in HwU.inc).
437+c starting at 1 and be increased for each plot.
438 c o) The second argument is a string that will apear above the
439 c histogram. Do not use brackets "(" or ")" inside this string.
440 c o) The third, forth and fifth arguments are the number of bis, the
441 c lower edge of the first bin and the upper edge of the last
442-c bin. There is a maximum of 100 bins per histogram.
443+c bin.
444 c o) When including scale and/or PDF uncertainties, declare a
445 c histogram for each weight, and compute the uncertainties from the
446 c final set of histograms
447
448=== removed symlink 'Template/NLO/MCatNLO/include/HwU.inc'
449=== target was u'../../FixedOrderAnalysis/HwU.inc'

Subscribers

People subscribed via source and target branches

to all changes: