Commit 9bebf7a9 authored by Juergen Reuter's avatar Juergen Reuter
Browse files

b41f63e Merge branch '431-powheg-events-with-real-partition-for-ee-jj' into 'master'

parent db0f029c
......@@ -39,7 +39,7 @@ scale = let int njet = count [jet] in
sqrts
endif
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 10 GeV
process ppzj_wp = pr, pr => Z, jet { nlo_calculation = real $restrictions="!W+:W-" }
......
......@@ -53,7 +53,7 @@ $lhapdf_file = "MSTW2008nlo68cl"
sqrts = 1.30000E+04
$rng_method = "rng_stream"
$integration_method = "vamp2"
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 1.00000E+01
$restrictions = "!W+:W-"
| Process library 'nlo_ppzj_real_partition_lib': recorded process 'ppzj_wp'
......
......@@ -22,7 +22,7 @@ mtop = 175 GeV
?alphas_is_fixed = false
?alphas_from_mz = true
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 10 GeV
$fks_mapping_type = "resonances"
......
......@@ -53,10 +53,13 @@ simulate (powheg_1_p1)
n_events = 2
simulate (powheg_1_p1)
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 5 GeV
process powheg_1_p2 = E1, e1 => t, T { nlo_calculation = full }
integrate (powheg_1_p2) { iterations = 1:100:"gw" }
expect (integral(powheg_1_p1) == integral(powheg_1_p2))
{tolerance = 3 * sqrt(error(powheg_1_p1)**2 + error(powheg_1_p2)**2)}
simulate (powheg_1_p2)
......@@ -24,7 +24,7 @@ process real_partition_1_p1 = e1, E1 => u, U { nlo_calculation = real }
integrate (real_partition_1_p1)
seed = 1
?nlo_use_real_partition = true
$real_partition_mode = "on"
# Large scale to reduce required statistics
real_partition_scale = 100 GeV
......
......@@ -132,7 +132,7 @@ n_events = 2
| Events: actual unweighting efficiency = 100.00 %
| Events: closing ASCII file 'powheg_1_p1.debug'
Warning: There have been 1 ( 11.11%) POWHEG grid excesses.
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 5.00000E+00
| Process library 'powheg_1_lib': unloading
| Process library 'powheg_1_lib': open
......
......@@ -132,7 +132,7 @@ n_events = 2
| Events: actual unweighting efficiency = 100.00 %
| Events: closing ASCII file 'powheg_1_p1.debug'
Warning: There have been 1 ( 11.11%) POWHEG grid excesses.
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 5.00000E+00
| Process library 'powheg_1_lib': unloading
| Process library 'powheg_1_lib': open
......
......@@ -132,7 +132,7 @@ n_events = 2
| Events: actual unweighting efficiency = 100.00 %
| Events: closing ASCII file 'powheg_1_p1.debug'
Warning: There have been 1 ( 11.11%) POWHEG grid excesses.
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 5.00000E+00
| Process library 'powheg_1_lib': unloading
| Process library 'powheg_1_lib': open
......
......@@ -9,7 +9,7 @@ SM.mtop => 1.75000E+02
?use_vamp_equivalences = false
?alphas_is_fixed = false
?alphas_from_mz = true
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 1.00000E+01
$fks_mapping_type = "resonances"
[user variable] pr = PDG(2, 1, -2, -1)
......
......@@ -539,7 +539,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc = Fortran-compiler
......
......@@ -76,7 +76,7 @@ Warning: No cuts have been defined.
1 0 2.292E+02 2.19E+01 9.56 0.00 20.1
|=============================================================================|
seed = 1
?nlo_use_real_partition = true
$real_partition_mode = "on"
real_partition_scale = 1.00000E+02
| Process library 'real_partition_1_lib': unloading
| Process library 'real_partition_1_lib': open
......
......@@ -190,6 +190,7 @@ $gosam_fc = ""
$dalitz_plot = ""
$nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
$real_partition_mode = "default"
$fc => Fortran-compiler
$fcflags => Fortran-flags
$fclibs => Fortran-libs
......@@ -343,7 +344,6 @@ $fclibs => Fortran-libs
?keep_failed_events = false
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
?nlo_reuse_amplitudes_fks = false
?rebuild_library = true
?recompile_library = false
......@@ -739,7 +739,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc => Fortran-compiler
......@@ -1173,7 +1173,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc => Fortran-compiler
......
......@@ -239,7 +239,6 @@ seed = 2
?keep_failed_events = false
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
?nlo_reuse_amplitudes_fks = false
?rebuild_library = true
?recompile_library = false
......@@ -338,6 +337,7 @@ $gosam_fc = ""
$dalitz_plot = ""
$nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
$real_partition_mode = "default"
$fc => Fortran-compiler
$fcflags => Fortran-flags
$fclibs => Fortran-libs
......
......@@ -393,7 +393,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc => "Fortran-compiler"
......
......@@ -407,7 +407,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc => "Fortran-compiler"
......
......@@ -422,7 +422,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc => "Fortran-compiler"
......@@ -942,7 +942,7 @@ $nlo_correction_type = "QCD"
$exclude_gauge_splittings = "c:b:t:e2:e3"
?nlo_use_born_scale = false
?nlo_cut_all_real_sqmes = false
?nlo_use_real_partition = false
$real_partition_mode = "default"
real_partition_scale = 1.000000000000E+01
?nlo_reuse_amplitudes_fks = false
$fc => "Fortran-compiler"
......
......@@ -6674,7 +6674,8 @@ end module real_subtraction
integer, parameter, public :: POWHEG = 2
@ %def real subtraction parameters
@
@ Translates the purpose into a string.
% However, purpose is never set and this routine is never used.
<<real subtraction: public>>=
public :: this_purpose
<<real subtraction: procedures>>=
......@@ -7699,6 +7700,7 @@ is true).
sqme_alr = rsub%evaluate_emitter_region (alr, this_emitter, i_phs, i_res, &
alpha_coupling)
if (rsub%purpose == INTEGRATION .or. rsub%purpose == FIXED_ORDER_EVENTS) then
!!! (PS-2021-09-13) This is always executed as purpose is never set.
i_con = rsub%get_i_contributor (alr)
sqme_alr = sqme_alr * rsub%get_phs_factor (i_con)
end if
......@@ -8334,6 +8336,7 @@ However, it still remains to be motivated analytically and to be validated.
case (INTEGRATION, FIXED_ORDER_EVENTS)
sqme = sqme * xi**2 / xi_tilde * real_kinematics%jac(i_phs)%jac(1)
case (POWHEG)
!!! (PS 2021-09-21) This has been implemented but is never executed so far.
if (.not. isr) then
s = real_kinematics%cms_energy2
sqme = sqme * real_kinematics%jac(i_phs)%jac(1) * s / (8 * twopi3) * xi
......
......@@ -3867,10 +3867,16 @@ but we assume that the generator cuts applied for parton showering are weak cuts
avoid divergent phase space regions, thus $\mathcal{B}$ is large and does not influence $N$
so we can skip these phase space points.
A more appropriate way to solve this problem of Born zeroes is given by the real partition.
With real partition, there are two possible cases if $\mathcal{B} = 0$: if the real phase
space point belongs to the real singular, the real kinematic is Born-like and very likely
also fails the cuts. If the real kinematic is not Born-like, it belongs to the real finite.
In both cases, it is $\mathcal{R} = \mathcal{B} = 0$ and we can savely skip the phase space
point.
In cases where just the real kinematic failed the cuts, $\mathcal{R}$ was artificially set
to zero and we have no way to recover it. These cases will have no influence on $N$.
It might be possible that this leads to an underestimate of $N$ in some bins.
The number of POWHEG excess events will be an indicator for the severity of this effect.
to zero and we have no way to recover it. We skip these phase space points too. If we apply
the same cuts during integration and event generation, these points are not relevant.
<<POWHEG matching: powheg matching hook: TBP>>=
procedure :: evaluate => powheg_matching_hook_evaluate
<<POWHEG matching: procedures>>=
......
......@@ -2692,7 +2692,7 @@ slash, so the default is [[.]].
cmd = ""
end if
end function workspace_cmd
@ %def workspace_prefix
@ %def workspace_path
@ %def workspace_cmd
......@@ -3412,12 +3412,12 @@ NOTE: temporarily made public.
class(test_writer_1_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_1_before_compile
subroutine test_writer_1_after_compile (writer, id)
class(test_writer_1_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_1_after_compile
@ %def test_writer_1_type_name
@ %def test_writer_1_mk test_writer_1_if
@ %def test_writer_1_md5sum test_writer_1_int_sub
......@@ -3538,12 +3538,12 @@ lines. Since all procedures are NOPASS, we can reuse two of the TBP.
class(test_writer_2_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_2_before_compile
subroutine test_writer_2_after_compile (writer, id)
class(test_writer_2_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_2_after_compile
@ %def test_writer_2_type_name test_writer_2_wr
@ %def test_writer_2_before_compile test_writer_2_after_compile
@
......@@ -3891,12 +3891,12 @@ We will reuse this later, therefore public.
class(test_writer_4_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_4_before_compile
subroutine test_writer_4_after_compile (writer, id)
class(test_writer_4_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_4_after_compile
@ %def test_writer_2_type_name test_writer_4_wr
@ %def test_writer_4_before_compile test_writer_4_after_compile
@
......@@ -4214,7 +4214,7 @@ so no wrapper is needed.
if (.not. verbose) then
write (unit, "(5A)") TAB // '@echo " FC " $@'
end if
write (unit, "(5A)") TAB, "$(LTFCOMPILE) $<"
write (unit, "(5A)") TAB, "$(LTFCOMPILE) $<"
end subroutine test_writer_5_mk
subroutine test_writer_5_src (writer, id)
......@@ -4248,12 +4248,12 @@ so no wrapper is needed.
class(test_writer_5_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_5_before_compile
subroutine test_writer_5_after_compile (writer, id)
class(test_writer_5_t), intent(in) :: writer
type(string_t), intent(in) :: id
end subroutine test_writer_5_after_compile
@ %def test_writer_5_type_name test_writer_5_mk
@ %def test_writer_5_if
@ %def test_writer_5_before_compile test_writer_5_after_compile
......@@ -4589,7 +4589,7 @@ Makefile writer.
write (unit, "(5A)") char (id), ".lo: ", char (id), ".c"
if (.not. verbose) then
write (unit, "(5A)") TAB // '@echo " FC " $@'
end if
end if
write (unit, "(5A)") TAB, "$(LTCCOMPILE) $<"
end subroutine test_writer_6_mk
......@@ -5417,7 +5417,7 @@ the provided templates.
@ %def process_component_def_show
@ Compute the MD5 sum of a process component. We reset the stored MD5
sum to the empty string (so a previous value is not included in the
calculation), the write a temporary file and calculate the MD5 sum of
calculation), then write a temporary file and calculate the MD5 sum of
that file.
This implies that all data that are displayed by the [[write]] method
......@@ -5708,7 +5708,9 @@ abstract type.
logical :: valid
i_skip = 0; if (present (i_skip_in)) i_skip = i_skip_in
n = count (component%associated_components /= 0) - 1
if (i_skip > 0) n = n - 1
if (i_skip > 0) then
if (component%associated_components(i_skip) > 0) n = n - 1
end if
allocate (list (n))
j = 1
do i = 1, size(component%associated_components)
......@@ -6569,7 +6571,7 @@ Return the number of processes supported by the library.
entry => current
end function process_def_list_get_process_def_ptr
@ %def process_def_list_get_process_def_ptr
@ %def process_def_list_get_process_def_ptr
@ Return true if a given process is in the library.
<<Process libraries: process def list: TBP>>=
procedure :: contains => process_def_list_contains
......@@ -7248,7 +7250,7 @@ The [[workspace]] optional argument puts any library code in a subdirectory.
end subroutine process_library_write_makefile
@ %def process_library_write_makefile
@
@
@ Write the driver source code for the library to file, if there are
external processes.
<<Process libraries: process library: TBP>>=
......@@ -7275,7 +7277,7 @@ external processes.
call msg_message ("Process library '" // char (lib%basename) &
// "': writing driver")
unit = free_unit ()
open (unit, &
open (unit, &
file = char (workspace_prefix (workspace) &
& // lib%driver%basename // ".f90"), &
status="replace", action="write")
......
......@@ -216,7 +216,7 @@ we introduce a MD5 sum calculator.
end subroutine mci_compute_md5sum
end interface
@ %def mci_compute_md5sum@
@ %def mci_compute_md5sum
@ Record the index of the MCI object within a process. For
multi-component processes with more than one integrator, the
integrator should know about its own index, so file names can be
......@@ -11510,6 +11510,34 @@ may contain extra passes, these are ignored.
end subroutine list_pass_update_from_ref
@ %def list_pass_update_from_ref
<<MCI vamp2: list pass: TBP>>=
procedure :: has_last_integral => list_pass_has_last_integral
procedure :: get_last_integral => list_pass_get_last_integral
<<MCI vamp2: procedures>>=
function list_pass_has_last_integral(self) result (flag)
class(list_pass_t), intent(in) :: self
logical :: flag
flag = associated(self%current)
if (flag) flag = self%current%integral_defined
end function list_pass_has_last_integral
subroutine list_pass_get_last_integral(self, integral, error, efficiency)
class(list_pass_t), intent(in) :: self
real(default), intent(out) :: integral
real(default), intent(out) :: error
real(default), intent(out) :: efficiency
if (self%has_last_integral()) then
integral = self%current%get_integral()
error = self%current%get_error()
efficiency = self%current%get_efficiency()
else
integral = 0
error = 0
efficiency = 0
end if
end subroutine list_pass_get_last_integral
@ %def list_pass_has_last_integral list_pass_get_last_integral
@ Output. Write the complete linked list to the specified unit.
<<MCI vamp2: list pass: TBP>>=
procedure :: write => list_pass_write
......@@ -12249,10 +12277,11 @@ that the user wants to perform and average all iterations in this pass.
select type (mci_ref)
type is (mci_vamp2_t)
call mci%list_pass%update_from_ref (mci_ref%list_pass, success)
if (mci%list_pass%current%integral_defined) then
mci%integral = mci%list_pass%current%get_integral ()
mci%error = mci%list_pass%current%get_error ()
mci%efficiency = mci%list_pass%current%get_efficiency ()
if (mci%list_pass%has_last_integral()) then
call mci%list_pass%get_last_integral( &
integral = mci%integral, &
error = mci%error, &
efficiency = mci%efficiency)
mci%integral_known = .true.
mci%error_known = .true.
mci%efficiency_known = .true.
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment