Index: trunk/src/recola/recola.nw =================================================================== --- trunk/src/recola/recola.nw (revision 8178) +++ trunk/src/recola/recola.nw (revision 8179) @@ -1,3078 +1,3225 @@ % -*- ess-noweb-default-code-mode: f90-mode; noweb-default-code-mode: f90-mode; -*- % WHIZARD code as NOWEB source: interface to Recola 1-loop library @ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \chapter{Recola Interface} \section{Recola wrapper} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <<[[recola_wrapper.f90]]>>= <> module recola_wrapper use recola !NODEP! use kinds <> use constants, only: zero use diagnostics, only: msg_fatal, msg_message, msg_debug, msg_debug2, D_ME_METHODS use io_units, only: given_output_unit <> <> <> <> <> contains <> end module recola_wrapper @ %def recola_wrapper @ <>= public :: rclwrap_is_active <>= logical, parameter :: rclwrap_is_active = .true. @ %def rclwrap_is_active @ Returns the particle string corresponding to a pdg code used in the Recola process definition <>= public :: get_recola_particle_string <>= elemental function get_recola_particle_string (pdg) result (name) type(string_t) :: name integer, intent(in) :: pdg select case (pdg) case (1) name = var_str ("d") case (-1) name = var_str ("d~") case (2) name = var_str ("u") case (-2) name = var_str ("u~") case (3) name = var_str ("s") case (-3) name = var_str ("s~") case (4) name = var_str ("c") case (-4) name = var_str ("c~") case (5) name = var_str ("b") case (-5) name = var_str ("b~") case (6) name = var_str ("t") case (-6) name = var_str ("t~") case (11) name = var_str ("e-") case (-11) name = var_str ("e+") case (12) name = var_str ("nu_e") case (-12) name = var_str ("nu_e~") case (13) name = var_str ("mu-") case (-13) name = var_str ("mu+") case (14) name = var_str ("nu_mu") case (-14) name = var_str ("nu_mu~") case (15) name = var_str ("tau-") case (-15) name = var_str ("tau+") case (16) name = var_str ("nu_tau") case (-16) name = var_str ("nu_tau~") case (21) name = var_str ("g") case (22) name = var_str ("A") case (23) name = var_str ("Z") case (24) name = var_str ("W+") case (-24) name = var_str ("W-") case (25) name = var_str ("H") end select end function get_recola_particle_string @ %def get_recola_particle_string @ -<>= - public :: rclwrap_define_process <>= subroutine rclwrap_define_process (id, process_string, order) integer, intent(in) :: id type(string_t), intent(in) :: process_string type(string_t), intent(in) :: order call msg_debug2 (D_ME_METHODS, "define_process_rcl") call define_process_rcl (id, char (process_string), char (order)) end subroutine rclwrap_define_process @ %def rclwrap_define_process +@ This defines a wrapper for the information required to define a RECOLA +process. It is used to collect the process definitions in an array. +<>= + type :: rcl_process_t + private + integer :: id + type(string_t) :: process_string + type(string_t) :: order + contains + <> + end type rcl_process_t + +@ %def rcl_process_t +@ +<>= + interface rcl_process_t + module procedure new_rcl_process_t + end interface + +@ %def rcl_process_t +@ +<>= + function new_rcl_process_t (id, process_string, order) + integer, intent(in) :: id + type(string_t), intent(in) :: process_string, order + type(rcl_process_t) :: new_rcl_process_t + new_rcl_process_t%id = id + new_rcl_process_t%process_string = process_string + new_rcl_process_t%order = order + end function new_rcl_process_t + +@ %def new_rcl_process_t +<>= + procedure :: get_params => rcl_process_get_params +<>= + subroutine rcl_process_get_params (prc, id, process_string, order) + class(rcl_process_t), intent(in) :: prc + integer, intent(out) :: id + type(string_t), intent(out) :: process_string + type(string_t), intent(out) :: order + id = prc%id + process_string = prc%process_string + order = prc%order + end subroutine rcl_process_get_params + +@ %def rcl_process_get_params +@ Output. +<>= + procedure :: write => rcl_process_write +<>= + subroutine rcl_process_write (object, unit) + class(rcl_process_t), intent(in) :: object + integer, intent(in), optional :: unit + integer :: u + u = given_output_unit (unit) + write (u, "(1x,A,I0,2(1x,A,1x))") "RECOLA process:", & + "id=", object%id, "process_string=", char(object%process_string), & + "order=", char(object%order) + end subroutine rcl_process_write + +@ %def rcl_process_write @ This defines a singleton object, located in this module only, that -controls RECOLA initialization. When WHIZARD compiles processes, it -should also run the RECOLA "`controller"', which actually initializes -RECOLA for integration. The main complication is that this has to be +controls RECOLA initialization and process management. When WHIZARD +compiles processes, it should also run the RECOLA "`controller"', +which actually initializes RECOLA for integration and manages process +information in an array. The main complication is that this has to be done after all processes have been registered, and cannot be redone. We could work with module variables directly, but the singleton pattern, e.g., allows us to work with multiple RECOLA instances, if this becomes possible in the future. Type and object can be private. <>= type :: rcl_controller_t private logical :: active = .false. + logical :: defined = .false. logical :: done = .false. integer :: recola_id = 0 + type(rcl_process_t), dimension (:), allocatable :: processes + integer :: n_processes = 0 contains <> end type rcl_controller_t @ %def rcl_controller_t <>= type(rcl_controller_t) :: rcl_controller @ %def rcl_controller +@ Add a RECOLA process to the controller. This will make sure that +processes can be redefined if additional definitions are to be made +after process generation. +<>= + procedure :: add_process => rcl_controller_add_process +<>= + subroutine rcl_controller_add_process (rcl, process) + class(rcl_controller_t), intent(inout) :: rcl + type(rcl_process_t), intent(in) :: process + type(rcl_process_t), dimension (:), allocatable :: temp + if (rcl%n_processes == size(rcl%processes)) then + allocate( temp(2 * rcl%n_processes) ) + temp(:rcl%n_processes) = rcl%processes + call move_alloc(temp, rcl%processes) + end if + rcl%processes(rcl%n_processes + 1) = process + rcl%n_processes = rcl%n_processes + 1 + end subroutine rcl_controller_add_process + +@ %def rcl_controller_add_process +@ Define all processes added to the controller, and only them. +If processes have been defined before, RECOLA is reset. +<>= + procedure :: define_processes => rcl_controller_define_processes +<>= + subroutine rcl_controller_define_processes (rcl) + class(rcl_controller_t), intent(inout) :: rcl + integer :: id, i + type(string_t) :: process_string + type(string_t) :: order + if (rcl%defined) then + if (.not. rcl%done) call rclwrap_generate_processes () + call msg_debug2 (D_ME_METHODS, "reset_recola_rcl") + call reset_recola_rcl () + end if + do i = 1, rcl%n_processes + call rcl%processes(i)%get_params(id, process_string, order) + call rclwrap_define_process (id, process_string, order) + end do + rcl%defined = .true. + rcl%done = .false. + end subroutine rcl_controller_define_processes + +@ %def rcl_controller_define_processes @ Revert to initial state. Also, reset RECOLA (only if it has already done something). <>= procedure :: reset => rcl_controller_reset <>= subroutine rcl_controller_reset (rcl) class(rcl_controller_t), intent(inout) :: rcl if (rcl%active .or. rcl%done) then call msg_debug2 (D_ME_METHODS, "reset_recola_rcl") + if (allocated (rcl%processes)) deallocate (rcl%processes) call reset_recola_rcl () end if rcl%active = .false. + rcl%defined = .false. rcl%done = .false. rcl%recola_id = 0 + rcl%n_processes = 0 end subroutine rcl_controller_reset @ %def rcl_controller_reset @ Output. <>= procedure :: write => rcl_controller_write <>= subroutine rcl_controller_write (object, unit) class(rcl_controller_t), intent(in) :: object integer, intent(in), optional :: unit integer :: u u = given_output_unit (unit) - write (u, "(1x,A,2(1x,A,L1),1(1x,A,I0))") "RECOLA controller:", & + write (u, "(1x,A,2(1x,A,L1),2(1x,A,I0))") "RECOLA controller:", & "active=", object%active, "done=", object%done, & - "id=", object%recola_id + "id=", object%recola_id, "n_processes=", object%n_processes end subroutine rcl_controller_write @ %def rcl_controller_write @ Return a new numeric process ID, incrementing the counter once. <>= procedure :: get_new_id => rcl_controller_get_new_id <>= subroutine rcl_controller_get_new_id (object, id) class(rcl_controller_t), intent(inout) :: object integer, intent(out) :: id object%recola_id = object%recola_id + 1 id = object%recola_id end subroutine rcl_controller_get_new_id @ %def rcl_controller_get_new_id @ Return the current numeric process ID without incrementing the counter. <>= procedure :: get_current_id => rcl_controller_get_current_id <>= subroutine rcl_controller_get_current_id (object, id) class(rcl_controller_t), intent(inout) :: object integer, intent(out) :: id id = object%recola_id end subroutine rcl_controller_get_current_id @ %def rcl_controller_get_current_id @ Do not allow activation if processes have been calculated previously. Otherwise set the flag. <>= procedure :: activate => rcl_controller_activate <>= subroutine rcl_controller_activate (rcl) class(rcl_controller_t), intent(inout) :: rcl - if (rcl_controller%done) then - call msg_fatal ("Recola interface: attempt at adding processes after initialization") - else - rcl_controller%active = .true. - end if + if ( .not. allocated(rcl%processes) ) allocate ( rcl%processes(10) ) + rcl_controller%active = .true. end subroutine rcl_controller_activate @ %def rcl_controller_activate @ Start process initialization by calling the RECOLA API. Do not allow this twice (skip silently), and skip anyway if there is no activation. <>= procedure :: generate_processes => rcl_controller_generate_processes <>= subroutine rcl_controller_generate_processes (rcl) class(rcl_controller_t), intent(inout) :: rcl if (rcl_controller%active) then if (.not. rcl_controller%done) then call msg_message ("Recola: preparing processes for integration") call generate_processes_rcl () rcl_controller%done = .true. end if end if end subroutine rcl_controller_generate_processes @ %def rcl_controller_generate_processes @ Return a new numeric RECOLA process ID. The singleton nature of the controller guarantees that the ID is unique. <>= public :: rclwrap_get_new_recola_id <>= subroutine rclwrap_get_new_recola_id (id) integer, intent(out) :: id call rcl_controller%get_new_id (id) end subroutine rclwrap_get_new_recola_id @ %def rclwrap_get_new_recola_id -@ Return the number of registered processes, which coincides with the -current numeric RECOLA process ID. +@ Return the current numeric RECOLA process ID. This coincides with the amount +of IDs currently in use. <>= - public :: rclwrap_get_n_processes + public :: rclwrap_get_current_recola_id <>= - function rclwrap_get_n_processes () result (n) + function rclwrap_get_current_recola_id () result (n) integer :: n call rcl_controller%get_current_id (n) - end function rclwrap_get_n_processes + end function rclwrap_get_current_recola_id -@ %def rclwrap_get_n_processes +@ %def rclwrap_get_current_recola_id @ This procedure records the fact that there is a recola process pending, so we will have to call [[generate_processes]] before we can calculate anything with Recola. <>= public :: rclwrap_request_generate_processes <>= subroutine rclwrap_request_generate_processes () call msg_debug2 (D_ME_METHODS, "request_generate_processes_rcl") call rcl_controller%activate () end subroutine rclwrap_request_generate_processes @ %def rclwrap_request_generate_processes -@ We call this after all processes have been registered, so RECOLA can -initialize itself for integration. +@ Add a process to be defined later +<>= + public :: rclwrap_add_process +<>= + subroutine rclwrap_add_process (id, process_string, order) + integer, intent(in) :: id + type(string_t), intent(in) :: process_string, order + type(rcl_process_t) :: prc + call msg_debug2 (D_ME_METHODS, "add_process_rcl: id", id) + prc = rcl_process_t (id, process_string, order) + call rcl_controller%add_process (prc) + end subroutine rclwrap_add_process + +@ %def rclwrap_add_process +@ Define all added processes. Reset if processes were already defined. +<>= + public :: rclwrap_define_processes +<>= + subroutine rclwrap_define_processes () + call msg_debug2 (D_ME_METHODS, "define_processes_rcl") + call rcl_controller%define_processes () + end subroutine rclwrap_define_processes + +@ %def rclwrap_define_processes +@ We call this after all processes have been added and defined, +so RECOLA can initialize itself for integration. <>= public :: rclwrap_generate_processes <>= subroutine rclwrap_generate_processes () call msg_debug2 (D_ME_METHODS, "generate_processes_rcl") call rcl_controller%generate_processes () end subroutine rclwrap_generate_processes @ %def rclwrap_generate_processes @ <>= public :: rclwrap_compute_process <>= subroutine rclwrap_compute_process (id, p, order, sqme) integer, intent(in) :: id real(double), intent(in), dimension(:,:) :: p character(len=*), intent(in) :: order real(double), intent(out), dimension(0:1), optional :: sqme call msg_debug2 (D_ME_METHODS, "compute_process_rcl") call compute_process_rcl (id, p, order, sqme) end subroutine rclwrap_compute_process @ %def rclwrap_compute_process @ <>= public :: rclwrap_get_amplitude <>= subroutine rclwrap_get_amplitude (id, g_power, order, col, hel, amp) integer, intent(in) :: id, g_power character(len=*), intent(in) :: order integer, dimension(:), intent(in) :: col, hel complex(double), intent(out) :: amp call msg_debug2 (D_ME_METHODS, "get_amplitude_rcl") call get_amplitude_rcl (id, g_power, order, col, hel, amp) end subroutine rclwrap_get_amplitude @ %def rclwrap_get_amplitude @ <>= public :: rclwrap_get_squared_amplitude <>= subroutine rclwrap_get_squared_amplitude (id, alphas_power, order, sqme) integer, intent(in) :: id, alphas_power character(len=*), intent(in) :: order real(double), intent(out) :: sqme call msg_debug2 (D_ME_METHODS, "get_squared_amplitude_rcl") call get_squared_amplitude_rcl (id, alphas_power, order, sqme) end subroutine rclwrap_get_squared_amplitude @ %def rclwrap_get_squared_amplitude @ <>= public :: rclwrap_set_pole_mass <>= subroutine rclwrap_set_pole_mass (pdg_id, mass, width) integer, intent(in) :: pdg_id real(double), intent(in) :: mass, width call msg_debug2 (D_ME_METHODS, "rclwrap_set_pole_mass of ", pdg_id) select case (abs(pdg_id)) case (11) if (width > zero) & call msg_fatal ("Recola pole mass: Attempting to set non-zero electron width!") call set_pole_mass_electron_rcl (mass) case (13) call set_pole_mass_muon_rcl (mass, width) case (15) call set_pole_mass_tau_rcl (mass, width) case (1) if (width > zero) & call msg_fatal ("Recola pole mass: Attempting to set non-zero down-quark width!") call set_pole_mass_down_rcl (mass) case (2) if (width > zero) & call msg_fatal ("Recola pole mass: Attempting to set non-zero up-quark width!") call set_pole_mass_up_rcl (mass) case (3) if (width > zero) & call msg_fatal ("Recola pole mass: Attempting to set non-zero strange-quark width!") call set_pole_mass_strange_rcl (mass) case (4) call set_pole_mass_charm_rcl (mass, width) case (5) call set_pole_mass_bottom_rcl (mass, width) case (6) call set_pole_mass_top_rcl (mass, width) case (23) call set_pole_mass_z_rcl (mass, width) case (24) call set_pole_mass_w_rcl (mass, width) case (25) call set_pole_mass_h_rcl (mass, width) case default call msg_fatal ("Recola pole mass: Unsupported particle") end select end subroutine rclwrap_set_pole_mass @ %def rclwrap_set_pole_mass @ <>= public :: rclwrap_set_onshell_mass <>= subroutine rclwrap_set_onshell_mass (pdg_id, mass, width) integer, intent(in) :: pdg_id real(double), intent(in) :: mass, width call msg_debug2 (D_ME_METHODS, "rclwrap_set_onshell_mass of ", pdg_id) select case (abs(pdg_id)) case (23) call set_onshell_mass_z_rcl (mass, width) case (24) call set_onshell_mass_w_rcl (mass, width) case default call msg_fatal ("Recola onshell mass: Only for W and Z") end select end subroutine rclwrap_set_onshell_mass @ %def rclwrap_set_onshell_mass @ <>= public :: rclwrap_use_gfermi_scheme <>= subroutine rclwrap_use_gfermi_scheme (gf) real(double), intent(in), optional :: gf call msg_debug2 (D_ME_METHODS, "use_gfermi_scheme_rcl", & real(gf, kind=default)) call use_gfermi_scheme_rcl (gf) end subroutine rclwrap_use_gfermi_scheme @ %def rclwrap_use_gfermi_scheme @ <>= public :: rclwrap_set_light_fermions <>= subroutine rclwrap_set_light_fermions (m) real(double), intent(in) :: m call msg_debug2 (D_ME_METHODS, "set_light_fermions_rcl", & real(m, kind=default)) call set_light_fermions_rcl (m) end subroutine rclwrap_set_light_fermions @ %def rclwrap_set_light_fermions @ <>= public :: rclwrap_set_light_fermion <>= subroutine rclwrap_set_light_fermion (pdg_id) integer, intent(in) :: pdg_id call msg_debug2 (D_ME_METHODS, "rclwrap_set_light_fermion", pdg_id) select case (abs(pdg_id)) case (1) call set_light_down_rcl () case (2) call set_light_up_rcl () case (3) call set_light_strange_rcl () case (4) call set_light_charm_rcl () case (5) call set_light_bottom_rcl () case (6) call set_light_top_rcl () case (11) call set_light_electron_rcl () case (13) call set_light_muon_rcl () case (15) call set_light_tau_rcl () end select end subroutine rclwrap_set_light_fermion @ %def rclwrap_set_light_fermion @ <>= public :: rclwrap_unset_light_fermion <>= subroutine rclwrap_unset_light_fermion (pdg_id) integer, intent(in) :: pdg_id call msg_debug2 (D_ME_METHODS, "rclwrap_unset_light_fermion", pdg_id) select case (abs(pdg_id)) case (1) call unset_light_down_rcl () case (2) call unset_light_up_rcl () case (3) call unset_light_strange_rcl () case (4) call unset_light_charm_rcl () case (5) call unset_light_bottom_rcl () case (6) call unset_light_top_rcl () case (11) call unset_light_electron_rcl () case (13) call unset_light_muon_rcl () case (15) call unset_light_tau_rcl () end select end subroutine rclwrap_unset_light_fermion @ %def rclwrap_unset_light_fermion @ <>= public :: rclwrap_set_onshell_scheme <>= subroutine rclwrap_set_onshell_scheme call msg_debug2 (D_ME_METHODS, "set_on_shell_scheme_rcl") call set_on_shell_scheme_rcl () end subroutine rclwrap_set_onshell_scheme @ %def rclwrap_set_onshell_scheme @ <>= public :: rclwrap_set_alpha_s <>= subroutine rclwrap_set_alpha_s (alpha_s, mu, nf) real(double), intent(in) :: alpha_s, mu integer, intent(in) :: nf call msg_debug2 (D_ME_METHODS, "set_alphas_rcl") call set_alphas_rcl (alpha_s, mu, nf) end subroutine rclwrap_set_alpha_s @ %def rclwrap_set_alpha_s @ <>= public :: rclwrap_get_alpha_s <>= function rclwrap_get_alpha_s () result (alpha_s) real(double) :: alpha_s call msg_debug2 (D_ME_METHODS, "get_alphas_rcl") call get_alphas_rcl (alpha_s) end function rclwrap_get_alpha_s @ %def rclwrap_get_alpha_s @ <>= public :: rclwrap_get_helicity_configurations <>= subroutine rclwrap_get_helicity_configurations (id, hel) integer, intent(in) :: id integer, intent(inout), dimension(:,:), allocatable :: hel call get_helicity_configurations_rcl (id, hel) end subroutine rclwrap_get_helicity_configurations @ %def rclwrap_get_helicity_configurations @ <>= public :: rclwrap_get_color_configurations <>= subroutine rclwrap_get_color_configurations (id, col) integer, intent(in) :: id integer, intent(out), dimension(:,:), allocatable :: col call get_colour_configurations_rcl (id, col) end subroutine rclwrap_get_color_configurations @ %def rclwrap_get_color_configurations @ Selects dimensional regularization for soft singularities. <>= public :: rclwrap_use_dim_reg_soft <>= subroutine rclwrap_use_dim_reg_soft () call msg_debug2 (D_ME_METHODS, "use_dim_reg_soft_rcl") call use_dim_reg_soft_rcl () end subroutine rclwrap_use_dim_reg_soft @ %def rclwrap_use_dim_reg_soft @ Selects mass regularization for soft singularities and sets the mass regulator in GeV to [[m]]. <>= public :: rclwrap_use_mass_reg_soft <>= subroutine rclwrap_use_mass_reg_soft (m) real(double), intent(in) :: m call msg_debug2 (D_ME_METHODS, "use_mass_reg_soft_rcl") call use_mass_reg_soft_rcl (m) end subroutine rclwrap_use_mass_reg_soft @ %def rclwrap_use_mass_reg_soft @ Sets the UV pole parameterization $\Delta_{UV}$. <>= public :: rclwrap_set_delta_uv <>= subroutine rclwrap_set_delta_uv (d) real(double), intent(in) :: d call msg_debug2 (D_ME_METHODS, "set_delta_uv_rcl") call set_delta_uv_rcl (d) end subroutine rclwrap_set_delta_uv @ %def rclwrap_set_delta_uv @ <>= public :: rclwrap_set_mu_uv <>= subroutine rclwrap_set_mu_uv (mu) real(double), intent(in) :: mu call msg_debug2 (D_ME_METHODS, "set_mu_uv_rcl") call set_mu_uv_rcl (mu) end subroutine rclwrap_set_mu_uv @ %def rclwrap_set_mu_uv @ Sets the IR pole parameterizations $\Delta_{IR}$ and $\Delta_2$. <>= public :: rclwrap_set_delta_ir <>= subroutine rclwrap_set_delta_ir (d, d2) real(double), intent(in) :: d, d2 call msg_debug2 (D_ME_METHODS, "set_delta_ir_rcl", & real(d, kind=default)) call msg_debug2 (D_ME_METHODS, "set_delta_ir_rcl", & real(d2, kind=default)) call set_delta_ir_rcl (d, d2) end subroutine rclwrap_set_delta_ir @ %def rclwrap_set_delta_ir @ <>= public :: rclwrap_set_mu_ir <>= subroutine rclwrap_set_mu_ir (mu) real(double), intent(in) :: mu call msg_debug2 (D_ME_METHODS, "set_mu_ir_rcl") call set_mu_ir_rcl (mu) end subroutine rclwrap_set_mu_ir @ %def rclwrap_set_mu_ir @ <>= public :: rclwrap_get_renormalization_scale <>= function rclwrap_get_renormalization_scale () result (mu) real(double) :: mu call msg_debug2 (D_ME_METHODS, "get_renormalization_scale_rcl") call get_renormalization_scale_rcl (mu) end function rclwrap_get_renormalization_scale @ %def rclwrap_get_renormalization_scale @ <>= public :: rclwrap_get_flavor_scheme <>= subroutine rclwrap_get_flavor_scheme (nf) integer, intent(out) :: nf call msg_debug2 (D_ME_METHODS, "get_flavour_scheme_rcl") call get_flavour_scheme_rcl (nf) end subroutine rclwrap_get_flavor_scheme @ %def rclwrap_get_flavor_scheme @ <>= public :: rclwrap_use_alpha0_scheme <>= subroutine rclwrap_use_alpha0_scheme (al0) real(double), intent(in), optional :: al0 call msg_debug2 (D_ME_METHODS, "use_alpha0_scheme_rcl") call use_alpha0_scheme_rcl (al0) end subroutine rclwrap_use_alpha0_scheme @ %def rclwrap_use_alpha0_scheme @ <>= public :: rclwrap_use_alphaz_scheme <>= subroutine rclwrap_use_alphaz_scheme (alz) real(double), intent(in), optional :: alz call msg_debug2 (D_ME_METHODS, "use_alphaz_scheme_rcl") call use_alphaz_scheme_rcl (alz) end subroutine rclwrap_use_alphaz_scheme @ %def rclwrap_use_alphaz_scheme @ <>= public :: rclwrap_set_complex_mass_scheme <>= subroutine rclwrap_set_complex_mass_scheme () call msg_debug2 (D_ME_METHODS, "set_complex_mass_scheme_rcl") call set_complex_mass_scheme_rcl () end subroutine rclwrap_set_complex_mass_scheme @ %def rclwrap_set_complex_mass_scheme @ <>= public :: rclwrap_set_resonant_particle <>= subroutine rclwrap_set_resonant_particle (pdg_id) integer, intent(in) :: pdg_id call msg_debug2 (D_ME_METHODS, "set_resonant_particle_rcl") call set_resonant_particle_rcl (char(get_recola_particle_string (pdg_id))) end subroutine rclwrap_set_resonant_particle @ %def rclwrap_set_resonant_particle @ <>= public :: rclwrap_switch_on_resonant_self_energies <>= subroutine rclwrap_switch_on_resonant_self_energies () call msg_debug2 (D_ME_METHODS, "switchon_resonant_selfenergies_rcl") call switchon_resonant_selfenergies_rcl () end subroutine rclwrap_switch_on_resonant_self_energies @ %def rclwrap_switch_on_resonant_self_energies @ <>= public :: rclwrap_switch_off_resonant_self_energies <>= subroutine rclwrap_switch_off_resonant_self_energies () call msg_debug2 (D_ME_METHODS, "switchoff_resonant_selfenergies_rcl") call switchoff_resonant_selfenergies_rcl () end subroutine rclwrap_switch_off_resonant_self_energies @ %def rclwrap_switch_off_resonant_self_energies @ <>= public :: rclwrap_set_draw_level_branches <>= subroutine rclwrap_set_draw_level_branches (n) integer, intent(in) :: n call msg_debug2 (D_ME_METHODS, "set_draw_level_branches_rcl") call set_draw_level_branches_rcl (n) end subroutine rclwrap_set_draw_level_branches @ %def rclwrap_set_draw_level_branches @ <>= public :: rclwrap_set_print_level_amplitude <>= subroutine rclwrap_set_print_level_amplitude (n) integer, intent(in) :: n call msg_debug2 (D_ME_METHODS, "set_print_level_amplitude_rcl") call set_print_level_amplitude_rcl (n) end subroutine rclwrap_set_print_level_amplitude @ %def rclwrap_set_print_level_amplitude @ <>= public :: rclwrap_set_print_level_squared_amplitude <>= subroutine rclwrap_set_print_level_squared_amplitude (n) integer, intent(in) :: n call msg_debug2 (D_ME_METHODS, "set_print_level_squared_amplitude_rcl") call set_print_level_squared_amplitude_rcl (n) end subroutine rclwrap_set_print_level_squared_amplitude @ %def rclwrap_set_print_level_squared_amplitude @ <>= public :: rclwrap_set_print_level_correlations <>= subroutine rclwrap_set_print_level_correlations (n) integer, intent(in) :: n call msg_debug2 (D_ME_METHODS, "set_print_level_correlations_rcl") call set_print_level_correlations_rcl (n) end subroutine rclwrap_set_print_level_correlations @ %def rclwrap_set_print_level_correlations @ <>= public :: rclwrap_set_print_level_RAM <>= subroutine rclwrap_set_print_level_RAM (n) integer, intent(in) :: n call msg_debug2 (D_ME_METHODS, "set_print_level_RAM_rcl") call set_print_level_RAM_rcl (n) end subroutine rclwrap_set_print_level_RAM @ %def rclwrap_set_print_level_RAM @ <>= public :: rclwrap_scale_coupling3 <>= subroutine rclwrap_scale_coupling3 (pdg_id1, pdg_id2, pdg_id3, factor) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3 complex(double), intent(in) :: factor call msg_debug2 (D_ME_METHODS, "scale_coupling3_rcl") call scale_coupling3_rcl (factor, char(get_recola_particle_string (pdg_id1)), & char(get_recola_particle_string (pdg_id2)), char(get_recola_particle_string (pdg_id3))) end subroutine rclwrap_scale_coupling3 @ %def rclwrap_scale_coupling3 @ <>= public :: rclwrap_scale_coupling4 <>= subroutine rclwrap_scale_coupling4 (pdg_id1, pdg_id2, pdg_id3, pdg_id4, factor) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3, pdg_id4 complex(double), intent(in) :: factor call msg_debug2 (D_ME_METHODS, "scale_coupling4_rcl") call scale_coupling4_rcl (factor, char(get_recola_particle_string (pdg_id1)), & char(get_recola_particle_string (pdg_id2)), char(get_recola_particle_string (pdg_id3)), & char(get_recola_particle_string (pdg_id4))) end subroutine rclwrap_scale_coupling4 @ %def rclwrap_scale_coupling4 @ <>= public :: rclwrap_switch_off_coupling3 <>= subroutine rclwrap_switch_off_coupling3 (pdg_id1, pdg_id2, pdg_id3) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3 call msg_debug2 (D_ME_METHODS, "switchoff_coupling3_rcl") call switchoff_coupling3_rcl (char(get_recola_particle_string (pdg_id1)), & char(get_recola_particle_string (pdg_id2)), char(get_recola_particle_string (pdg_id3))) end subroutine rclwrap_switch_off_coupling3 @ %def rclwrap_switch_off_coupling3 @ <>= public :: rclwrap_switch_off_coupling4 <>= subroutine rclwrap_switch_off_coupling4 (pdg_id1, pdg_id2, pdg_id3, pdg_id4) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3, pdg_id4 call msg_debug2 (D_ME_METHODS, "switchoff_coupling4_rcl") call switchoff_coupling4_rcl & (char(get_recola_particle_string (pdg_id1)), & char(get_recola_particle_string (pdg_id2)), & char(get_recola_particle_string (pdg_id3)), & char(get_recola_particle_string (pdg_id4))) end subroutine rclwrap_switch_off_coupling4 @ %def rclwrap_switch_off_coupling4 @ <>= public :: rclwrap_set_ifail <>= subroutine rclwrap_set_ifail (i) integer, intent(in) :: i call msg_debug2 (D_ME_METHODS, "set_ifail_rcl") call set_ifail_rcl (i) end subroutine rclwrap_set_ifail @ %def rclwrap_set_ifail @ <>= public :: rclwrap_get_ifail <>= subroutine rclwrap_get_ifail (i) integer, intent(out) :: i call msg_debug2 (D_ME_METHODS, "get_ifail_rcl") call get_ifail_rcl (i) end subroutine rclwrap_get_ifail @ %def rclwrap_get_ifail @ <>= public :: rclwrap_set_output_file <>= subroutine rclwrap_set_output_file (filename) character(len=*), intent(in) :: filename call msg_debug2 (D_ME_METHODS, "set_output_file_rcl") call set_output_file_rcl (filename) end subroutine rclwrap_set_output_file @ %def rclwrap_set_output_file @ <>= public :: rclwrap_set_gs_power <>= subroutine rclwrap_set_gs_power (id, gs_array) integer, intent(in) :: id integer, dimension(:,:), intent(in) :: gs_array call msg_debug2 (D_ME_METHODS, "set_gs_power_rcl") call set_gs_power_rcl (id, gs_array) end subroutine rclwrap_set_gs_power @ %def rclwrap_set_gs_power @ <>= public :: rclwrap_select_gs_power_born_amp <>= subroutine rclwrap_select_gs_power_born_amp (id, gs_power) integer, intent(in) :: id, gs_power call msg_debug2 (D_ME_METHODS, "select_gs_power_BornAmpl_rcl") call select_gs_power_BornAmpl_rcl (id, gs_power) end subroutine rclwrap_select_gs_power_born_amp @ %def rclwrap_select_gs_power_born_amp @ <>= public :: rclwrap_unselect_gs_power_born_amp <>= subroutine rclwrap_unselect_gs_power_born_amp (id, gs_power) integer, intent(in) :: id, gs_power call msg_debug2 (D_ME_METHODS, "unselect_gs_power_BornAmpl_rcl") call unselect_gs_power_BornAmpl_rcl (id, gs_power) end subroutine rclwrap_unselect_gs_power_born_amp @ %def rclwrap_unselect_gs_power_born_amp @ <>= public :: rclwrap_select_gs_power_loop_amp <>= subroutine rclwrap_select_gs_power_loop_amp (id, gs_power) integer, intent(in) :: id, gs_power call msg_debug2 (D_ME_METHODS, "select_gs_power_LoopAmpl_rcl") call select_gs_power_LoopAmpl_rcl (id, gs_power) end subroutine rclwrap_select_gs_power_loop_amp @ %def rclwrap_select_gs_power_loop_amp @ <>= public :: rclwrap_unselect_gs_power_loop_amp <>= subroutine rclwrap_unselect_gs_power_loop_amp (id, gs_power) integer, intent(in) :: id, gs_power call msg_debug2 (D_ME_METHODS, "unselect_gs_power_LoopAmpl_rcl") call unselect_gs_power_LoopAmpl_rcl (id, gs_power) end subroutine rclwrap_unselect_gs_power_loop_amp @ %def rclwrap_unselect_gs_power_loop_amp @ <>= public :: rclwrap_select_all_gs_powers_born_amp <>= subroutine rclwrap_select_all_gs_powers_born_amp (id) integer, intent(in) :: id call msg_debug2 (D_ME_METHODS, "select_all_gs_powers_BornAmpl_rcl") call select_all_gs_powers_BornAmpl_rcl (id) end subroutine rclwrap_select_all_gs_powers_born_amp @ %def rclwrap_select_all_gs_powers_born_amp @ <>= public :: rclwrap_unselect_all_gs_powers_loop_amp <>= subroutine rclwrap_unselect_all_gs_powers_loop_amp (id) integer, intent(in) :: id call msg_debug2 (D_ME_METHODS, "unselect_all_gs_powers_BornAmpl_rcl") call unselect_all_gs_powers_BornAmpl_rcl (id) end subroutine rclwrap_unselect_all_gs_powers_loop_amp @ %def rclwrap_unselect_all_gs_powers_loop_amp @ <>= public :: rclwrap_select_all_gs_powers_loop_amp <>= subroutine rclwrap_select_all_gs_powers_loop_amp (id) integer, intent(in) :: id call msg_debug2 (D_ME_METHODS, "select_all_gs_powers_LoopAmpl_rcl") call select_all_gs_powers_LoopAmpl_rcl (id) end subroutine rclwrap_select_all_gs_powers_loop_amp @ %def rclwrap_select_all_gs_powers_loop_amp @ <>= public :: rclwrap_unselect_all_gs_powers_born_amp <>= subroutine rclwrap_unselect_all_gs_powers_born_amp (id) integer, intent(in) :: id call msg_debug2 (D_ME_METHODS, "unselect_all_gs_powers_LoopAmpl_rcl") call unselect_all_gs_powers_LoopAmpl_rcl (id) end subroutine rclwrap_unselect_all_gs_powers_born_amp @ %def rclwrap_unselect_all_gs_powers_born_amp @ <>= public :: rclwrap_set_resonant_squared_momentum <>= subroutine rclwrap_set_resonant_squared_momentum (id, i_res, p2) integer, intent(in) :: id, i_res real(double), intent(in) :: p2 call msg_debug2 (D_ME_METHODS, "set_resonant_squared_momentum_rcl") call set_resonant_squared_momentum_rcl (id, i_res, p2) end subroutine rclwrap_set_resonant_squared_momentum @ %def rclwrap_set_resonant_squared_momentum @ <>= public :: rclwrap_compute_running_alpha_s <>= subroutine rclwrap_compute_running_alpha_s (Q, nf, n_loops) real(double), intent(in) :: Q integer, intent(in) :: nf, n_loops call msg_debug2 (D_ME_METHODS, "compute_running_alphas_rcl") call compute_running_alphas_rcl (Q, nf, n_loops) end subroutine rclwrap_compute_running_alpha_s @ %def rclwrap_compute_running_alpha_s @ <>= public :: rclwrap_set_dynamic_settings <>= subroutine rclwrap_set_dynamic_settings () call msg_debug2 (D_ME_METHODS, "set_dynamic_settings_rcl") call set_dynamic_settings_rcl (1) end subroutine rclwrap_set_dynamic_settings @ %def rclwrap_set_dynamic_settings @ <>= public :: rclwrap_rescale_process <>= subroutine rclwrap_rescale_process (id, order, sqme) integer, intent(in) :: id character(len=*), intent(in) :: order real(double), dimension(0:1), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "rescale_process_rcl") call rescale_process_rcl (id, order, sqme) end subroutine rclwrap_rescale_process @ %def rclwrap_rescale_process @ <>= public :: rclwrap_get_polarized_squared_amplitude <>= subroutine rclwrap_get_polarized_squared_amplitude (id, & alphas_power, order, hel, sqme) integer, intent(in) :: id, alphas_power character(len=*), intent(in) :: order integer, dimension(:), intent(in) :: hel real(double), intent(out) :: sqme call msg_debug2 (D_ME_METHODS, "get_polarized_squared_amplitude_rcl") call get_polarized_squared_amplitude_rcl (id, alphas_power, & order, hel, sqme) end subroutine rclwrap_get_polarized_squared_amplitude @ %def rclwrap_get_polarized_squared_amplitude @ <>= public :: rclwrap_compute_color_correlation <>= subroutine rclwrap_compute_color_correlation (id, p, & i1, i2, sqme) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p integer, intent(in) :: i1, i2 real(double), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "compute_colour_correlation_rcl") call compute_colour_correlation_rcl (id, p, i1, i2, sqme) end subroutine rclwrap_compute_color_correlation @ %def rclwrap_compute_color_correlation @ <>= public :: rclwrap_compute_all_color_correlations <>= subroutine rclwrap_compute_all_color_correlations (id, p) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p call msg_debug2 (D_ME_METHODS, "compute_all_colour_correlations_rcl") call compute_all_colour_correlations_rcl (id, p) end subroutine rclwrap_compute_all_color_correlations @ %def rclwrap_compute_all_color_correlations @ <>= public :: rclwrap_rescale_color_correlation <>= subroutine rclwrap_rescale_color_correlation (id, i1, i2, sqme) integer, intent(in) :: id, i1, i2 real(double), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "rescale_colour_correlation_rcl") call rescale_colour_correlation_rcl (id, i1, i2, sqme) end subroutine rclwrap_rescale_color_correlation @ %def rclwrap_rescale_color_correlation @ <>= public :: rclwrap_rescale_all_color_correlations <>= subroutine rclwrap_rescale_all_color_correlations (id) integer, intent(in) :: id call msg_debug2 (D_ME_METHODS, "rescale_all_colour_correlations_rcl") call rescale_all_colour_correlations_rcl (id) end subroutine rclwrap_rescale_all_color_correlations @ %def rclwrap_rescale_all_color_correlations @ <>= public :: rclwrap_get_color_correlation <>= subroutine rclwrap_get_color_correlation (id, alphas_power, i1, i2, sqme) integer, intent(in) :: id, alphas_power, i1, i2 real(double), intent(out) :: sqme call msg_debug2 (D_ME_METHODS, "get_colour_correlation_rcl") call get_colour_correlation_rcl (id, alphas_power, i1, i2, sqme) end subroutine rclwrap_get_color_correlation @ %def rclwrap_get_color_correlation @ <>= public :: rclwrap_compute_spin_correlation <>= subroutine rclwrap_compute_spin_correlation (id, p, i_photon, pol, sqme) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p integer, intent(in) :: i_photon complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "compute_spin_correlation_rcl") call compute_spin_correlation_rcl (id, p, i_photon, pol, sqme) end subroutine rclwrap_compute_spin_correlation @ %def rclwrap_compute_spin_correlation @ <>= public :: rclwrap_rescale_spin_correlation <>= subroutine rclwrap_rescale_spin_correlation (id, i_photon, pol, sqme) integer, intent(in) :: id, i_photon complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "rescale_spin_correlation_rcl") call rescale_spin_correlation_rcl (id, i_photon, pol, sqme) end subroutine rclwrap_rescale_spin_correlation @ %def rclwrap_rescale_spin_correlation @ <>= public :: rclwrap_get_spin_correlation <>= subroutine rclwrap_get_spin_correlation (id, alphas_power, sqme) integer, intent(in) :: id, alphas_power real(double), intent(out) :: sqme call msg_debug2 (D_ME_METHODS, "get_spin_correlation_rcl") call get_spin_correlation_rcl (id, alphas_power, sqme) end subroutine rclwrap_get_spin_correlation @ %def rclwrap_get_spin_correlation @ <>= public :: rclwrap_compute_spin_color_correlation <>= subroutine rclwrap_compute_spin_color_correlation (id, p, & i_gluon, i_spectator, pol, sqme) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p integer, intent(in) :: i_gluon, i_spectator complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "compute_spin_colour_correlation_rcl") call compute_spin_colour_correlation_rcl (id, p, & i_gluon, i_spectator, pol, sqme) end subroutine rclwrap_compute_spin_color_correlation @ %def rclwrap_compute_spin_color_correlation @ <>= public :: rclwrap_rescale_spin_color_correlation <>= subroutine rclwrap_rescale_spin_color_correlation (id, i_gluon, & i_spectator, pol, sqme) integer, intent(in) :: id, i_gluon, i_spectator complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme call msg_debug2 (D_ME_METHODS, "rescale_spin_colour_correlation_rcl") call rescale_spin_colour_correlation_rcl (id, i_gluon, & i_spectator, pol, sqme) end subroutine rclwrap_rescale_spin_color_correlation @ %def rclwrap_rescale_spin_color_correlation @ <>= public :: rclwrap_get_spin_color_correlation <>= subroutine rclwrap_get_spin_color_correlation (id, alphas_power, & i_gluon, i_spectator, sqme) integer, intent(in) :: id, alphas_power, i_gluon, i_spectator real(double), intent(out) :: sqme call msg_debug2 (D_ME_METHODS, "get_spin_colour_correlation_rcl") call get_spin_colour_correlation_rcl (id, alphas_power, & i_gluon, i_spectator, sqme) end subroutine rclwrap_get_spin_color_correlation @ %def rclwrap_get_spin_color_correlation @ <>= public :: rclwrap_get_momenta <>= subroutine rclwrap_get_momenta (id, p) integer, intent(in) :: id real(double), dimension(:,:), intent(out) :: p call msg_debug2 (D_ME_METHODS, "get_momenta_rcl") call get_momenta_rcl (id, p) end subroutine rclwrap_get_momenta @ %def rclwrap_get_momenta @ The reset routine is essential. But note that it doesn't reset the Recola parameters, just the processes. For LOL, Recola's reset routine crashes the program if there was no process before. So, rather reset indirectly via the controller. <>= public :: rclwrap_reset_recola <>= subroutine rclwrap_reset_recola call msg_debug (D_ME_METHODS, "rclwrap_reset_recola") call rcl_controller%reset () end subroutine rclwrap_reset_recola @ %def rclwrap_reset_recola @ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \section{Recola dummy replacement module} <<[[recola_wrapper_dummy.f90]]>>= <> module recola_wrapper use kinds <> <> <> <> contains <> end module recola_wrapper @ %def recola_wrapper_dummy @ <>= public :: rclwrap_is_active <>= logical, parameter :: rclwrap_is_active = .false. @ %def rclwrap_is_active @ <>= public :: get_recola_particle_string <>= elemental function get_recola_particle_string (pdg) result (name) type(string_t) :: name integer, intent(in) :: pdg name = var_str ("?") end function get_recola_particle_string @ %def get_recola_paritcle_string @ <>= public :: rclwrap_get_new_recola_id <>= subroutine rclwrap_get_new_recola_id (id) integer, intent(out) :: id id = 0 end subroutine rclwrap_get_new_recola_id @ %def rclwrap_get_new_recola_id @ <>= - public :: rclwrap_get_n_processes + public :: rclwrap_get_current_recola_id <>= - function rclwrap_get_n_processes () result (n) + function rclwrap_get_current_recola_id () result (n) integer :: n n = 0 - end function rclwrap_get_n_processes + end function rclwrap_get_current_recola_id -@ %def rclwrap_get_n_processes +@ %def rclwrap_get_current_recola_id @ <>= - public :: rclwrap_define_process -<>= - subroutine rclwrap_define_process (id, process_string, order) - integer, intent(in) :: id - type(string_t), intent(in) :: process_string - type(string_t), intent(in) :: order - end subroutine rclwrap_define_process - -@ %def rclwrap_define_process -@ -<>= public :: rclwrap_request_generate_processes <>= subroutine rclwrap_request_generate_processes () end subroutine rclwrap_request_generate_processes @ %def rclwrap_request_generate_processes @ <>= + public :: rclwrap_add_process +<>= + subroutine rclwrap_add_process (id, process_string, order) + integer, intent(in) :: id + type(string_t), intent(in) :: process_string, order + end subroutine rclwrap_add_process + +@ %def rclwrap_add_process +@ +<>= + public :: rclwrap_define_processes +<>= + subroutine rclwrap_define_processes () + end subroutine rclwrap_define_processes + +@ %def rclwrap_define_processes +@ +<>= public :: rclwrap_generate_processes <>= subroutine rclwrap_generate_processes () end subroutine rclwrap_generate_processes @ %def rclwrap_generate_processes @ <>= public :: rclwrap_compute_process <>= subroutine rclwrap_compute_process (id, p, order, sqme) integer, intent(in) :: id real(double), intent(in), dimension(:,:) :: p character(len=*), intent(in) :: order real(double), intent(out), dimension(0:1), optional :: sqme end subroutine rclwrap_compute_process @ %def rclwrap_compute_process @ <>= public :: rclwrap_get_amplitude <>= subroutine rclwrap_get_amplitude (id, g_power, order, col, hel, amp) integer, intent(in) :: id, g_power character(len=*), intent(in) :: order integer, dimension(:), intent(in) :: col, hel complex(double), intent(out) :: amp end subroutine rclwrap_get_amplitude @ %def rclwrap_get_amplitude @ <>= public :: rclwrap_get_squared_amplitude <>= subroutine rclwrap_get_squared_amplitude (id, alphas_power, order, sqme) integer, intent(in) :: id, alphas_power character(len=*), intent(in) :: order real(double), intent(out) :: sqme end subroutine rclwrap_get_squared_amplitude @ %def rclwrap_get_squared_amplitude @ <>= public :: rclwrap_set_pole_mass <>= subroutine rclwrap_set_pole_mass (pdg_id, mass, width) integer, intent(in) :: pdg_id real(double), intent(in) :: mass, width end subroutine rclwrap_set_pole_mass @ %def rclwrap_set_pole_mass @ <>= public :: rclwrap_set_onshell_mass <>= subroutine rclwrap_set_onshell_mass (pdg_id, mass, width) integer, intent(in) :: pdg_id real(double), intent(in) :: mass, width end subroutine rclwrap_set_onshell_mass @ %def rclwrap_set_onshell_mass @ <>= public :: rclwrap_use_gfermi_scheme <>= subroutine rclwrap_use_gfermi_scheme (gf) real(double), intent(in), optional :: gf end subroutine rclwrap_use_gfermi_scheme @ %def rclwrap_use_gfermi_scheme @ <>= public :: rclwrap_set_light_fermions <>= subroutine rclwrap_set_light_fermions (m) real(double), intent(in) :: m end subroutine rclwrap_set_light_fermions @ %def rclwrap_set_light_fermions @ <>= public :: rclwrap_set_light_fermion <>= subroutine rclwrap_set_light_fermion (pdg_id) integer, intent(in) :: pdg_id end subroutine rclwrap_set_light_fermion @ %def rclwrap_set_light_fermion @ <>= public :: rclwrap_unset_light_fermion <>= subroutine rclwrap_unset_light_fermion (pdg_id) integer, intent(in) :: pdg_id end subroutine rclwrap_unset_light_fermion @ %def rclwrap_unset_light_fermion @ <>= public :: rclwrap_set_onshell_scheme <>= subroutine rclwrap_set_onshell_scheme end subroutine rclwrap_set_onshell_scheme @ %def rclwrap_set_onshell_scheme @ <>= public :: rclwrap_set_alpha_s <>= subroutine rclwrap_set_alpha_s (alpha_s, mu, nf) real(double), intent(in) :: alpha_s, mu integer, intent(in) :: nf end subroutine rclwrap_set_alpha_s @ %def rclwrap_set_alpha_s @ <>= public :: rclwrap_get_alpha_s <>= function rclwrap_get_alpha_s () result (alpha_s) real(double) :: alpha_s end function rclwrap_get_alpha_s @ %def rclwrap_get_alpha_s @ <>= public :: rclwrap_get_helicity_configurations <>= subroutine rclwrap_get_helicity_configurations (id, hel) integer, intent(in) :: id integer, intent(inout), dimension(:,:), allocatable :: hel end subroutine rclwrap_get_helicity_configurations @ %def rclwrap_get_helicity_configurations @ <>= public :: rclwrap_get_color_configurations <>= subroutine rclwrap_get_color_configurations (id, col) integer, intent(in) :: id integer, intent(out), dimension(:,:), allocatable :: col end subroutine rclwrap_get_color_configurations @ %def rclwrap_get_color_configurations @ <>= public :: rclwrap_use_dim_reg_soft <>= subroutine rclwrap_use_dim_reg_soft () end subroutine rclwrap_use_dim_reg_soft @ %def rclwrap_use_dim_reg_soft @ <>= public :: rclwrap_use_mass_reg_soft <>= subroutine rclwrap_use_mass_reg_soft (m) real(double), intent(in) :: m end subroutine rclwrap_use_mass_reg_soft @ %def rclwrap_use_mass_reg_soft @ <>= public :: rclwrap_set_delta_uv <>= subroutine rclwrap_set_delta_uv (d) real(double), intent(in) :: d end subroutine rclwrap_set_delta_uv @ %def rclwrap_set_delta_uv @ <>= public :: rclwrap_set_mu_uv <>= subroutine rclwrap_set_mu_uv (mu) real(double), intent(in) :: mu end subroutine rclwrap_set_mu_uv @ %def rclwrap_set_mu_uv @ <>= public :: rclwrap_set_delta_ir <>= subroutine rclwrap_set_delta_ir (d, d2) real(double), intent(in) :: d, d2 end subroutine rclwrap_set_delta_ir @ %def rclwrap_set_delta_ir @ <>= public :: rclwrap_set_mu_ir <>= subroutine rclwrap_set_mu_ir (mu) real(double), intent(in) :: mu end subroutine rclwrap_set_mu_ir @ %def rclwrap_set_mu_ir @ <>= public :: rclwrap_get_renormalization_scale <>= function rclwrap_get_renormalization_scale () result (mu) real(double) :: mu end function rclwrap_get_renormalization_scale @ %def rclwrap_get_renormalization_scale @ <>= public :: rclwrap_get_flavor_scheme <>= subroutine rclwrap_get_flavor_scheme (nf) integer, intent(out) :: nf end subroutine rclwrap_get_flavor_scheme @ %def rclwrap_get_flavor_scheme @ <>= public :: rclwrap_use_alpha0_scheme <>= subroutine rclwrap_use_alpha0_scheme (al0) real(double), intent(in), optional :: al0 end subroutine rclwrap_use_alpha0_scheme @ %def rclwrap_use_alpha0_scheme @ <>= public :: rclwrap_use_alphaz_scheme <>= subroutine rclwrap_use_alphaz_scheme (alz) real(double), intent(in), optional :: alz end subroutine rclwrap_use_alphaz_scheme @ %def rclwrap_use_alphaz_scheme @ <>= public :: rclwrap_set_complex_mass_scheme <>= subroutine rclwrap_set_complex_mass_scheme () end subroutine rclwrap_set_complex_mass_scheme @ %def rclwrap_set_complex_mass_scheme @ <>= public :: rclwrap_set_resonant_particle <>= subroutine rclwrap_set_resonant_particle (pdg_id) integer, intent(in) :: pdg_id end subroutine rclwrap_set_resonant_particle @ %def rclwrap_set_resonant_particle @ <>= public :: rclwrap_switch_on_resonant_self_energies <>= subroutine rclwrap_switch_on_resonant_self_energies () end subroutine rclwrap_switch_on_resonant_self_energies @ %def rclwrap_switch_on_resonant_self_energies @ <>= public :: rclwrap_switch_off_resonant_self_energies <>= subroutine rclwrap_switch_off_resonant_self_energies () end subroutine rclwrap_switch_off_resonant_self_energies @ %def rclwrap_switch_off_resonant_self_energies @ <>= public :: rclwrap_set_draw_level_branches <>= subroutine rclwrap_set_draw_level_branches (n) integer, intent(in) :: n end subroutine rclwrap_set_draw_level_branches @ %def rclwrap_set_draw_level_branches @ <>= public :: rclwrap_set_print_level_amplitude <>= subroutine rclwrap_set_print_level_amplitude (n) integer, intent(in) :: n end subroutine rclwrap_set_print_level_amplitude @ %def rclwrap_set_print_level_amplitude @ <>= public :: rclwrap_set_print_level_squared_amplitude <>= subroutine rclwrap_set_print_level_squared_amplitude (n) integer, intent(in) :: n end subroutine rclwrap_set_print_level_squared_amplitude @ %def rclwrap_set_print_level_squared_amplitude @ <>= public :: rclwrap_set_print_level_correlations <>= subroutine rclwrap_set_print_level_correlations (n) integer, intent(in) :: n end subroutine rclwrap_set_print_level_correlations @ %def rclwrap_set_print_level_correlations @ <>= public :: rclwrap_set_print_level_RAM <>= subroutine rclwrap_set_print_level_RAM (n) integer, intent(in) :: n end subroutine rclwrap_set_print_level_RAM @ %def rclwrap_set_print_level_RAM @ <>= public :: rclwrap_scale_coupling3 <>= subroutine rclwrap_scale_coupling3 (pdg_id1, pdg_id2, pdg_id3, factor) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3 complex(double), intent(in) :: factor end subroutine rclwrap_scale_coupling3 @ %def rclwrap_scale_coupling3 @ <>= public :: rclwrap_scale_coupling4 <>= subroutine rclwrap_scale_coupling4 (pdg_id1, pdg_id2, pdg_id3, pdg_id4, factor) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3, pdg_id4 complex(double), intent(in) :: factor end subroutine rclwrap_scale_coupling4 @ %def rclwrap_scale_coupling4 @ <>= public :: rclwrap_switch_off_coupling3 <>= subroutine rclwrap_switch_off_coupling3 (pdg_id1, pdg_id2, pdg_id3) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3 end subroutine rclwrap_switch_off_coupling3 @ %def rclwrap_switch_off_coupling3 @ <>= public :: rclwrap_switch_off_coupling4 <>= subroutine rclwrap_switch_off_coupling4 (pdg_id1, pdg_id2, pdg_id3, pdg_id4) integer, intent(in) :: pdg_id1, pdg_id2, pdg_id3, pdg_id4 end subroutine rclwrap_switch_off_coupling4 @ %def rclwrap_switch_off_coupling4 @ <>= public :: rclwrap_set_ifail <>= subroutine rclwrap_set_ifail (i) integer, intent(in) :: i end subroutine rclwrap_set_ifail @ %def rclwrap_set_ifail @ <>= public :: rclwrap_get_ifail <>= subroutine rclwrap_get_ifail (i) integer, intent(out) :: i end subroutine rclwrap_get_ifail @ %def rclwrap_get_ifail @ <>= public :: rclwrap_set_output_file <>= subroutine rclwrap_set_output_file (filename) character(len=*), intent(in) :: filename end subroutine rclwrap_set_output_file @ %def rclwrap_set_output_file @ <>= public :: rclwrap_set_gs_power <>= subroutine rclwrap_set_gs_power (id, gs_array) integer, intent(in) :: id integer, dimension(:,:), intent(in) :: gs_array end subroutine rclwrap_set_gs_power @ %def rclwrap_set_gs_power @ <>= public :: rclwrap_select_gs_power_born_amp <>= subroutine rclwrap_select_gs_power_born_amp (id, gs_power) integer, intent(in) :: id, gs_power end subroutine rclwrap_select_gs_power_born_amp @ %def rclwrap_select_gs_power_born_amp @ <>= public :: rclwrap_unselect_gs_power_born_amp <>= subroutine rclwrap_unselect_gs_power_born_amp (id, gs_power) integer, intent(in) :: id, gs_power end subroutine rclwrap_unselect_gs_power_born_amp @ %def rclwrap_unselect_gs_power_born_amp @ <>= public :: rclwrap_select_gs_power_loop_amp <>= subroutine rclwrap_select_gs_power_loop_amp (id, gs_power) integer, intent(in) :: id, gs_power end subroutine rclwrap_select_gs_power_loop_amp @ %def rclwrap_select_gs_power_loop_amp @ <>= public :: rclwrap_unselect_gs_power_loop_amp <>= subroutine rclwrap_unselect_gs_power_loop_amp (id, gs_power) integer, intent(in) :: id, gs_power end subroutine rclwrap_unselect_gs_power_loop_amp @ %def rclwrap_unselect_gs_power_loop_amp @ <>= public :: rclwrap_select_all_gs_powers_born_amp <>= subroutine rclwrap_select_all_gs_powers_born_amp (id) integer, intent(in) :: id end subroutine rclwrap_select_all_gs_powers_born_amp @ %def rclwrap_select_all_gs_powers_born_amp @ <>= public :: rclwrap_unselect_all_gs_powers_loop_amp <>= subroutine rclwrap_unselect_all_gs_powers_loop_amp (id) integer, intent(in) :: id end subroutine rclwrap_unselect_all_gs_powers_loop_amp @ %def rclwrap_unselect_all_gs_powers_loop_amp @ <>= public :: rclwrap_select_all_gs_powers_loop_amp <>= subroutine rclwrap_select_all_gs_powers_loop_amp (id) integer, intent(in) :: id end subroutine rclwrap_select_all_gs_powers_loop_amp @ %def rclwrap_select_all_gs_powers_loop_amp @ <>= public :: rclwrap_unselect_all_gs_powers_born_amp <>= subroutine rclwrap_unselect_all_gs_powers_born_amp (id) integer, intent(in) :: id end subroutine rclwrap_unselect_all_gs_powers_born_amp @ %def rclwrap_unselect_all_gs_powers_born_amp @ <>= public :: rclwrap_set_resonant_squared_momentum <>= subroutine rclwrap_set_resonant_squared_momentum (id, i_res, p2) integer, intent(in) :: id, i_res real(double), intent(in) :: p2 end subroutine rclwrap_set_resonant_squared_momentum @ %def rclwrap_set_resonant_squared_momentum @ <>= public :: rclwrap_compute_running_alpha_s <>= subroutine rclwrap_compute_running_alpha_s (Q, nf, n_loops) real(double), intent(in) :: Q integer, intent(in) :: nf, n_loops end subroutine rclwrap_compute_running_alpha_s @ %def rclwrap_compute_running_alpha_s @ <>= public :: rclwrap_set_dynamic_settings <>= subroutine rclwrap_set_dynamic_settings () end subroutine rclwrap_set_dynamic_settings @ %def rclwrap_set_dynamic_settings @ <>= public :: rclwrap_rescale_process <>= subroutine rclwrap_rescale_process (id, order, sqme) integer, intent(in) :: id character(len=*), intent(in) :: order real(double), dimension(0:1), intent(out), optional :: sqme end subroutine rclwrap_rescale_process @ %def rclwrap_rescale_process @ <>= public :: rclwrap_get_polarized_squared_amplitude <>= subroutine rclwrap_get_polarized_squared_amplitude (id, & alphas_power, order, hel, sqme) integer, intent(in) :: id, alphas_power character(len=*), intent(in) :: order integer, dimension(:), intent(in) :: hel real(double), intent(out) :: sqme end subroutine rclwrap_get_polarized_squared_amplitude @ %def rclwrap_get_polarized_squared_amplitude @ <>= public :: rclwrap_compute_color_correlation <>= subroutine rclwrap_compute_color_correlation (id, p, & i1, i2, sqme) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p integer, intent(in) :: i1, i2 real(double), intent(out), optional :: sqme end subroutine rclwrap_compute_color_correlation @ %def rclwrap_compute_color_correlation @ <>= public :: rclwrap_compute_all_color_correlations <>= subroutine rclwrap_compute_all_color_correlations (id, p) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p end subroutine rclwrap_compute_all_color_correlations @ %def rclwrap_compute_all_color_correlations @ <>= public :: rclwrap_rescale_color_correlation <>= subroutine rclwrap_rescale_color_correlation (id, i1, i2, sqme) integer, intent(in) :: id, i1, i2 real(double), intent(out), optional :: sqme end subroutine rclwrap_rescale_color_correlation @ %def rclwrap_rescale_color_correlation @ <>= public :: rclwrap_rescale_all_color_correlations <>= subroutine rclwrap_rescale_all_color_correlations (id) integer, intent(in) :: id end subroutine rclwrap_rescale_all_color_correlations @ %def rclwrap_rescale_all_color_correlations @ <>= public :: rclwrap_get_color_correlation <>= subroutine rclwrap_get_color_correlation (id, alphas_power, i1, i2, sqme) integer, intent(in) :: id, alphas_power, i1, i2 real(double), intent(out) :: sqme end subroutine rclwrap_get_color_correlation @ %def rclwrap_get_color_correlation @ <>= public :: rclwrap_compute_spin_correlation <>= subroutine rclwrap_compute_spin_correlation (id, p, i_photon, pol, sqme) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p integer, intent(in) :: i_photon complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme end subroutine rclwrap_compute_spin_correlation @ %def rclwrap_compute_spin_correlation @ <>= public :: rclwrap_rescale_spin_correlation <>= subroutine rclwrap_rescale_spin_correlation (id, i_photon, pol, sqme) integer, intent(in) :: id, i_photon complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme end subroutine rclwrap_rescale_spin_correlation @ %def rclwrap_rescale_spin_correlation @ <>= public :: rclwrap_get_spin_correlation <>= subroutine rclwrap_get_spin_correlation (id, alphas_power, sqme) integer, intent(in) :: id, alphas_power real(double), intent(out) :: sqme end subroutine rclwrap_get_spin_correlation @ %def rclwrap_get_spin_correlation @ <>= public :: rclwrap_compute_spin_color_correlation <>= subroutine rclwrap_compute_spin_color_correlation (id, p, & i_gluon, i_spectator, pol, sqme) integer, intent(in) :: id real(double), dimension(:,:), intent(in) :: p integer, intent(in) :: i_gluon, i_spectator complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme end subroutine rclwrap_compute_spin_color_correlation @ %def rclwrap_compute_spin_color_correlation @ <>= public :: rclwrap_rescale_spin_color_correlation <>= subroutine rclwrap_rescale_spin_color_correlation (id, i_gluon, & i_spectator, pol, sqme) integer, intent(in) :: id, i_gluon, i_spectator complex(double), dimension(:), intent(in) :: pol real(double), intent(out), optional :: sqme end subroutine rclwrap_rescale_spin_color_correlation @ %def rclwrap_rescale_spin_color_correlation @ <>= public :: rclwrap_get_spin_color_correlation <>= subroutine rclwrap_get_spin_color_correlation (id, alphas_power, & i_gluon, i_spectator, sqme) integer, intent(in) :: id, alphas_power, i_gluon, i_spectator real(double), intent(out) :: sqme end subroutine rclwrap_get_spin_color_correlation @ %def rclwrap_get_spin_color_correlation @ <>= public :: rclwrap_get_momenta <>= subroutine rclwrap_get_momenta (id, p) integer, intent(in) :: id real(double), dimension(:,:), intent(out) :: p end subroutine rclwrap_get_momenta @ %def rclwrap_get_momenta @ <>= public :: rclwrap_reset_recola <>= subroutine rclwrap_reset_recola end subroutine rclwrap_reset_recola @ %def rclwrap_reset_recola @ \clearpage %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \section{Recola Core} The recola core object and auxiliary types and objects. <<[[prc_recola.f90]]>>= <> module prc_recola use recola_wrapper !NODEP! use kinds use constants, only: pi, zero <> use string_utils, only: str use system_defs, only: TAB use diagnostics use io_units use lorentz use physics_defs use variables, only: var_list_t use os_interface, only: os_data_t use sm_qcd, only: qcd_t use model_data, only: model_data_t use prc_core, only: prc_core_state_t use prc_core_def, only: prc_core_driver_t, prc_core_def_t use prc_user_defined use process_libraries, only: process_library_t <> <> <> <> <> contains <> end module prc_recola @ %def prc_recola @ \subsection{Sanity check} Checks the [[rclwrap_is_active]] flag and aborts the program if the dummy is used. <>= public :: abort_if_recola_not_active <>= subroutine abort_if_recola_not_active () if (.not. rclwrap_is_active) call msg_fatal ("You want to use Recola, ", & [var_str("but either the compiler with which Whizard has been build "), & var_str("is not supported by it, or you have not linked Recola "), & var_str("correctly to Whizard. Either reconfigure Whizard with a path to "), & var_str("a valid Recola installation (for details consult the manual), "), & var_str("or choose a different matrix-element method.")]) end subroutine abort_if_recola_not_active @ %def abort_if_recola_not_active @ \subsection{Process definition} When defining a RECOLA process, we store the process-specific flags and parameters. <>= public :: recola_def_t <>= type, extends (user_defined_def_t) :: recola_def_t type(string_t) :: suffix type(string_t) :: order integer :: alpha_power = 0 integer :: alphas_power = 0 contains <> end type recola_def_t @ %def recola_def_t @ <>= procedure, nopass :: type_string => recola_def_type_string <>= function recola_def_type_string () result (string) type(string_t) :: string string = "recola" end function recola_def_type_string @ %def recola_def_type_string @ Not implemented yet. <>= procedure :: write => recola_def_write <>= subroutine recola_def_write (object, unit) class(recola_def_t), intent(in) :: object integer, intent(in) :: unit end subroutine recola_def_write @ %def recola_def_write @ <>= procedure :: read => recola_def_read <>= subroutine recola_def_read (object, unit) class(recola_def_t), intent(out) :: object integer, intent(in) :: unit end subroutine recola_def_read @ %def recola_def_read @ The initializer has the responsibility to store all process- and method-specific parameters, such that they can be used later by the writer and by the driver for this process. Also, it allocates the writer. For RECOLA, the writer (i) creates full-fledged \oMega\ matrix element code which we need for the interface. (TODO: I guess that we don't need the matrix element, just the management part. Maybe modify the \oMega\ call such that no matrix-element code is written.) (ii) registers the process definition with the RECOLA library which has been linked. The latter task does not involve external code. Note that all management stuff is taken care of by the base type(s) methods. Here, we introduce only RECOLA-specific procedures, in addition. The NLO flag is true only for virtual matrix elements. <>= procedure :: init => recola_def_init <>= subroutine recola_def_init (object, basename, model_name, & prt_in, prt_out, nlo_type, alpha_power, alphas_power) class(recola_def_t), intent(inout) :: object type(string_t), intent(in) :: basename, model_name type(string_t), dimension(:), intent(in) :: prt_in, prt_out integer, intent(in) :: nlo_type integer, intent(in) :: alpha_power integer, intent(in) :: alphas_power call msg_debug (D_ME_METHODS, "recola_def_init: " & // char (basename) // ", nlo_type", nlo_type) object%basename = basename object%alpha_power = alpha_power object%alphas_power = alphas_power allocate (recola_writer_t :: object%writer) select case (nlo_type) case (BORN) object%suffix = '_BORN' object%order = "LO" case (NLO_REAL) object%suffix = '_REAL' object%order = "LO" case (NLO_VIRTUAL) object%suffix = '_LOOP' object%order = "NLO" case (NLO_SUBTRACTION) object%suffix = '_SUB' object%order = "LO" case (NLO_MISMATCH) object%suffix = '_MISMATCH' object%order = "LO" case (NLO_DGLAP) object%suffix = '_DGLAP' object%order = "LO" end select select type (writer => object%writer) class is (recola_writer_t) call writer%init (model_name, prt_in, prt_out) - call writer%set_id (basename) + call writer%set_id (basename // object%suffix) call writer%set_order (object%order) call writer%set_coupling_powers (object%alpha_power, object%alphas_power) end select end subroutine recola_def_init @ %def recola_def_init @ \subsection{Writer object} The RECOLA writer takes the additional resposibility of transferring process information to RECOLA. <>= type, extends (prc_user_defined_writer_t) :: recola_writer_t private type(string_t) :: id type(string_t) :: order integer :: alpha_power = 0 integer :: alphas_power = 0 contains <> end type recola_writer_t @ %def recola_writer_t @ <>= procedure, nopass :: type_name => recola_writer_type_name <>= function recola_writer_type_name () result (string) type(string_t) :: string string = "recola" end function recola_writer_type_name @ %def recola_writer_type_name @ Set the process ID string as used by WHIZARD. <>= procedure :: set_id => recola_writer_set_id <>= subroutine recola_writer_set_id (writer, id) class(recola_writer_t), intent(inout) :: writer type(string_t), intent(in) :: id call msg_debug2 (D_ME_METHODS, "Recola writer: id = " // char (id)) writer%id = id end subroutine recola_writer_set_id @ %def recola_writer_set_id @ Set the NLO flag. <>= procedure :: set_order => recola_writer_set_order <>= subroutine recola_writer_set_order (writer, order) class(recola_writer_t), intent(inout) :: writer type(string_t), intent(in) :: order call msg_debug2 (D_ME_METHODS, "Recola writer: order = " // char (order)) writer%order = order end subroutine recola_writer_set_order @ %def recola_writer_set_order @ Set coupling powers. <>= procedure :: set_coupling_powers => recola_writer_set_coupling_powers <>= subroutine recola_writer_set_coupling_powers (writer, alpha_power, alphas_power) class(recola_writer_t), intent(inout) :: writer integer, intent(in) :: alpha_power integer, intent(in) :: alphas_power call msg_debug2 (D_ME_METHODS, "Recola writer: alphas_power", alphas_power) call msg_debug2 (D_ME_METHODS, "Recola writer: alpha_power", alpha_power) writer%alpha_power = alpha_power writer%alphas_power = alphas_power end subroutine recola_writer_set_coupling_powers @ %def recola_writer_set_coupling_powers @ The Makefile code contains all of the code that the [[user_defined]] base method generates, plus an extra clause that extracts a shorthand listing of all flavor combinations for the current process. This list is required by [[make source]], so it can be read and used for declaring the RECOLA processes. There is one glitch here: we use the component-specific source file but write a flavor list for the process, without component extension. That is, we must not have more than one component at this stage. NB: We might actually extend \oMega\ to produce this shorthand listing. <>= procedure :: write_makefile_code => recola_writer_write_makefile_code <>= function flv_file_name (id) type(string_t), intent(in) :: id type(string_t) :: flv_file_name flv_file_name = id // ".flv.dat" end function flv_file_name subroutine recola_writer_write_makefile_code & (writer, unit, id, os_data, verbose, testflag) class(recola_writer_t), intent(in) :: writer integer, intent(in) :: unit type(string_t), intent(in) :: id type(os_data_t), intent(in) :: os_data logical, intent(in) :: verbose logical, intent(in), optional :: testflag type(string_t) :: src_file type(string_t) :: flv_file call writer%base_write_makefile_code (unit, id, os_data, verbose, testflag) src_file = trim (char(id)) // ".f90" flv_file = flv_file_name (writer%id) write (unit, *) write (unit, "(5A)") "# Flavor state listing for RECOLA process generation" write (unit, "(5A)") char (flv_file), ": ", char (src_file) if (verbose) then write (unit, "(5A)", advance="no") TAB else write (unit, "(5A)") TAB, '@echo " MAKE ', char (flv_file), '"' write (unit, "(5A)", advance="no") TAB, "@" end if write (unit, "(5A)") & "grep 'data table_flavor_states' $< ", & "| sed -e 's/.*\/\(.*\)\/.*/\1/' -e 's/,//g' > $@" write (unit, "(5A)") "SOURCES += ", char (flv_file) write (unit, "(5A)") "CLEAN_SOURCES += ", char (flv_file) end subroutine recola_writer_write_makefile_code @ %def recola_writer_write_makefile_code @ To communicate the process definition to RECOLA, we must know the following: the process definition, expanded in terms of flavor states, and the process order (LO/NLO). We will ask for a new numeric ID, create a process string using RECOLA conventions, and define the process. The [[request_generate_processes]] enables the RECOLA internal process compiler, which can be called only after all processes have been defined. <>= procedure :: register_processes => prc_recola_register_processes <>= - subroutine prc_recola_register_processes (writer) + subroutine prc_recola_register_processes (writer, recola_ids) class(recola_writer_t), intent(in) :: writer integer :: recola_id + integer :: i_flv integer :: n_tot integer :: unit, iostat + integer, dimension (:), intent(inout) :: recola_ids integer, dimension(:), allocatable :: pdg type(string_t), dimension(:), allocatable :: particle_names type(string_t) :: process_string integer :: i_part !!! TODO (cw-2016-08-08): Include helicities call msg_message ("Recola: registering processes for '" // char (writer%id) // "'") + i_flv = 0 n_tot = writer%n_in + writer%n_out allocate (pdg (n_tot)) allocate (particle_names (n_tot)) call open_flv_list (writer%id, unit) + call rclwrap_request_generate_processes () SCAN_FLV_LIST: do read (unit, *, iostat = iostat) pdg if (iostat < 0) then exit SCAN_FLV_LIST else if (iostat > 0) then call err_flv_list (writer%id) end if + i_flv = i_flv + 1 call rclwrap_get_new_recola_id (recola_id) + recola_ids(i_flv) = recola_id particle_names(:) = get_recola_particle_string (pdg) process_string = var_str ("") do i_part = 1, n_tot process_string = process_string // & particle_names (i_part) // var_str (" ") if (i_part == writer%n_in) then process_string = process_string // var_str ("-> ") end if end do call msg_message ("Recola: " & // "process #" // char (str (recola_id)) & // ": " // char (process_string) & // "(" // char (writer%order) // ")") - call rclwrap_define_process (recola_id, process_string, writer%order) + call rclwrap_add_process (recola_id, process_string, writer%order) + call rclwrap_define_processes () end do SCAN_FLV_LIST call close_flv_list (unit) - call rclwrap_request_generate_processes () call msg_debug (D_ME_METHODS, "RECOLA: processes for '" & // char (writer%id) // "' registered") end subroutine prc_recola_register_processes @ %def prc_recola_register_processes @ Manage the list of flavor combinations for the current process. We rely on this being created along with the \oMega\ call. <>= subroutine open_flv_list (id, unit) type(string_t), intent(in) :: id integer, intent(out) :: unit type(string_t) :: flv_file integer :: iostat flv_file = flv_file_name (id) open (file = char (flv_file), newunit = unit, & status = "old", action = "read", & iostat = iostat) if (iostat /= 0) then call msg_fatal ("Recola: attempt to open flavor-list file '" & // char (flv_file) // "' failed") end if end subroutine open_flv_list subroutine err_flv_list (id) type(string_t), intent(in) :: id type(string_t) :: flv_file flv_file = flv_file_name (id) call msg_fatal ("Recola: error while reading from flavor-list file '" & // char (flv_file) // "'") end subroutine err_flv_list subroutine close_flv_list (unit) integer, intent(in) :: unit close (unit) end subroutine close_flv_list @ %def open_flv_list @ %def err_flv_list @ %def close_flv_list @ \subsection{Driver object} A core driver is required by design. However, we are not going to load any external dynamical libraries, so this is a dummy. <>= type, extends (user_defined_driver_t) :: recola_driver_t contains <> end type recola_driver_t @ %def recola_driver_t @ <>= procedure :: allocate_driver => recola_def_allocate_driver <>= subroutine recola_def_allocate_driver (object, driver, basename) class(recola_def_t), intent(in) :: object class(prc_core_driver_t), intent(out), allocatable :: driver type(string_t), intent(in) :: basename call msg_debug2 (D_ME_METHODS, "recola_def_allocate_driver") allocate (recola_driver_t :: driver) end subroutine recola_def_allocate_driver @ %def recola_def_allocate_driver @ <>= procedure, nopass :: type_name => recola_driver_type_name <>= function recola_driver_type_name () result (type) type(string_t) :: type type = "Recola" end function recola_driver_type_name @ %def recola_driver_type_name @ \subsection{Process object} We create [[prc_recola_t]] as an extension of the [[prc_user_defined_base_t]], which in turn inherits from [[prc_core_t]]. This way, we can use a lot of the existing interfaces in the actual code. However, we have to stick to the rules and implement the deferred type-bound procedures of [[prc_core_t]]. <>= public :: prc_recola_t <>= type, extends (prc_user_defined_base_t) :: prc_recola_t - integer :: recola_id = 0 - integer, dimension(:,:), allocatable :: color_state + integer, dimension (:), allocatable :: recola_ids + integer, dimension (:,:), allocatable :: color_state integer :: n_f = 0 logical :: helicity_and_color_arrays_are_replaced = .false. contains <> end type prc_recola_t @ %def prc_recola_t @ <>= procedure :: write_name => prc_recola_write_name <>= subroutine prc_recola_write_name (object, unit) class(prc_recola_t), intent(in) :: object integer, intent(in), optional :: unit integer :: u u = given_output_unit (unit) write (u,"(1x,A)") "Core: Recola" end subroutine prc_recola_write_name @ %def prc_recola_write_name @ <>= procedure :: has_matrix_element => prc_recola_has_matrix_element <>= function prc_recola_has_matrix_element (object) result (flag) logical :: flag class(prc_recola_t), intent(in) :: object flag = .true. end function prc_recola_has_matrix_element @ %def prc_recola_has_matrix_element @ Not implemented yet. <>= procedure :: write => prc_recola_write <>= subroutine prc_recola_write (object, unit) class(prc_recola_t), intent(in) :: object integer, intent(in), optional :: unit end subroutine prc_recola_write @ %def prc_recola_write @ \subsection{Accompanying state object} This must be implemented, but is unused. <>= type, extends (user_defined_state_t) :: recola_state_t contains <> end type recola_state_t @ %def recola_state_t @ <>= procedure :: write => recola_state_write <>= subroutine recola_state_write (object, unit) class(recola_state_t), intent(in) :: object integer, intent(in), optional :: unit end subroutine recola_state_write @ %def recola_state_write @ <>= procedure :: allocate_workspace => prc_recola_allocate_workspace <>= subroutine prc_recola_allocate_workspace (object, core_state) class(prc_recola_t), intent(in) :: object class(prc_core_state_t), intent(inout), allocatable :: core_state allocate (recola_state_t :: core_state) end subroutine prc_recola_allocate_workspace @ %def prc_recola_allocate_workspace @ \subsection{Recola process data} This information is stored in the associated [[def]] object. To obtain it, we need a type cast. <>= procedure :: get_alpha_power => prc_recola_get_alpha_power procedure :: get_alphas_power => prc_recola_get_alphas_power <>= function prc_recola_get_alpha_power (object) result (p) class(prc_recola_t), intent(in) :: object integer :: p p = 0 if (associated (object%def)) then select type (def => object%def) type is (recola_def_t) p = def%alpha_power end select end if end function prc_recola_get_alpha_power function prc_recola_get_alphas_power (object) result (p) class(prc_recola_t), intent(in) :: object integer :: p p = 0 if (associated (object%def)) then select type (def => object%def) type is (recola_def_t) p = def%alphas_power end select end if end function prc_recola_get_alphas_power @ %def prc_recola_get_alpha_power @ %def prc_recola_get_alphas_power @ <>= procedure :: compute_alpha_s => prc_recola_compute_alpha_s <>= subroutine prc_recola_compute_alpha_s (object, core_state, ren_scale) class(prc_recola_t), intent(in) :: object class(user_defined_state_t), intent(inout) :: core_state real(default), intent(in) :: ren_scale core_state%alpha_qcd = object%qcd%alpha%get (ren_scale) end subroutine prc_recola_compute_alpha_s @ %def prc_recola_compute_alpha_s @ <>= procedure :: includes_polarization => prc_recola_includes_polarization <>= function prc_recola_includes_polarization (object) result (polarized) logical :: polarized class(prc_recola_t), intent(in) :: object polarized = .false. end function prc_recola_includes_polarization @ %def prc_recola_includes_polarization @ \subsection{Prepare for process evaluation} This has become obsolete and is empty. <>= procedure :: create_and_load_extra_libraries => & prc_recola_create_and_load_extra_libraries <>= subroutine prc_recola_create_and_load_extra_libraries & (core, flv_states, var_list, os_data, libname, model, i_core, is_nlo) class(prc_recola_t), intent(inout) :: core integer, intent(in), dimension(:,:), allocatable :: flv_states type(var_list_t), intent(in) :: var_list type(os_data_t), intent(in) :: os_data type(string_t), intent(in) :: libname type(model_data_t), intent(in), target :: model integer, intent(in) :: i_core logical, intent(in) :: is_nlo call msg_debug (D_ME_METHODS, "prc_recola_create_and_load_extra_libraries (no-op)") end subroutine prc_recola_create_and_load_extra_libraries @ %def prc_recola_create_and_load_extra_libraries @ Set all Recola parameters to their correct values. We use the model object for masses and such. Note that the QCD object provides the [[n_f]] parameter which affects $\alpha_s$ evaluation. Note that this is executed before the [[init]] method below, which defines and prepares the Recola process objects. This is in line with the Recola workflow, however. <>= procedure :: set_parameters => prc_recola_set_parameters <>= subroutine prc_recola_set_parameters (object, qcd, model) class(prc_recola_t), intent(inout) :: object type(qcd_t), intent(in) :: qcd class(model_data_t), intent(in), target, optional :: model call msg_debug (D_ME_METHODS, "RECOLA: set_parameters") object%qcd = qcd call rclwrap_set_dynamic_settings () - call rclwrap_set_pole_mass & (11, dble(model%get_real (var_str ('me'))), 0._double) call rclwrap_set_pole_mass & (13, dble(model%get_real (var_str ('mmu'))), 0._double) call rclwrap_set_pole_mass & (15, dble(model%get_real (var_str ('mtau'))), 0._double) call rclwrap_set_pole_mass (1, 0._double, 0._double) call rclwrap_set_pole_mass (2, 0._double, 0._double) call rclwrap_set_pole_mass (3, dble(model%get_real (var_str ('ms'))), 0._double) call rclwrap_set_pole_mass (4, dble(model%get_real (var_str ('mc'))), 0._double) call rclwrap_set_pole_mass (5, dble(model%get_real (var_str ('mb'))), 0._double) call rclwrap_set_pole_mass (6, dble(model%get_real (var_str ('mtop'))), & dble(model%get_real (var_str ('wtop')))) call rclwrap_set_pole_mass (23, dble(model%get_real (var_str ('mZ'))), & dble(model%get_real (var_str ('wZ')))) call rclwrap_set_pole_mass (24, dble(model%get_real (var_str ('mW'))), & dble(model%get_real (var_str ('wW')))) call rclwrap_set_pole_mass (25, dble(model%get_real (var_str ('mH'))), & dble(model%get_real (var_str ('wH')))) call rclwrap_use_gfermi_scheme (dble(model%get_real (var_str ('GF')))) call rclwrap_set_light_fermions (0._double) call rclwrap_set_delta_ir (0._double, dble(pi**2 / 6)) end subroutine prc_recola_set_parameters @ %def prc_recola_set_parameters @ <>= procedure :: set_mu_ir => prc_recola_set_mu_ir <>= subroutine prc_recola_set_mu_ir (object, mu) class(prc_recola_t), intent(inout) :: object real(default), intent(in) :: mu call rclwrap_set_mu_ir (dble(mu)) end subroutine prc_recola_set_mu_ir @ %def prc_recola_set_mu_ir @ Extend the base-type initialization method by Recola-specific initialization. We take the process definitions from the [[def]] object, which has been filled before. The [[writer]] component of the process-definition object can now complete its task and prepare the Recola processes. Sadly, we have to completely reset Recola first, since Recola does not allow to modify \emph{anything} after process definition. Also, we cannot really make use of Recola's multi-process capability without violating the Whizard convention that the parameter settings at process integration time apply, not at process definition time. Each new process (i.e., process-integration) object will thus trigger a complete new Recola instance. TODO: There is support just for one process, Recola ID = 1. This should be expanded. <>= procedure :: init => prc_recola_init <>= subroutine prc_recola_init (object, def, lib, id, i_component) class(prc_recola_t), intent(inout) :: object class(prc_core_def_t), intent(in), target :: def type(process_library_t), intent(in), target :: lib type(string_t), intent(in) :: id integer, intent(in) :: i_component + integer :: n_flv call msg_debug (D_ME_METHODS, "RECOLA: init process object") call object%base_init (def, lib, id, i_component) - call rclwrap_reset_recola () + n_flv = object%get_n_flvs (1) + allocate (object%recola_ids(n_flv)) select type (writer => object%def%writer) type is (recola_writer_t) - call writer%register_processes () + call writer%register_processes (object%recola_ids) end select call rclwrap_generate_processes () - object%recola_id = rclwrap_get_n_processes () call object%replace_helicity_and_color_arrays () end subroutine prc_recola_init @ %def prc_recola_init @ Recola can compute dressed amplitudes, but it needs helicity and color to be in its own format to do so. <>= procedure :: replace_helicity_and_color_arrays => & prc_recola_replace_helicity_and_color_arrays <>= subroutine prc_recola_replace_helicity_and_color_arrays (object) + ! TODO: Adjust routine for multiple recola ids class(prc_recola_t), intent(inout) :: object integer, dimension(:,:), allocatable :: col_recola integer :: i call msg_debug (D_ME_METHODS, "RECOLA: replace_helicity_and_color_arrays") deallocate (object%data%hel_state) call rclwrap_get_helicity_configurations & - (object%recola_id, object%data%hel_state) - call rclwrap_get_color_configurations (object%recola_id, col_recola) + (object%recola_ids(1), object%data%hel_state) + call rclwrap_get_color_configurations (object%recola_ids(1), col_recola) allocate (object%color_state (object%data%n_in + object%data%n_out, & size (col_recola, dim = 2))) do i = 1, size (col_recola, dim = 2) object%color_state (:, i) = col_recola (:, i) end do object%data%n_hel = size (object%data%hel_state, dim = 1) end subroutine prc_recola_replace_helicity_and_color_arrays @ %def prc_recola_replace_helicity_and_color_arrays @ \subsection{Compute matrix element} Computes the amplitude as a function of the phase space point, the flavor, helicity and color index. It is currently only used in the form by [[prc_omega_t]], all the other ones use different interfaces. H With RECOLA, we might be able to use this, too. The current implementation can fail due to missing helicity initialization. <>= procedure :: compute_amplitude => prc_recola_compute_amplitude <>= function prc_recola_compute_amplitude & (object, j, p, f, h, c, fac_scale, ren_scale, alpha_qcd_forced, & core_state) result (amp) complex(default) :: amp class(prc_recola_t), intent(in) :: object integer, intent(in) :: j type(vector4_t), intent(in), dimension(:) :: p integer, intent(in) :: f, h, c real(default), intent(in) :: fac_scale, ren_scale real(default), intent(in), allocatable :: alpha_qcd_forced class(prc_core_state_t), intent(inout), allocatable, optional :: & core_state real(double), dimension(0:3, object%data%n_in + object%data%n_out) :: & p_recola integer :: i logical :: new_event complex(double) :: amp_dble call msg_debug2 (D_ME_METHODS, "prc_recola_compute_amplitude") if (present (core_state)) then if (allocated (core_state)) then select type (core_state) type is (recola_state_t) new_event = core_state%new_kinematics core_state%new_kinematics = .false. end select end if end if if (new_event) then do i = 1, object%data%n_in + object%data%n_out p_recola(:, i) = dble(p(i)%p) end do - call rclwrap_compute_process (object%recola_id, p_recola, 'LO') + call rclwrap_compute_process (object%recola_ids(f), p_recola, 'LO') end if - call rclwrap_get_amplitude (object%recola_id, 0, 'LO', & + call rclwrap_get_amplitude (object%recola_ids(f), 0, 'LO', & object%color_state (:, c), object%data%hel_state (h, :), amp_dble) amp = amp_dble end function prc_recola_compute_amplitude @ %def prc_recola_compute_amplitude @ <>= procedure :: compute_sqme => prc_recola_compute_sqme <>= subroutine prc_recola_compute_sqme (object, i_flv, i_hel, p, & ren_scale, sqme, bad_point) class(prc_recola_t), intent(in) :: object integer, intent(in) :: i_flv, i_hel type(vector4_t), dimension(:), intent(in) :: p real(default), intent(in) :: ren_scale real(default), intent(out) :: sqme logical, intent(out) :: bad_point real(double) :: sqme_dble real(double), dimension(0:3, object%data%n_in + object%data%n_out) :: & p_recola real(default) :: alpha_s integer :: i integer :: alphas_power ! TODO sbrass: Helicity for RECOLA call msg_debug2 (D_ME_METHODS, "prc_recola_compute_sqme") do i = 1, object%data%n_in + object%data%n_out p_recola(:, i) = dble(p(i)%p) end do call rclwrap_set_mu_ir (dble (ren_scale)) alpha_s = object%qcd%alpha%get (ren_scale) call msg_debug2 (D_ME_METHODS, "alpha_s", alpha_s) call msg_debug2 (D_ME_METHODS, "ren_scale", ren_scale) call rclwrap_set_alpha_s (dble (alpha_s), dble (ren_scale), object%qcd%n_f) - call rclwrap_compute_process (i_flv, p_recola, 'LO') + call rclwrap_compute_process (object%recola_ids(i_flv), p_recola, 'LO') call rclwrap_get_squared_amplitude & - (i_flv, object%get_alphas_power (), 'LO', sqme_dble) + (object%recola_ids(i_flv), object%get_alphas_power (), 'LO', sqme_dble) sqme = real(sqme_dble, kind=default) bad_point = .false. end subroutine prc_recola_compute_sqme @ %def prc_recola_compute_sqme @ <>= procedure :: compute_sqme_virt => prc_recola_compute_sqme_virt <>= subroutine prc_recola_compute_sqme_virt (object, i_flv, i_hel, & p, ren_scale, sqme, bad_point) class(prc_recola_t), intent(in) :: object integer, intent(in) :: i_flv, i_hel type(vector4_t), dimension(:), intent(in) :: p real(default), intent(in) :: ren_scale real(default), dimension(4), intent(out) :: sqme real(default) :: amp logical, intent(out) :: bad_point real(double), dimension(0:3, object%data%n_in + object%data%n_out) :: & p_recola real(double) :: sqme_dble real(default) :: alpha_s integer :: i ! TODO sbrass Helicity for RECOLA call msg_debug2 (D_ME_METHODS, "prc_recola_compute_sqme_virt") sqme = zero do i = 1, object%data%n_in + object%data%n_out p_recola(:, i) = dble(p(i)%p) end do call rclwrap_set_mu_ir (dble (ren_scale)) alpha_s = object%qcd%alpha%get (ren_scale) call rclwrap_set_alpha_s (dble (alpha_s), dble (ren_scale), object%qcd%n_f) - call rclwrap_compute_process (object%recola_id, p_recola, 'NLO') + call rclwrap_compute_process (object%recola_ids(i_flv), p_recola, 'NLO') call rclwrap_get_squared_amplitude & - (object%recola_id, object%get_alphas_power () + 1, 'NLO', sqme_dble) + (object%recola_ids(i_flv), object%get_alphas_power () + 1, 'NLO', sqme_dble) sqme(3) = sqme_dble call rclwrap_get_squared_amplitude & - (object%recola_id, object%get_alphas_power (), 'LO', sqme_dble) + (object%recola_ids(i_flv), object%get_alphas_power (), 'LO', sqme_dble) sqme(4) = sqme_dble bad_point = .false. end subroutine prc_recola_compute_sqme_virt @ %def prc_recola_compute_sqme_virt @ \subsection{Unit tests} <<[[prc_recola_ut.f90]]>>= <> module prc_recola_ut use unit_tests use prc_recola_uti <> <> contains <> end module prc_recola_ut @ %def prc_recola_ut @ <<[[prc_recola_uti.f90]]>>= <> module prc_recola_uti use recola_wrapper !NODEP! use, intrinsic :: iso_c_binding !NODEP! use kinds <> use constants use format_utils, only: write_separator use numeric_utils, only: assert_equal use os_interface use particle_specifiers, only: new_prt_spec use prc_core_def use process_constants use process_libraries use prc_core use prc_omega <> <> contains <> <> end module prc_recola_uti @ %def prc_recola_uti @ <>= public :: prc_recola_test <>= subroutine prc_recola_test (u, results) integer, intent(in) :: u type(test_results_t), intent(inout) :: results <> end subroutine prc_recola_test @ %def prc_recola_test @ \subsubsection{Testing a fixed flavor matrix element computation} <>= function get_omega_parameter_array () result (par) real(default), dimension(25) :: par par = zero par(1) = 1.16637d-5 ! gf par(2) = 91.153480619182744_default ! mZ par(3) = 80.357973609877547_default ! mW par(4) = 125._default ! mH par(5) = rclwrap_get_alpha_s () ! alpha_s par(12) = 173.2_default ! mt par(14) = 2.4942663787728243_default ! wZ par(15) = 2.0842989982782196_default ! wW par(22) = one / sqrt (sqrt (two) * par(1)) ! par%v - Higgs expectation value par(23) = par(3) / par(2) ! par%cw par(24) = sqrt (one - par(23)**2) ! par%sw par(25) = two * par(24) * par(3) / par(22) end function get_omega_parameter_array @ %def get_omega_parameter_array @ <>= call test (prc_recola_1, "prc_recola_1", & "Registering a RECOLA process and computing the amplitude", & u, results) <>= public :: prc_recola_1 <>= subroutine prc_recola_1 (u) integer, intent(in) :: u real(double) :: p(0:3,1:4) real(double) :: sqrts = 500._double real(double) :: m_e = 0._double real(double) :: m_mu = 0._double real(double) :: p_x_out, p_y_out, p_z_out, p_z_in integer :: h_e_p, h_e_m, h_mu_p, h_mu_m, counter real(double) :: sqme integer :: i integer, dimension(:), allocatable :: col_recola, hel_recola complex(double) :: amp_recola complex(default) :: amp_recola_default real(default), parameter :: ee = 0.3 !!! Electromagnetic coupling type(process_library_t) :: lib class(prc_core_def_t), allocatable :: def type(process_def_entry_t), pointer :: entry type(string_t), dimension(:), allocatable :: prt_in, prt_out type(os_data_t) :: os_data type(process_constants_t) :: data class(prc_core_driver_t), allocatable :: driver complex(default) :: amp integer, dimension(:,:), allocatable :: helicities write (u, "(A)") "* Test output: prc_recola_1" write (u, "(A)") "* Purpose: Test interface to RECOLA and compare matrix elements with O'Mega" write (u, "(A)") p_z_in = sqrt ((sqrts / 2)**2 - m_e**2) p_z_out = 0._double p_y_out = sqrts / 10._default p_x_out = sqrt ((sqrts / 2)**2 - p_y_out**2 - p_z_out**2 - m_mu**2) p(:,1) = [sqrts / 2, 0._double, 0._double, p_z_in] p(:,2) = [sqrts / 2, 0._double, 0._double, -p_z_in] p(:,3) = [sqrts / 2, p_x_out, p_y_out, p_z_out] p(:,4) = [sqrts / 2, -p_x_out, -p_y_out, -p_z_out] write (u, "(A)") "Use phase-space point: " do i = 1, 4 write (u, "(4(F12.3,1x))") p(:,1) end do write (u, "(A)") call write_separator (u) write (u, "(A)") write (u, "(A)") "* RECOLA: Evaluate process" counter = 1 + call rclwrap_request_generate_processes () write (u, "(A)") "* RECOLA: Define process e+ e- -> mu+ mu- at leading order" - call rclwrap_define_process (counter, var_str ('e+ e- -> mu+ mu-'), var_str ('LO')) + call rclwrap_add_process (counter, var_str ('e+ e- -> mu+ mu-'), var_str ('LO')) + call rclwrap_define_processes () write (u, "(A)") "* RECOLA: generate process" - call rclwrap_request_generate_processes () call rclwrap_generate_processes () call rclwrap_compute_process (1, p, 'LO') call rclwrap_get_helicity_configurations (1, helicities) allocate (hel_recola (4), col_recola (4)) col_recola = [0,0,0,0] write (u, "(A)") "* Setting up Omega to compute the same amplitude" call lib%init (var_str ("omega1")) allocate (prt_in (2), prt_out (2)) prt_in = [var_str ("e+"), var_str ("e-")] prt_out = [var_str ("mu+"), var_str ("mu-")] allocate (omega_def_t :: def) select type (def) type is (omega_def_t) call def%init (var_str ("SM"), prt_in, prt_out, & ufo = .false., ovm = .false., cms_scheme = .true.) end select allocate (entry) call entry%init (var_str ("omega1_a"), model_name = var_str ("SM"), & n_in = 2, n_components = 1) call entry%import_component (1, n_out = 2, & prt_in = new_prt_spec (prt_in), & prt_out = new_prt_spec (prt_out), & method = var_str ("omega"), & variant = def) call lib%append (entry) call os_data_init (os_data) call lib%configure (os_data) call lib%write_makefile (os_data, force = .true., verbose = .false.) call lib%clean (os_data, distclean = .false.) call lib%write_driver (force = .true.) call lib%load (os_data) call lib%connect_process (var_str ("omega1_a"), 1, data, driver) select type (driver) type is (omega_driver_t) call driver%init (get_omega_parameter_array (), 3) call driver%new_event (real(p, kind = default)) do i = 1, 6 call rclwrap_get_amplitude (1, 0, 'LO', col_recola, helicities (:,i), amp_recola) end do do i = 1, 16 call rclwrap_get_amplitude (1, 0, 'LO', col_recola, data%hel_state (:,i), amp_recola) amp_recola = amp_recola * cmplx (0, -1, double) amp_recola_default = amp_recola call driver%get_amplitude (1, i, 1, amp) write(u,"(A,4(I2),A)") "Helicity: [",data%hel_state (:,i),"]" call assert_equal (u, amp, amp_recola_default, rel_smallness = 1.E-7_default) end do end select call rclwrap_reset_recola () write (u, "(A)") write (u, "(A)") "* End of test output: prc_recola_1" end subroutine prc_recola_1 @ %def prc_recola_1 @ \subsubsection{Testing a fixed flavor matrix element computation for 2->3} <>= call test (prc_recola_2, "prc_recola_2", & "Registering a RECOLA process and computing the amplitude for 2->3 process", & u, results) <>= public :: prc_recola_2 <>= subroutine prc_recola_2 (u) integer, intent(in) :: u real(double) :: p(0:3,1:5) real(double) :: sqrts = 700._double real(double) :: m_e = 0._double real(double) :: m_mu = 0._double real(double) :: p_x_out, p_y_out, p_z_out, p_z_in real(double) :: sqme integer :: i integer, dimension(:), allocatable :: col_recola, hel_recola integer, dimension(:,:), allocatable :: helicities complex(double) :: amp_recola complex(default) :: amp_recola_default real(default), parameter :: ee = 0.3 !!! Electromagnetic coupling type(process_library_t) :: lib class(prc_core_def_t), allocatable :: def type(process_def_entry_t), pointer :: entry type(string_t), dimension(:), allocatable :: prt_in, prt_out type(os_data_t) :: os_data type(process_constants_t) :: data class(prc_core_driver_t), allocatable :: driver complex(default) :: amp integer :: n_allowed write (u, "(A)") "* Test output: prc_recola_2" write (u, "(A)") "* Purpose: Test interface to RECOLA and compare matrix elements with O'Mega for 2->3 process" write (u, "(A)") p_z_in = sqrt ((sqrts / 2)**2 - m_e**2) p(:,1) = [sqrts / 2, 0._double, 0._double, p_z_in] p(:,2) = [sqrts / 2, 0._double, 0._double, -p_z_in] p(:,3) = [243.49323116_double, -141.69619338_double, -108.30640321_double, 165.77353656_double] p(:,4) = [337.53250628_double, 143.95931207_double, 110.19717026_double, -284.71124482_double] p(:,5) = [118.97426257_double, -2.2631186860_double, -1.8907670459_double, 118.93770827_double] write (u, "(A)") "Use phase-space point: " do i = 1, 5 write (u, "(4(F12.3,1x))") p(:,1) end do write (u, "(A)") call write_separator (u) write (u, "(A)") write (u, "(A)") "* RECOLA: Evaluate process" + call rclwrap_request_generate_processes () write (u, "(A)") "* RECOLA: Define process e+ e- -> mu+ mu- A at leading order" - call rclwrap_define_process (2, var_str ('e+ e- -> mu+ mu- A'), var_str ('LO')) + call rclwrap_add_process (2, var_str ('e+ e- -> mu+ mu- A'), var_str ('LO')) + call rclwrap_define_processes () write (u, "(A)") "* RECOLA: generate process" - call rclwrap_request_generate_processes () call rclwrap_generate_processes () call rclwrap_compute_process (2, p, 'LO') call rclwrap_get_helicity_configurations (2, helicities) allocate (hel_recola (5), col_recola (5)) col_recola = [0,0,0,0,0] write (u, "(A)") "* Setting up Omega to compute the same amplitude" call lib%init (var_str ("omega2")) allocate (prt_in (2), prt_out (3)) prt_in = [var_str ("e+"), var_str ("e-")] prt_out = [var_str ("mu+"), var_str ("mu-"), var_str("A")] allocate (omega_def_t :: def) select type (def) type is (omega_def_t) call def%init (var_str ("SM"), prt_in, prt_out, & ufo = .false., ovm = .false.) end select allocate (entry) call entry%init (var_str ("omega2_a"), model_name = var_str ("SM"), & n_in = 2, n_components = 1) call entry%import_component (1, n_out = 3, & prt_in = new_prt_spec (prt_in), & prt_out = new_prt_spec (prt_out), & method = var_str ("omega"), & variant = def) call lib%append (entry) call os_data_init (os_data) call lib%configure (os_data) call lib%write_makefile (os_data, force = .true., verbose = .false.) call lib%clean (os_data, distclean = .false.) call lib%write_driver (force = .true.) call lib%load (os_data) call lib%connect_process (var_str ("omega2_a"), 1, data, driver) select type (driver) type is (omega_driver_t) call driver%init (get_omega_parameter_array (), 3) call driver%new_event (real(p, kind = default)) do i = 1, 32 call rclwrap_get_amplitude & (2, 0, 'LO', col_recola, data%hel_state (:,i), amp_recola) if (data%hel_state(3,i) * data%hel_state(4,i) * & data%hel_state(5,i) == -1) then amp_recola = amp_recola * cmplx (0, -1, double) else amp_recola = amp_recola * cmplx (0, 1, double) end if amp_recola_default = amp_recola call driver%get_amplitude (1, i, 1, amp) write(u,"(A,5(I2),A)") "Helicity: [", data%hel_state (:,i),"]" write(u,"(A,2(F12.7,1x),A,2(F12.7,1x))") "RECOLA:", & amp_recola,", O'MEGA:", amp call assert_equal & (u, amp, amp_recola_default, rel_smallness = 1.E-6_default) end do end select call rclwrap_reset_recola () write (u, "(A)") write (u, "(A)") "* End of test output: prc_recola_2" end subroutine prc_recola_2 @ %def prc_recola_2 @ Index: trunk/share/tests/functional_tests/ref-output/recola_2.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_2.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_2.ref (revision 8179) @@ -1,62 +1,62 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 | Process library 'recola_2_lib': recorded process 'recola_2_p1' seed = 1 sqrts = 5.00000E+02 | Integrate: current process library needs compilation | Process library 'recola_2_lib': compiling ... | Process library 'recola_2_lib': writing makefile | Process library 'recola_2_lib': removing old files | Process library 'recola_2_lib': writing driver | Process library 'recola_2_lib': creating source code | Process library 'recola_2_lib': compiling sources | Process library 'recola_2_lib': linking | Process library 'recola_2_lib': loading | Process library 'recola_2_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_2_p1: -| Recola: registering processes for 'recola_2_p1' +| Recola: registering processes for 'recola_2_p1_BORN' | Recola: process #1: e+ e- -> t t~ (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | e+ (mass = 5.1099700E-04 GeV) | e- (mass = 5.1099700E-04 GeV) | sqrts = 5.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_2_p1.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_2_p1' | Library name = 'recola_2_lib' | Process index = 1 | Process components: | 1: 'recola_2_p1_i1': e+, e- => t, tbar [recola] | ------------------------------------------------------------------------ | Phase space: 1 channels, 2 dimensions | Phase space: found 1 channel, collected in 1 grove. | Phase space: Using 1 equivalence between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_2_p1' | Integrate: iterations = 1:100 | Integrator: 1 chains, 1 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 5.068E+02 2.49E+01 4.92 0.49 46.3 |-----------------------------------------------------------------------------| 1 100 5.068E+02 2.49E+01 4.92 0.49 46.3 |=============================================================================| | There were no errors and 1 warning(s). | WHIZARD run finished. |=============================================================================| Index: trunk/share/tests/functional_tests/ref-output/recola_4.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_4.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_4.ref (revision 8179) @@ -1,63 +1,63 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 alphas_power = 2 | Process library 'recola_4_lib': recorded process 'recola_4_p1' seed = 1 sqrts = 5.00000E+02 | Integrate: current process library needs compilation | Process library 'recola_4_lib': compiling ... | Process library 'recola_4_lib': writing makefile | Process library 'recola_4_lib': removing old files | Process library 'recola_4_lib': writing driver | Process library 'recola_4_lib': creating source code | Process library 'recola_4_lib': compiling sources | Process library 'recola_4_lib': linking | Process library 'recola_4_lib': loading | Process library 'recola_4_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_4_p1: -| Recola: registering processes for 'recola_4_p1' +| Recola: registering processes for 'recola_4_p1_BORN' | Recola: process #1: u u~ -> t t~ (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | u (mass = 0.0000000E+00 GeV) | ubar (mass = 0.0000000E+00 GeV) | sqrts = 5.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_4_p1.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_4_p1' | Library name = 'recola_4_lib' | Process index = 1 | Process components: | 1: 'recola_4_p1_i1': u, ubar => t, tbar [recola] | ------------------------------------------------------------------------ | Phase space: 1 channels, 2 dimensions | Phase space: found 1 channel, collected in 1 grove. | Phase space: Using 1 equivalence between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_4_p1' | Integrate: iterations = 1:100 | Integrator: 1 chains, 1 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 1.791E+04 1.76E+02 0.98 0.10 83.0 |-----------------------------------------------------------------------------| 1 100 1.791E+04 1.76E+02 0.98 0.10 83.0 |=============================================================================| | There were no errors and 1 warning(s). | WHIZARD run finished. |=============================================================================| Index: trunk/share/tests/functional_tests/ref-output/recola_6.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_6.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_6.ref (revision 8179) @@ -1,63 +1,64 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 | Process library 'recola_6_lib': recorded process 'recola_6_p1' seed = 1 sqrts = 5.00000E+02 | Integrate: current process library needs compilation | Process library 'recola_6_lib': compiling ... | Process library 'recola_6_lib': writing makefile | Process library 'recola_6_lib': removing old files | Process library 'recola_6_lib': writing driver | Process library 'recola_6_lib': creating source code | Process library 'recola_6_lib': compiling sources | Process library 'recola_6_lib': linking | Process library 'recola_6_lib': loading | Process library 'recola_6_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_6_p1: -| Recola: registering processes for 'recola_6_p1' +| Recola: registering processes for 'recola_6_p1_BORN' | Recola: process #1: e+ e- -> u u~ (LO) | Recola: process #2: e+ e- -> d d~ (LO) | Recola: preparing processes for integration +| Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | e+ (mass = 5.1099700E-04 GeV) | e- (mass = 5.1099700E-04 GeV) | sqrts = 5.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_6_p1.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_6_p1' | Library name = 'recola_6_lib' | Process index = 1 | Process components: | 1: 'recola_6_p1_i1': e+, e- => u:d, ubar:dbar [recola] | ------------------------------------------------------------------------ | Phase space: 2 channels, 2 dimensions | Phase space: found 2 channels, collected in 2 groves. | Phase space: Using 2 equivalences between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_6_p1' | Integrate: iterations = 1:100 | Integrator: 2 chains, 2 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 9.584E+02 7.55E+01 7.87 0.79 34.7 |-----------------------------------------------------------------------------| 1 100 9.584E+02 7.55E+01 7.87 0.79 34.7 |=============================================================================| | There were no errors and 1 warning(s). | WHIZARD run finished. |=============================================================================| Index: trunk/share/tests/functional_tests/ref-output/recola_8.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_8.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_8.ref (revision 8179) @@ -1,60 +1,60 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 alphas_power = 0 | Process library 'recola_8_lib': recorded process 'recola_8_p1' seed = 1 | Integrate: current process library needs compilation | Process library 'recola_8_lib': compiling ... | Process library 'recola_8_lib': writing makefile | Process library 'recola_8_lib': removing old files | Process library 'recola_8_lib': writing driver | Process library 'recola_8_lib': creating source code | Process library 'recola_8_lib': compiling sources | Process library 'recola_8_lib': linking | Process library 'recola_8_lib': loading | Process library 'recola_8_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_8_p1: -| Recola: registering processes for 'recola_8_p1' +| Recola: registering processes for 'recola_8_p1_BORN' | Recola: process #1: t -> W+ b (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (decay): | t (mass = 1.7310000E+02 GeV) | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_8_p1.i1.phs' | ------------------------------------------------------------------------ | Process [decay]: 'recola_8_p1' | Library name = 'recola_8_lib' | Process index = 1 | Process components: | 1: 'recola_8_p1_i1': t => W+, b [recola] | ------------------------------------------------------------------------ | Phase space: 1 channels, 2 dimensions | Phase space: found 1 channel, collected in 1 grove. | Phase space: Using 1 equivalence between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_8_p1' | Integrate: iterations = 1:100 | Integrator: 1 chains, 1 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[GeV] Error[GeV] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 1.495E+00 0.00E+00 0.00 0.00 100.0 |-----------------------------------------------------------------------------| 1 100 1.495E+00 0.00E+00 0.00 0.00 100.0 |=============================================================================| | There were no errors and 1 warning(s). | WHIZARD run finished. |=============================================================================| Index: trunk/share/tests/functional_tests/ref-output/recola_1.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_1.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_1.ref (revision 8179) @@ -1,229 +1,229 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 | Process library 'recola_1_lib': recorded process 'recola_1_p1' seed = 1 sqrts = 1.00000E+02 | Integrate: current process library needs compilation | Process library 'recola_1_lib': compiling ... | Process library 'recola_1_lib': writing makefile | Process library 'recola_1_lib': removing old files | Process library 'recola_1_lib': writing driver | Process library 'recola_1_lib': creating source code | Process library 'recola_1_lib': compiling sources | Process library 'recola_1_lib': linking | Process library 'recola_1_lib': loading | Process library 'recola_1_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_1_p1: -| Recola: registering processes for 'recola_1_p1' +| Recola: registering processes for 'recola_1_p1_BORN' | Recola: process #1: e+ e- -> mu+ mu- (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | e+ (mass = 5.1099700E-04 GeV) | e- (mass = 5.1099700E-04 GeV) | sqrts = 1.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_1_p1.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_1_p1' | Library name = 'recola_1_lib' | Process index = 1 | Process components: | 1: 'recola_1_p1_i1': e+, e- => mu+, mu- [recola] | ------------------------------------------------------------------------ | Phase space: 2 channels, 2 dimensions | Phase space: found 2 channels, collected in 2 groves. | Phase space: Using 2 equivalences between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_1_p1' | Integrate: iterations = 1:100 | Integrator: 2 chains, 2 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 5.750E+04 3.95E+03 6.87 0.69 41.5 |-----------------------------------------------------------------------------| 1 100 5.750E+04 3.95E+03 6.87 0.69 41.5 |=============================================================================| ?sample_pacify = true ?debug_decay = false ?debug_process = false ?debug_verbose = false ?write_raw = false n_events = 2 | Starting simulation for process 'recola_1_p1' | Simulate: using integration grids from file 'recola_1_p1.m1.vg' | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 2 | Simulation: requested number of events = 2 | corr. to luminosity [fb-1] = 3.4784E-05 | Events: writing to ASCII file 'recola_1_p1.debug' | Events: generating 2 unweighted, unpolarized events ... | Events: event normalization mode '1' | ... event sample complete. | Events: actual unweighting efficiency = 66.67 % | Events: closing ASCII file 'recola_1_p1.debug' | There were no errors and 1 warning(s). | WHIZARD run finished. |=============================================================================| Contents of recola_1.debug: ======================================================================== Event #1 ------------------------------------------------------------------------ Unweighted = T Normalization = '1' Helicity handling = drop Keep correlations = F ------------------------------------------------------------------------ Squared matrix el. (ref) = 4.86324E-02 Squared matrix el. (prc) = 4.86324E-02 Event weight (ref) = 1.00000E+00 Event weight (prc) = 1.00000E+00 ------------------------------------------------------------------------ Selected MCI group = 1 Selected term = 1 Selected channel = 1 ------------------------------------------------------------------------ Passed selection = T Reweighting factor = 1.00000E+00 Analysis flag = T ======================================================================== Event transform: trivial (hard process) ------------------------------------------------------------------------ Associated process: 'recola_1_p1' TAO random-number generator: seed = 65538 calls = 3 Number of tries = 1 ------------------------------------------------------------------------ Particle set: ------------------------------------------------------------------------ Particle 1 [i] f(-11) E = 5.000000E+01 P = 0.000000E+00 0.000000E+00 5.000000E+01 T = 2.611179340E-07 Children: 3 4 Particle 2 [i] f(11) E = 5.000000E+01 P = 0.000000E+00 0.000000E+00 -5.000000E+01 T = 2.611179340E-07 Children: 3 4 Particle 3 [o] f(-13) E = 5.000000E+01 P = -2.613124E+01 4.259687E+01 -1.629005E+00 T = 1.116369517E-02 Parents: 1 2 Particle 4 [o] f(13) E = 5.000000E+01 P = 2.613124E+01 -4.259687E+01 1.629005E+00 T = 1.116369517E-02 Parents: 1 2 ======================================================================== Local variables: ------------------------------------------------------------------------ sqrts* = 1.00000E+02 sqrts_hat* => 1.00000E+02 n_in* => 2 n_out* => 2 n_tot* => 4 $process_id* => "recola_1_p1" process_num_id* => [unknown integer] sqme* => 4.86324E-02 sqme_ref* => 4.86324E-02 event_index* => 1 event_weight* => 1.00000E+00 event_weight_ref* => 1.00000E+00 event_excess* => 0.00000E+00 ------------------------------------------------------------------------ subevent: 1 prt(i:-11|-5.0000000E+01; 0.0000000E+00, 0.0000000E+00,-5.0000000E+01| 2.6111793E-07| 1) 2 prt(i:11|-5.0000000E+01; 0.0000000E+00, 0.0000000E+00, 5.0000000E+01| 2.6111793E-07| 2) 3 prt(o:-13| 5.0000000E+01;-2.6131239E+01, 4.2596872E+01,-1.6290046E+00| 1.1163695E-02| 3) 4 prt(o:13| 5.0000000E+01; 2.6131239E+01,-4.2596872E+01, 1.6290046E+00| 1.1163695E-02| 4) ======================================================================== ======================================================================== Event #2 ------------------------------------------------------------------------ Unweighted = T Normalization = '1' Helicity handling = drop Keep correlations = F ------------------------------------------------------------------------ Squared matrix el. (ref) = 1.60258E-01 Squared matrix el. (prc) = 1.60258E-01 Event weight (ref) = 1.00000E+00 Event weight (prc) = 1.00000E+00 ------------------------------------------------------------------------ Selected MCI group = 1 Selected term = 1 Selected channel = 2 ------------------------------------------------------------------------ Passed selection = T Reweighting factor = 1.00000E+00 Analysis flag = T ======================================================================== Event transform: trivial (hard process) ------------------------------------------------------------------------ Associated process: 'recola_1_p1' TAO random-number generator: seed = 65538 calls = 6 Number of tries = 1 ------------------------------------------------------------------------ Particle set: ------------------------------------------------------------------------ Particle 1 [i] f(-11) E = 5.000000E+01 P = 0.000000E+00 0.000000E+00 5.000000E+01 T = 2.611179340E-07 Children: 3 4 Particle 2 [i] f(11) E = 5.000000E+01 P = 0.000000E+00 0.000000E+00 -5.000000E+01 T = 2.611179340E-07 Children: 3 4 Particle 3 [o] f(-13) E = 5.000000E+01 P = -4.528467E+00 -2.379808E+01 4.373938E+01 T = 1.116369517E-02 Parents: 1 2 Particle 4 [o] f(13) E = 5.000000E+01 P = 4.528467E+00 2.379808E+01 -4.373938E+01 T = 1.116369517E-02 Parents: 1 2 ======================================================================== Local variables: ------------------------------------------------------------------------ sqrts* = 1.00000E+02 sqrts_hat* => 1.00000E+02 n_in* => 2 n_out* => 2 n_tot* => 4 $process_id* => "recola_1_p1" process_num_id* => [unknown integer] sqme* => 1.60258E-01 sqme_ref* => 1.60258E-01 event_index* => 2 event_weight* => 1.00000E+00 event_weight_ref* => 1.00000E+00 event_excess* => 0.00000E+00 ------------------------------------------------------------------------ subevent: 1 prt(i:-11|-5.0000000E+01; 0.0000000E+00, 0.0000000E+00,-5.0000000E+01| 2.6111793E-07| 1) 2 prt(i:11|-5.0000000E+01; 0.0000000E+00, 0.0000000E+00, 5.0000000E+01| 2.6111793E-07| 2) 3 prt(o:-13| 5.0000000E+01;-4.5284672E+00,-2.3798084E+01, 4.3739376E+01| 1.1163695E-02| 3) 4 prt(o:13| 5.0000000E+01; 4.5284672E+00, 2.3798084E+01,-4.3739376E+01| 1.1163695E-02| 4) ======================================================================== Index: trunk/share/tests/functional_tests/ref-output/recola_3.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_3.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_3.ref (revision 8179) @@ -1,62 +1,62 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 alphas_power = 2 | Process library 'recola_3_lib': recorded process 'recola_3_p1' seed = 1 sqrts = 5.00000E+01 | Integrate: current process library needs compilation | Process library 'recola_3_lib': compiling ... | Process library 'recola_3_lib': writing makefile | Process library 'recola_3_lib': removing old files | Process library 'recola_3_lib': writing driver | Process library 'recola_3_lib': creating source code | Process library 'recola_3_lib': compiling sources | Process library 'recola_3_lib': linking | Process library 'recola_3_lib': loading | Process library 'recola_3_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_3_p1: -| Recola: registering processes for 'recola_3_p1' +| Recola: registering processes for 'recola_3_p1_BORN' | Recola: process #1: g g -> g g (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | gl (mass = 0.0000000E+00 GeV) | gl (mass = 0.0000000E+00 GeV) | sqrts = 5.000000000000E+01 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_3_p1.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_3_p1' | Library name = 'recola_3_lib' | Process index = 1 | Process components: | 1: 'recola_3_p1_i1': gl, gl => gl, gl [recola] | ------------------------------------------------------------------------ | Phase space: 5 channels, 2 dimensions | Phase space: found 5 channels, collected in 2 groves. | Phase space: Using 18 equivalences between channels. | Phase space: wood | Applying user-defined cuts. | Starting integration for process 'recola_3_p1' | Integrate: iterations = 1:100 | Integrator: 2 chains, 5 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 4.416E+08 4.70E+07 10.64 1.06 28.1 |-----------------------------------------------------------------------------| 1 100 4.416E+08 4.70E+07 10.64 1.06 28.1 |=============================================================================| | WHIZARD run finished. |=============================================================================| Index: trunk/share/tests/functional_tests/ref-output/recola_5.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_5.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_5.ref (revision 8179) @@ -1,102 +1,102 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true $method = "recola" openmp_num_threads = 1 | Process library 'recola_5_lib': recorded process 'recola_5_p1' | Process library 'recola_5_lib': recorded process 'recola_5_p2' seed = 1 sqrts = 1.00000E+02 | Integrate: current process library needs compilation | Process library 'recola_5_lib': compiling ... | Process library 'recola_5_lib': writing makefile | Process library 'recola_5_lib': removing old files | Process library 'recola_5_lib': writing driver | Process library 'recola_5_lib': creating source code | Process library 'recola_5_lib': compiling sources | Process library 'recola_5_lib': linking | Process library 'recola_5_lib': loading | Process library 'recola_5_lib': ... success. | Integrate: compilation done | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 1 | Initializing integration for process recola_5_p1: -| Recola: registering processes for 'recola_5_p1' +| Recola: registering processes for 'recola_5_p1_BORN' | Recola: process #1: e+ e- -> mu+ mu- (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | e+ (mass = 5.1099700E-04 GeV) | e- (mass = 5.1099700E-04 GeV) | sqrts = 1.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_5_p1.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_5_p1' | Library name = 'recola_5_lib' | Process index = 1 | Process components: | 1: 'recola_5_p1_i1': e+, e- => mu+, mu- [recola] | ------------------------------------------------------------------------ | Phase space: 2 channels, 2 dimensions | Phase space: found 2 channels, collected in 2 groves. | Phase space: Using 2 equivalences between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_5_p1' | Integrate: iterations = 1:100 | Integrator: 2 chains, 2 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 5.750E+04 3.95E+03 6.87 0.69 41.5 |-----------------------------------------------------------------------------| 1 100 5.750E+04 3.95E+03 6.87 0.69 41.5 |=============================================================================| | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 2 | Initializing integration for process recola_5_p2: -| Recola: registering processes for 'recola_5_p2' -| Recola: process #1: e+ e- -> tau+ tau- (LO) +| Recola: registering processes for 'recola_5_p2_BORN' +| Recola: process #2: e+ e- -> tau+ tau- (LO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | e+ (mass = 5.1099700E-04 GeV) | e- (mass = 5.1099700E-04 GeV) | sqrts = 1.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_5_p2.i1.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_5_p2' | Library name = 'recola_5_lib' | Process index = 2 | Process components: | 1: 'recola_5_p2_i1': e+, e- => tau+, tau- [recola] | ------------------------------------------------------------------------ | Phase space: 2 channels, 2 dimensions | Phase space: found 2 channels, collected in 2 groves. | Phase space: Using 2 equivalences between channels. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_5_p2' | Integrate: iterations = 1:100 | Integrator: 2 chains, 2 channels, 2 dimensions | Integrator: Using VAMP channel equivalences | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 5.410E+04 4.10E+03 7.57 0.76 39.0 |-----------------------------------------------------------------------------| 1 100 5.410E+04 4.10E+03 7.57 0.76 39.0 |=============================================================================| | There were no errors and 2 warning(s). | WHIZARD run finished. |=============================================================================| Index: trunk/share/tests/functional_tests/ref-output/recola_7.ref =================================================================== --- trunk/share/tests/functional_tests/ref-output/recola_7.ref (revision 8178) +++ trunk/share/tests/functional_tests/ref-output/recola_7.ref (revision 8179) @@ -1,129 +1,129 @@ ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true SM.mW => 8.03760E+01 SM.mZ => 9.11876E+01 SM.GF => 1.16637E-05 SM.wZ => 2.49520E+00 SM.wW => 2.12400E+00 SM.mmu => 0.00000E+00 SM.me => 0.00000E+00 SM.mc => 0.00000E+00 SM.ms => 0.00000E+00 SM.wtop => 0.00000E+00 SM.mtop => 1.75000E+02 ?logging => true ?openmp_logging = false ?vis_history = false ?integration_timer = false ?pacify = true ?use_vamp_equivalences = false $loop_me_method = "recola" ?alphas_is_fixed = false ?alphas_from_mz = true ?alphas_from_lambda_qcd = false alpha_power = 2 alphas_power = 0 | Process library 'recola_7_lib': recorded process 'recola_7_p1' seed = 7777 sqrts = 5.00000E+02 | Integrate: current process library needs compilation | Process library 'recola_7_lib': compiling ... | Process library 'recola_7_lib': writing makefile | Process library 'recola_7_lib': removing old files | Process library 'recola_7_lib': writing driver | Process library 'recola_7_lib': creating source code | Process library 'recola_7_lib': compiling sources | Process library 'recola_7_lib': linking | Process library 'recola_7_lib': loading | Process library 'recola_7_lib': ... success. | Integrate: compilation done | QCD alpha: using a running strong coupling | RNG: Initializing TAO random-number generator | RNG: Setting seed for random-number generator to 7777 | Initializing integration for process recola_7_p1: -| Recola: registering processes for 'recola_7_p1' +| Recola: registering processes for 'recola_7_p1_LOOP' | Recola: process #1: e+ e- -> t t~ (NLO) | Recola: preparing processes for integration | Beam structure: [any particles] | Beam data (collision): | e+ (mass = 0.0000000E+00 GeV) | e- (mass = 0.0000000E+00 GeV) | sqrts = 5.000000000000E+02 GeV | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_7_p1.i1.phs' | Phase space: generating configuration ... | Phase space: ... success. | Phase space: writing configuration file 'recola_7_p1.i3.phs' | ------------------------------------------------------------------------ | Process [scattering]: 'recola_7_p1' | Library name = 'recola_7_lib' | Process index = 1 | Process components: | 1: 'recola_7_p1_i1': e+, e- => t, tbar [omega] | 2: 'recola_7_p1_i2': e+, e- => t, tbar, gl [omega], [real] | 3: 'recola_7_p1_i3': e+, e- => t, tbar [recola], [virtual] | 4: 'recola_7_p1_i4': e+, e- => t, tbar [inactive], [subtraction] | ------------------------------------------------------------------------ | Phase space: 1 channels, 2 dimensions | Phase space: found 1 channel, collected in 1 grove. | Phase space: no equivalences between channels used. | Phase space: wood | Phase space: 1 channels, 5 dimensions | Phase space: found 1 channel, collected in 1 grove. | Phase space: no equivalences between channels used. | Phase space: wood | Phase space: 1 channels, 2 dimensions | Phase space: found 1 channel, collected in 1 grove. | Phase space: no equivalences between channels used. | Phase space: wood Warning: No cuts have been defined. | Starting integration for process 'recola_7_p1' part 'born' | Integrate: iterations = 1:100 | Integrator: 1 chains, 1 channels, 2 dimensions | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 5.072E+02 2.64E+01 5.21 0.52 45.6 |-----------------------------------------------------------------------------| 1 100 5.072E+02 2.64E+01 5.21 0.52 45.6 |=============================================================================| | Starting integration for process 'recola_7_p1' part 'real' | Integrate: iterations = 1:100 | Integrator: 1 chains, 1 channels, 5 dimensions | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 -7.761E+01 4.10E+00 5.28 0.53 43.0 |-----------------------------------------------------------------------------| 1 100 -7.761E+01 4.10E+00 5.28 0.53 43.0 |=============================================================================| | Starting integration for process 'recola_7_p1' part 'virtual' | Integrate: iterations = 1:100 | Integrator: 1 chains, 1 channels, 2 dimensions | Integrator: 100 initial calls, 20 bins, stratified = T | Integrator: VAMP |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 100 1.388E+02 7.64E+00 5.51 0.55 45.4 |-----------------------------------------------------------------------------| 1 100 1.388E+02 7.64E+00 5.51 0.55 45.4 |=============================================================================| | Integrate: sum of all components |=============================================================================| | It Calls Integral[fb] Error[fb] Err[%] Acc Eff[%] Chi2 N[It] | |=============================================================================| 1 0 5.683E+02 2.78E+01 4.90 0.00 40.1 | NLO Correction: [O(alpha_s+1)/O(alpha_s)] | ( 12.0614 +- 1.82153 ) % |=============================================================================| | There were no errors and 1 warning(s). | WHIZARD run finished. |=============================================================================|