crosssections.h

Go to the documentation of this file.

/*
 *
 *    Copyright (c) 2015-2018
 *      SMASH Team
 *
 *    GNU General Public License (GPLv3 or later)
 *
 */

#ifndef SRC_INCLUDE_CROSSSECTIONS_H_
#define SRC_INCLUDE_CROSSSECTIONS_H_

#include "forwarddeclarations.h"
#include "isoparticletype.h"
#include "particles.h"
#include "potential_globals.h"
#include "processstring.h"

namespace smash {

namespace transit_high_energy {
const std::array<double, 2> sqrts_range_Npi = {1.9, 2.2};
const std::array<double, 2> sqrts_range_NN = {4., 5.};
const double sqrts_add_lower = 0.9;
const double sqrts_range = 1.0;

const double pipi_offset = 1.12;

const double KN_offset = 15.15;
}  // namespace transit_high_energy

class CrossSections {
 public:
  CrossSections(const ParticleList& incoming_particles, const double sqrt_s,
                const std::pair<FourVector, FourVector> potentials);

  CollisionBranchList generate_collision_list(
      double elastic_parameter, bool two_to_one_switch,
      ReactionsBitSet included_2to2, double low_snn_cut, bool strings_switch,
      bool use_AQM, bool strings_with_probability,
      NNbarTreatment nnbar_treatment, StringProcess* string_process) const;

  CollisionBranchPtr elastic(double elast_par, bool use_AQM) const;

  CollisionBranchList two_to_one() const;

  double formation(const ParticleType& type_resonance,
                   double cm_momentum_sqr) const;

  CollisionBranchList rare_two_to_two() const;

  CollisionBranchList two_to_two(ReactionsBitSet included_2to2) const;

  CollisionBranchList string_excitation(double total_string_xs,
                                        StringProcess* string_process,
                                        bool use_AQM) const;

  CollisionBranchPtr NNbar_annihilation(const double current_xs) const;

  CollisionBranchList NNbar_creation() const;

  double high_energy() const;

  double string_probability(bool strings_switch,
                            bool use_transition_probability, bool use_AQM,
                            bool treat_nnbar_with_strings) const;

  double probability_transit_high(const double region_lower,
                                  const double region_upper) const;

 private:
  double elastic_parametrization(bool use_AQM) const;

  double nn_el() const;

  double npi_el() const;

  double nk_el() const;

  CollisionBranchList npi_yk() const;

  CollisionBranchList bb_xx_except_nn(ReactionsBitSet included_2to2) const;

  CollisionBranchList nn_xx(ReactionsBitSet included_2to2) const;

  CollisionBranchList nk_xx(ReactionsBitSet included_2to2) const;

  CollisionBranchList deltak_xx(ReactionsBitSet included_2to2) const;

  CollisionBranchList ypi_xx(ReactionsBitSet included_2to2) const;

  CollisionBranchList dpi_xx(ReactionsBitSet included_2to2) const;

  CollisionBranchList dn_xx(ReactionsBitSet included_2to2) const;

  double string_hard_cross_section() const;

  CollisionBranchList bar_bar_to_nuc_nuc(const bool is_anti_particles) const;

  static double nn_to_resonance_matrix_element(double sqrts,
                                               const ParticleType& type_a,
                                               const ParticleType& type_b,
                                               const int twoI);

  template <class IntegrationMethod>
  CollisionBranchList find_nn_xsection_from_type(
      const ParticleTypePtrList& type_res_1,
      const ParticleTypePtrList& type_res_2,
      const IntegrationMethod integrator) const;

  double cm_momentum() const {
    const double m1 = incoming_particles_[0].effective_mass();
    const double m2 = incoming_particles_[1].effective_mass();
    return pCM(sqrt_s_, m1, m2);
  }

  const ParticleList incoming_particles_;

  const double sqrt_s_;

  const std::pair<FourVector, FourVector> potentials_;

  const bool is_BBbar_pair_;

  template <typename F>
  void add_channel(CollisionBranchList& process_list, F&& get_xsection,
                   double sqrts, const ParticleType& type_a,
                   const ParticleType& type_b) const {
    const double sqrt_s_min =
        type_a.min_mass_spectral() + type_b.min_mass_spectral();
    /* Determine wether the process is below the threshold. */
    double scale_B = 0.0;
    double scale_I3 = 0.0;
    bool is_below_threshold;
    FourVector incoming_momentum = FourVector();
    if (pot_pointer != nullptr) {
      for (const auto p : incoming_particles_) {
        incoming_momentum += p.momentum();
        scale_B += pot_pointer->force_scale(p.type()).first;
        scale_I3 +=
            pot_pointer->force_scale(p.type()).second * p.type().isospin3_rel();
      }
      scale_B -= pot_pointer->force_scale(type_a).first;
      scale_I3 -=
          pot_pointer->force_scale(type_a).second * type_a.isospin3_rel();
      scale_B -= pot_pointer->force_scale(type_b).first;
      scale_I3 -=
          pot_pointer->force_scale(type_b).second * type_b.isospin3_rel();
      is_below_threshold = (incoming_momentum + potentials_.first * scale_B +
                            potentials_.second * scale_I3)
                               .abs() <= sqrt_s_min;
    } else {
      is_below_threshold = (sqrts <= sqrt_s_min);
    }
    if (is_below_threshold) {
      return;
    }
    const auto xsection = get_xsection();
    if (xsection > really_small) {
      process_list.push_back(make_unique<CollisionBranch>(
          type_a, type_b, xsection, ProcessType::TwoToTwo));
    }
  }
};

}  // namespace smash

#endif  // SRC_INCLUDE_CROSSSECTIONS_H_