gibuu is hosted by Hepforge, IPPP Durham
GiBUU

TABLE OF CONTENTS


/XsectionRatios [ Modules ]

[ Top ] [ Modules ]

NAME

module XsectionRatios

NOTES

Computes and stores the ratios of the in-medium and vacuum BB -> BB + mesons cross sections (see M. Wagner et al., PRC 71, 034910 (2005).


XsectionRatios/flagScreen [ Global module-variables ]

[ Top ] [ XsectionRatios ] [ Global module-variables ]

SOURCE

  logical, save :: flagScreen = .false.

PURPOSE

If .true. -- in-medium screening is applied to the input cross section. If .false. -- no cross section modification.


XsectionRatios/ScreenMode [ Global module-variables ]

[ Top ] [ XsectionRatios ] [ Global module-variables ]

SOURCE

  integer, save :: ScreenMode = 1    ! 1 -- in-medium screening of NN total cross section
                                     !      according to Li and Machleidt, PRC 48 (1993) 1702
                                     !      and PRC 49 (1994) 566
                                     ! 2 -- in-medium screening of BB total cross section
                                     !      according to P. Daniewlewicz, NPA 673, 375 (2000);
                                     !      Acta. Phys. Pol. B 33, 45 (2002)
  ! NOTE
  ! relevant when  flagScreen = .true.

XsectionRatios/flagInMedium [ Global module-variables ]

[ Top ] [ XsectionRatios ] [ Global module-variables ]

SOURCE

  logical, save :: flagInMedium = .false.

PURPOSE

If .true. -- In-medium ratios are used to decide whether an event is accepted or not. If .false. -- The event is always accepted


XsectionRatios/InMediumMode [ Global module-variables ]

[ Top ] [ XsectionRatios ] [ Global module-variables ]

SOURCE

  integer, save :: InMediumMode = 1  ! 1 -- all events of the type BB -> BB (+ mesons) are subject to
                                     !      in-medium reduction following Eqs.(194),(195) of GiBUU review paper
                                     !      currently works in RMF mode only
                                     ! 2 -- only  N N <-> N Delta events are subject to in-medium reduction
                                     !      according to Song/Ko, arXiv:1403.7363
                                     !      works in all modes (Skyrme, RMF, cascade)
  ! NOTE
  ! relevant when flagInMedium = .true.

XsectionRatios/InMediumMode [ Global module-variables ]

[ Top ] [ XsectionRatios ] [ Global module-variables ]

SOURCE

  real, save :: alpha = 1.2

PURPOSE

Parameter which controls the density dependence of the NN <-> N Delta cross section. for the density dependence from: Song/Ko, arXiv:1403.7363 (InMediumMode=2)


XsectionRatios/shift0 [ Global module-variables ]

[ Top ] [ XsectionRatios ] [ Global module-variables ]

SOURCE

  real, save :: shift0 = 0.

PURPOSE

Mass shift m-m^* (GeV) for using in elementary particle collision mode.


XsectionRatios/getSigmaScreened [ Functions ]

[ Top ] [ XsectionRatios ] [ Functions ]

NAME

real function getSigmaScreened(pair,sigma,srts,mediumAtColl)

PURPOSE

Returns in-medium screened/reduced cross section (mbarn). See Li and Machleidt, PRC 48 (1993) 1702 and PRC 49 (1994) 566, P. Daniewlewicz, NPA 673, 375 (2000); Acta. Phys. Pol. B 33, 45 (2002) INPUT:

  • type(particle), dimension(1:2), intent(in) :: pair ! incoming particles
  • real, intent(in) :: sigma ! not screened cross section
  • real, intent(in) :: srts ! c.m. energy (needed for ScreenMode=1 only)
  • type(medium), intent(in) :: mediumAtColl ! medium info (needed for ScreenMode=1 only)


XsectionRatios/getSigmaTotal [ Functions ]

[ Top ] [ XsectionRatios ] [ Functions ]

NAME

real function getSigmaTotal(pair)

PURPOSE

Returns the total in-medium cross section (mbarn). INPUT:

  • type(particle), dimension(1:2), intent(in) :: pair ! incoming particles

NOTES: The total in-medium pp cross section is used for all types of colliding baryons/antibaryons. Returns zero if there is at least one incoming meson.


XsectionRatios/accept_event [ Functions ]

[ Top ] [ XsectionRatios ] [ Functions ]

NAME

logical function accept_event(pair,finalState)

PURPOSE

Accepts or rejects the event randomly, using the in-medium cross section ratios * returns .true. --- if the event is accepted * returns .false. --- if the event is rejected INPUT:

  • type(particle), dimension(1:2), intent(in) :: pair ! incoming particles
  • type(particle), dimension(:), intent(in) :: finalState ! outgoing particles


XsectionRatios/init [ Subroutines ]

[ Top ] [ XsectionRatios ] [ Subroutines ]

NAME

subroutine init

PURPOSE

Reads in input switches from namelist XsectionRatios_input.


XsectionRatios/XsectionRatios_input [ Namelists ]

[ Top ] [ XsectionRatios ] [ Namelists ]

NAME

NAMELIST /XsectionRatios_input/

PURPOSE

Includes the switches:


XsectionRatios/Ratio [ Functions ]

[ Top ] [ XsectionRatios ] [ Functions ]

NAME

real function Ratio(srtsStar,mStar,mass)

PURPOSE

Compute the ratio of the in-medium and vacuum B_1 + B_2 -> B_3 + B_4 + M_5 + M_6 + M_7 + M_8 cross sections INPUT:

  • real, intent(in) :: srtsStar ! c.m. energy (GeV),
  • real, intent(in), dimension(:) :: mStar, mass ! Dirac and vacuum masses of all particles (GeV)

OUTPUT:

  • real :: Ratio ! ratio of the in-medium and vacuum cross sections

NOTES: The particles 1 and 2 are incoming baryons, the particles 3 and 4 are outgoing baryons, the particles M_5, ..., M_8 are outgoing mesons. The routine will compute one of the ratios: B_1 + B_2 -> B_3 + B_4, B_1 + B_2 -> B_3 + B_4 + M_5, ..., B_1 + B_2 -> B_3 + B_4 + M_5 + ...+ M_8 depending on the the size n of the arrays mStar and mass, which must be in the interval [4:8].


XsectionRatios/Ratio_BaB [ Functions ]

[ Top ] [ XsectionRatios ] [ Functions ]

NAME

real function Ratio_BaB(pair,finalState)

PURPOSE

Compute the ratio of the in-medium and vacuum cross sections baryon + antibaryon -> mesons INPUT:

  • type(particle), dimension(1:2), intent(in) :: pair ! incoming particles
  • type(particle), dimension(:), intent(in) :: finalState ! outgoing particles

OUTPUT:

  • real :: Ratio_BaB ! ratio of the in-medium and vacuum cross sections

NOTES: Maximum 6 mesons can be treated (otherwice returns Ratio_BaB=1.).