Version: SMASH-3.0
smash::Experiment< Modus > Class Template Reference

#include <experiment.h>

template<typename Modus>
class smash::Experiment< Modus >

The main class, where the simulation of an experiment is executed.

The Experiment class owns all data (maybe indirectly) relevant for the execution of the experiment simulation. The experiment can be conducted in different running modi. Since the abstraction of these differences should not incur any overhead, the design is built around the Policy pattern.

The Policy pattern was defined by Andrei Alexandrescu in his book "Modern C++ Design: Generic Programming and Design Patterns Applied". Addison-Wesley:

A policy defines a class interface or a class template interface. The interface consists of one or all of the following: inner type definitions, member functions, and member variables.

The policy pattern can also be understood as a compile-time variant of the strategy pattern.

The Modus template parameter defines the "policy" of the Experiment class. It determines several aspects of the experiment execution at compile time. The original strategy pattern would select these differences at run time, thus incurring an overhead. This overhead becomes severe in cases where calls to strategy/policy functions are done very frequently. Using the policy pattern, the compiler can fully optimize: It creates a new instance of all functions in Experiment for all different Modus types.

Definition at line 187 of file experiment.h.

Inheritance diagram for smash::Experiment< Modus >:
smash::ExperimentBase

Public Member Functions

void run () override
 Runs the experiment. More...
 
 Experiment (Configuration &config, const std::filesystem::path &output_path)
 Create a new Experiment. More...
 
void initialize_new_event ()
 This is called in the beginning of each event. More...
 
void run_time_evolution (const double t_end, ParticleList &&add_plist={}, ParticleList &&remove_plist={})
 Runs the time evolution of an event with fixed-size time steps or without timesteps, from action to actions. More...
 
void do_final_decays ()
 Performs the final decays of an event. More...
 
void final_output ()
 Output at the end of an event. More...
 
Particlesfirst_ensemble ()
 Provides external access to SMASH particles. More...
 
std::vector< Particles > * all_ensembles ()
 Getter for all ensembles. More...
 
Modus * modus ()
 Provides external access to SMASH calculation modus. More...
 
void increase_event_number ()
 Increases the event number by one. More...
 
- Public Member Functions inherited from smash::ExperimentBase
 ExperimentBase ()=default
 
virtual ~ExperimentBase ()=default
 The virtual destructor avoids undefined behavior when destroying derived objects. More...
 

Private Member Functions

bool perform_action (Action &action, int i_ensemble, bool include_pauli_blocking=true)
 Perform the given action. More...
 
void create_output (const std::string &format, const std::string &content, const std::filesystem::path &output_path, const OutputParameters &par)
 Create a list of output files. More...
 
void propagate_and_shine (double to_time, Particles &particles)
 Propagate all particles until time to_time without any interactions and shine dileptons. More...
 
void run_time_evolution_timestepless (Actions &actions, int i_ensemble, const double end_time_propagation)
 Performs all the propagations and actions during a certain time interval neglecting the influence of the potentials. More...
 
void intermediate_output ()
 Intermediate output during an event. More...
 
void update_potentials ()
 Recompute potentials on lattices if necessary. More...
 
double compute_min_cell_length (double dt) const
 Calculate the minimal size for the grid cells such that the ScatterActionsFinder will find all collisions within the maximal transverse distance (which is determined by the maximal cross section). More...
 
double next_output_time () const
 Shortcut for next output time. More...
 
void count_nonempty_ensembles ()
 Counts the number of ensembles in wich interactions took place at the end of an event. More...
 
bool is_finished ()
 Checks wether the desired number events have been calculated. More...
 

Private Attributes

ExperimentParameters parameters_
 Struct of several member variables. More...
 
DensityParameters density_param_
 Structure to precalculate and hold parameters for density computations. More...
 
Modus modus_
 Instance of the Modus template parameter. More...
 
std::vector< Particlesensembles_
 Complete particle list, all ensembles in one vector. More...
 
std::unique_ptr< Potentialspotentials_
 An instance of potentials class, that stores parameters of potentials, calculates them and their gradients. More...
 
std::unique_ptr< PauliBlockerpauli_blocker_
 An instance of PauliBlocker class that stores parameters needed for Pauli blocking calculations and computes phase-space density. More...
 
OutputsList outputs_
 A list of output formaters. More...
 
OutputPtr dilepton_output_
 The Dilepton output. More...
 
OutputPtr photon_output_
 The Photon output. More...
 
std::vector< bool > projectile_target_interact_
 Whether the projectile and the target collided. More...
 
std::vector< FourVectorbeam_momentum_ = {}
 The initial nucleons in the ColliderModus propagate with beam_momentum_, if Fermi motion is frozen. More...
 
std::vector< std::unique_ptr< ActionFinderInterface > > action_finders_
 The Action finder objects. More...
 
std::unique_ptr< DecayActionsFinderDileptondilepton_finder_
 The Dilepton Action Finder. More...
 
std::unique_ptr< ActionFinderInterfacephoton_finder_
 The (Scatter) Actions Finder for Direct Photons. More...
 
int n_fractional_photons_
 Number of fractional photons produced per single reaction. More...
 
std::unique_ptr< DensityLatticej_QBS_lat_
 4-current for j_QBS lattice output More...
 
std::unique_ptr< DensityLatticejmu_B_lat_
 Baryon density on the lattice. More...
 
std::unique_ptr< DensityLatticejmu_I3_lat_
 Isospin projection density on the lattice. More...
 
std::unique_ptr< DensityLatticejmu_el_lat_
 Electric charge density on the lattice. More...
 
std::unique_ptr< FieldsLatticefields_lat_
 Mean-field A^mu on the lattice. More...
 
std::unique_ptr< DensityLatticejmu_custom_lat_
 Custom density on the lattices. More...
 
DensityType dens_type_lattice_printout_ = DensityType::None
 Type of density for lattice printout. More...
 
std::unique_ptr< RectangularLattice< FourVector > > UB_lat_ = nullptr
 Lattices for Skyrme or VDF potentials (evaluated in the local rest frame) times the baryon flow 4-velocity. More...
 
std::unique_ptr< RectangularLattice< FourVector > > UI3_lat_ = nullptr
 Lattices for symmetry potentials (evaluated in the local rest frame) times the isospin flow 4-velocity. More...
 
std::unique_ptr< RectangularLattice< std::pair< ThreeVector, ThreeVector > > > FB_lat_
 Lattices for the electric and magnetic components of the Skyrme or VDF force. More...
 
std::unique_ptr< RectangularLattice< std::pair< ThreeVector, ThreeVector > > > FI3_lat_
 Lattices for the electric and magnetic component of the symmetry force. More...
 
std::unique_ptr< RectangularLattice< std::pair< ThreeVector, ThreeVector > > > EM_lat_
 Lattices for electric and magnetic field in fm^-2. More...
 
std::unique_ptr< RectangularLattice< EnergyMomentumTensor > > Tmn_
 Lattices of energy-momentum tensors for printout. More...
 
std::unique_ptr< RectangularLattice< FourVector > > old_jmu_auxiliary_
 Auxiliary lattice for values of jmu at a time step t0. More...
 
std::unique_ptr< RectangularLattice< FourVector > > new_jmu_auxiliary_
 Auxiliary lattice for values of jmu at a time step t0 + dt. More...
 
std::unique_ptr< RectangularLattice< std::array< FourVector, 4 > > > four_gradient_auxiliary_
 Auxiliary lattice for calculating the four-gradient of jmu. More...
 
std::unique_ptr< RectangularLattice< FourVector > > old_fields_auxiliary_
 Auxiliary lattice for values of Amu at a time step t0. More...
 
std::unique_ptr< RectangularLattice< FourVector > > new_fields_auxiliary_
 Auxiliary lattice for values of Amu at a time step t0 + dt. More...
 
std::unique_ptr< RectangularLattice< std::array< FourVector, 4 > > > fields_four_gradient_auxiliary_
 Auxiliary lattice for calculating the four-gradient of Amu. More...
 
bool printout_tmn_ = false
 Whether to print the energy-momentum tensor. More...
 
bool printout_tmn_landau_ = false
 Whether to print the energy-momentum tensor in Landau frame. More...
 
bool printout_v_landau_ = false
 Whether to print the 4-velocity in Landau frame. More...
 
bool printout_j_QBS_ = false
 Whether to print the Q, B, S 4-currents. More...
 
bool printout_lattice_td_ = false
 Whether to print the thermodynamics quantities evaluated on the lattices. More...
 
bool printout_full_lattice_ascii_td_ = false
 Whether to print the thermodynamics quantities evaluated on the lattices, point by point, in ASCII format. More...
 
bool printout_full_lattice_binary_td_ = false
 Whether to print the thermodynamics quantities evaluated on the lattices, point by point, in Binary format. More...
 
bool printout_full_lattice_any_td_ = false
 Whether to print the thermodynamics quantities evaluated on the lattices, point by point, in any format. More...
 
std::unique_ptr< GrandCanThermalizerthermalizer_
 Instance of class used for forced thermalization. More...
 
StringProcessprocess_string_ptr_
 Pointer to the string process class object, which is used to set the random seed for PYTHIA objects in each event. More...
 
const int nevents_ = 0
 Number of events. More...
 
int minimum_nonempty_ensembles_ = 0
 The number of ensembles, in which interactions take place, to be calculated. More...
 
EventCounting event_counting_ = EventCounting::Invalid
 The way in which the number of calculated events is specified. More...
 
int event_ = 0
 Current event. More...
 
int nonempty_ensembles_ = 0
 Number of ensembles containing an interaction. More...
 
int max_events_ = 0
 Maximum number of events to be calculated in order obtain the desired number of non-empty events using the MinimumNonemptyEnsembles option. More...
 
const double end_time_
 simulation time at which the evolution is stopped. More...
 
const double delta_time_startup_
 The clock's timestep size at start up. More...
 
const bool force_decays_
 This indicates whether we force all resonances to decay in the last timestep. More...
 
const bool use_grid_
 This indicates whether to use the grid. More...
 
const ExpansionProperties metric_
 This struct contains information on the metric to be used. More...
 
const bool dileptons_switch_
 This indicates whether dileptons are switched on. More...
 
const bool photons_switch_
 This indicates whether photons are switched on. More...
 
const bool bremsstrahlung_switch_
 This indicates whether bremsstrahlung is switched on. More...
 
const bool IC_output_switch_
 This indicates whether the IC output is enabled. More...
 
const TimeStepMode time_step_mode_
 This indicates whether to use time steps. More...
 
double max_transverse_distance_sqr_ = std::numeric_limits<double>::max()
 Maximal distance at which particles can interact in case of the geometric criterion, squared. More...
 
QuantumNumbers conserved_initial_
 The conserved quantities of the system. More...
 
double initial_mean_field_energy_
 The initial total mean field energy in the system. More...
 
SystemTimePoint time_start_ = SystemClock::now()
 system starting time of the simulation More...
 
DensityType dens_type_ = DensityType::None
 Type of density to be written to collision headers. More...
 
uint64_t interactions_total_ = 0
 Total number of interactions for current timestep. More...
 
uint64_t previous_interactions_total_ = 0
 Total number of interactions for previous timestep. More...
 
uint64_t wall_actions_total_ = 0
 Total number of wall-crossings for current timestep. More...
 
uint64_t previous_wall_actions_total_ = 0
 Total number of wall-crossings for previous timestep. More...
 
uint64_t total_pauli_blocked_ = 0
 Total number of Pauli-blockings for current timestep. More...
 
uint64_t total_hypersurface_crossing_actions_ = 0
 Total number of particles removed from the evolution in hypersurface crossing actions. More...
 
uint64_t discarded_interactions_total_ = 0
 Total number of discarded interactions, because they were invalidated before they could be performed. More...
 
double total_energy_removed_ = 0.0
 Total energy removed from the system in hypersurface crossing actions. More...
 
double total_energy_violated_by_Pythia_ = 0.0
 Total energy violation introduced by Pythia. More...
 
bool kinematic_cuts_for_IC_output_ = false
 This indicates whether kinematic cuts are enabled for the IC output. More...
 
int64_t seed_ = -1
 random seed for the next event. More...
 

Friends

class ExperimentBase
 
std::ostream & operator<< (std::ostream &out, const Experiment &e)
 Writes the initial state for the Experiment to the output stream. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from smash::ExperimentBase
static std::unique_ptr< ExperimentBasecreate (Configuration &config, const std::filesystem::path &output_path)
 Factory method that creates and initializes a new Experiment<Modus>. More...
 

Constructor & Destructor Documentation

◆ Experiment()

template<typename Modus >
smash::Experiment< Modus >::Experiment ( Configuration config,
const std::filesystem::path &  output_path 
)
explicit

Create a new Experiment.

This constructor is only called from the ExperimentBase::create factory method.

Parameters
[in,out]configThe Configuration object contains all initial setup of the experiment. It is forwarded to the constructors of member variables as needed. Note that the object is passed by non-const reference. This is only necessary for bookkeeping: Values are not only read, but actually taken out of the object. Thus, all values that remain were not used.
[in]output_pathThe directory where the output files are written.

Member Function Documentation

◆ run()

template<typename Modus >
void smash::Experiment< Modus >::run
overridevirtual

Runs the experiment.

The constructor does the setup of the experiment. The run function executes the complete experiment.

Implements smash::ExperimentBase.

Definition at line 2974 of file experiment.h.

2974  {
2975  const auto &mainlog = logg[LMain];
2976  for (event_ = 0; !is_finished(); event_++) {
2977  mainlog.info() << "Event " << event_;
2978 
2979  // Sample initial particles, start clock, some printout and book-keeping
2981 
2983 
2984  if (force_decays_) {
2985  do_final_decays();
2986  }
2987 
2988  // Output at event end
2989  final_output();
2990  }
2991 }
const bool force_decays_
This indicates whether we force all resonances to decay in the last timestep.
Definition: experiment.h:599
void run_time_evolution(const double t_end, ParticleList &&add_plist={}, ParticleList &&remove_plist={})
Runs the time evolution of an event with fixed-size time steps or without timesteps,...
Definition: experiment.h:2224
void initialize_new_event()
This is called in the beginning of each event.
Definition: experiment.h:1860
bool is_finished()
Checks wether the desired number events have been calculated.
Definition: experiment.h:2943
int event_
Current event.
Definition: experiment.h:574
void do_final_decays()
Performs the final decays of an event.
Definition: experiment.h:2733
void final_output()
Output at the end of an event.
Definition: experiment.h:2782
const double end_time_
simulation time at which the evolution is stopped.
Definition: experiment.h:586
std::array< einhard::Logger<>, std::tuple_size< LogArea::AreaTuple >::value > logg
An array that stores all pre-configured Logger objects.
Definition: logging.cc:39
static constexpr int LMain
Definition: experiment.h:88

◆ initialize_new_event()

template<typename Modus >
void smash::Experiment< Modus >::initialize_new_event

This is called in the beginning of each event.

It initializes particles according to selected modus, resets the clock and saves the initial conserved quantities for subsequent sanity checks.

Definition at line 1860 of file experiment.h.

1860  {
1862  logg[LExperiment].info() << "random number seed: " << seed_;
1863  /* Set seed for the next event. It has to be positive, so it can be entered
1864  * in the config.
1865  *
1866  * We have to be careful about the minimal integer, whose absolute value
1867  * cannot be represented. */
1868  int64_t r = random::advance();
1869  while (r == INT64_MIN) {
1870  r = random::advance();
1871  }
1872  seed_ = std::abs(r);
1873  /* Set the random seed used in PYTHIA hadronization
1874  * to be same with the SMASH one.
1875  * In this way we ensure that the results are reproducible
1876  * for every event if one knows SMASH random seed. */
1877  if (process_string_ptr_ != NULL) {
1879  }
1880 
1881  for (Particles &particles : ensembles_) {
1882  particles.reset();
1883  }
1884 
1885  // Sample particles according to the initial conditions
1886  double start_time = -1.0;
1887 
1888  // Sample impact parameter only once per all ensembles
1889  // It should be the same for all ensembles
1890  if (modus_.is_collider()) {
1891  modus_.sample_impact();
1892  logg[LExperiment].info("Impact parameter = ", modus_.impact_parameter(),
1893  " fm");
1894  }
1895  for (Particles &particles : ensembles_) {
1896  start_time = modus_.initial_conditions(&particles, parameters_);
1897  }
1898  /* For box modus make sure that particles are in the box. In principle, after
1899  * a correct initialization they should be, so this is just playing it safe.
1900  */
1901  for (Particles &particles : ensembles_) {
1902  modus_.impose_boundary_conditions(&particles, outputs_);
1903  }
1904  // Reset the simulation clock
1905  double timestep = delta_time_startup_;
1906 
1907  switch (time_step_mode_) {
1908  case TimeStepMode::Fixed:
1909  break;
1910  case TimeStepMode::None:
1911  timestep = end_time_ - start_time;
1912  // Take care of the box modus + timestepless propagation
1913  const double max_dt = modus_.max_timestep(max_transverse_distance_sqr_);
1914  if (max_dt > 0. && max_dt < timestep) {
1915  timestep = max_dt;
1916  }
1917  break;
1918  }
1919  std::unique_ptr<UniformClock> clock_for_this_event;
1920  if (modus_.is_list() && (timestep < 0.0)) {
1921  throw std::runtime_error(
1922  "Timestep for the given event is negative. \n"
1923  "This might happen if the formation times of the input particles are "
1924  "larger than the specified end time of the simulation.");
1925  }
1926  clock_for_this_event =
1927  std::make_unique<UniformClock>(start_time, timestep, end_time_);
1928  parameters_.labclock = std::move(clock_for_this_event);
1929 
1930  // Reset the output clock
1931  parameters_.outputclock->reset(start_time, true);
1932  // remove time before starting time in case of custom output times.
1933  parameters_.outputclock->remove_times_in_past(start_time);
1934 
1935  logg[LExperiment].debug(
1936  "Lab clock: t_start = ", parameters_.labclock->current_time(),
1937  ", dt = ", parameters_.labclock->timestep_duration());
1938 
1939  /* Save the initial conserved quantum numbers and total momentum in
1940  * the system for conservation checks */
1941  conserved_initial_ = QuantumNumbers(ensembles_);
1942  wall_actions_total_ = 0;
1944  interactions_total_ = 0;
1950  total_energy_removed_ = 0.0;
1952  // Print output headers
1953  logg[LExperiment].info() << hline;
1954  logg[LExperiment].info() << "Time[fm] Ekin[GeV] E_MF[GeV] ETotal[GeV] "
1955  << "ETot/N[GeV] D(ETot/N)[GeV] Scatt&Decays "
1956  << "Particles Comp.Time";
1957  logg[LExperiment].info() << hline;
1958  double E_mean_field = 0.0;
1959  if (potentials_) {
1960  // update_potentials();
1961  // if (parameters.outputclock->current_time() == 0.0 )
1962  // using the lattice is necessary
1963  if ((jmu_B_lat_ != nullptr)) {
1968  parameters_.labclock->timestep_duration(), true);
1969  // Because there was no lattice at t=-Delta_t, the time derivatives
1970  // drho_dt and dj^mu/dt at t=0 are huge, while they shouldn't be; we
1971  // overwrite the time derivative to zero by hand.
1972  for (auto &node : *jmu_B_lat_) {
1973  node.overwrite_drho_dt_to_zero();
1974  node.overwrite_djmu_dt_to_zero();
1975  }
1977  EM_lat_.get(), parameters_);
1978  }
1979  }
1980  initial_mean_field_energy_ = E_mean_field;
1983  parameters_.labclock->current_time(), E_mean_field,
1985 
1986  // Output at event start
1987  for (const auto &output : outputs_) {
1988  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
1989  auto event_info = fill_event_info(
1990  ensembles_, E_mean_field, modus_.impact_parameter(), parameters_,
1992  output->at_eventstart(ensembles_[i_ens],
1993  // Pretend each ensemble is an independent event
1994  event_ * parameters_.n_ensembles + i_ens,
1995  event_info);
1996  }
1997  // For thermodynamic output
1998  output->at_eventstart(ensembles_, event_);
1999  // For thermodynamic lattice output
2001  switch (dens_type_lattice_printout_) {
2002  case DensityType::Baryon:
2003  output->at_eventstart(event_, ThermodynamicQuantity::EckartDensity,
2005  break;
2007  output->at_eventstart(event_, ThermodynamicQuantity::EckartDensity,
2009  break;
2010  case DensityType::None:
2011  break;
2012  default:
2013  output->at_eventstart(event_, ThermodynamicQuantity::EckartDensity,
2015  }
2016  if (printout_tmn_) {
2017  output->at_eventstart(event_, ThermodynamicQuantity::Tmn,
2019  }
2020  if (printout_tmn_landau_) {
2021  output->at_eventstart(event_, ThermodynamicQuantity::TmnLandau,
2023  }
2024  if (printout_v_landau_) {
2025  output->at_eventstart(event_, ThermodynamicQuantity::LandauVelocity,
2027  }
2028  if (printout_j_QBS_) {
2029  output->at_eventstart(event_, ThermodynamicQuantity::j_QBS,
2031  }
2032  }
2033  }
2034 
2035  /* In the ColliderModus, if Fermi motion is frozen, assign the beam momenta
2036  * to the nucleons in both the projectile and the target. Every ensemble
2037  * gets the same beam momenta, so no need to create beam_momenta_ vector
2038  * for every ensemble.
2039  */
2040  if (modus_.is_collider() && modus_.fermi_motion() == FermiMotion::Frozen) {
2041  for (ParticleData &particle : ensembles_[0]) {
2042  const double m = particle.effective_mass();
2043  double v_beam = 0.0;
2044  if (particle.belongs_to() == BelongsTo::Projectile) {
2045  v_beam = modus_.velocity_projectile();
2046  } else if (particle.belongs_to() == BelongsTo::Target) {
2047  v_beam = modus_.velocity_target();
2048  }
2049  const double gamma = 1.0 / std::sqrt(1.0 - v_beam * v_beam);
2050  beam_momentum_.emplace_back(
2051  FourVector(gamma * m, 0.0, 0.0, gamma * v_beam * m));
2052  } // loop over particles
2053  }
2054 }
double initial_mean_field_energy_
The initial total mean field energy in the system.
Definition: experiment.h:641
bool printout_tmn_
Whether to print the energy-momentum tensor.
Definition: experiment.h:504
QuantumNumbers conserved_initial_
The conserved quantities of the system.
Definition: experiment.h:635
DensityParameters density_param_
Structure to precalculate and hold parameters for density computations.
Definition: experiment.h:370
double total_energy_removed_
Total energy removed from the system in hypersurface crossing actions.
Definition: experiment.h:694
DensityType dens_type_lattice_printout_
Type of density for lattice printout.
Definition: experiment.h:455
double max_transverse_distance_sqr_
Maximal distance at which particles can interact in case of the geometric criterion,...
Definition: experiment.h:626
bool printout_j_QBS_
Whether to print the Q, B, S 4-currents.
Definition: experiment.h:513
const TimeStepMode time_step_mode_
This indicates whether to use time steps.
Definition: experiment.h:620
std::unique_ptr< DensityLattice > jmu_custom_lat_
Custom density on the lattices.
Definition: experiment.h:452
bool printout_full_lattice_any_td_
Whether to print the thermodynamics quantities evaluated on the lattices, point by point,...
Definition: experiment.h:528
std::unique_ptr< DensityLattice > jmu_B_lat_
Baryon density on the lattice.
Definition: experiment.h:434
std::vector< FourVector > beam_momentum_
The initial nucleons in the ColliderModus propagate with beam_momentum_, if Fermi motion is frozen.
Definition: experiment.h:416
std::unique_ptr< RectangularLattice< FourVector > > new_jmu_auxiliary_
Auxiliary lattice for values of jmu at a time step t0 + dt.
Definition: experiment.h:490
const double delta_time_startup_
The clock's timestep size at start up.
Definition: experiment.h:593
std::unique_ptr< RectangularLattice< EnergyMomentumTensor > > Tmn_
Lattices of energy-momentum tensors for printout.
Definition: experiment.h:485
std::vector< Particles > ensembles_
Complete particle list, all ensembles in one vector.
Definition: experiment.h:379
SystemTimePoint time_start_
system starting time of the simulation
Definition: experiment.h:644
uint64_t previous_wall_actions_total_
Total number of wall-crossings for previous timestep.
Definition: experiment.h:671
OutputsList outputs_
A list of output formaters.
Definition: experiment.h:397
bool printout_v_landau_
Whether to print the 4-velocity in Landau frame.
Definition: experiment.h:510
std::unique_ptr< RectangularLattice< std::pair< ThreeVector, ThreeVector > > > EM_lat_
Lattices for electric and magnetic field in fm^-2.
Definition: experiment.h:482
Modus modus_
Instance of the Modus template parameter.
Definition: experiment.h:376
std::unique_ptr< RectangularLattice< FourVector > > old_jmu_auxiliary_
Auxiliary lattice for values of jmu at a time step t0.
Definition: experiment.h:488
std::unique_ptr< DensityLattice > j_QBS_lat_
4-current for j_QBS lattice output
Definition: experiment.h:431
bool printout_tmn_landau_
Whether to print the energy-momentum tensor in Landau frame.
Definition: experiment.h:507
std::unique_ptr< RectangularLattice< std::array< FourVector, 4 > > > four_gradient_auxiliary_
Auxiliary lattice for calculating the four-gradient of jmu.
Definition: experiment.h:493
StringProcess * process_string_ptr_
Pointer to the string process class object, which is used to set the random seed for PYTHIA objects i...
Definition: experiment.h:537
ExperimentParameters parameters_
Struct of several member variables.
Definition: experiment.h:367
std::unique_ptr< Potentials > potentials_
An instance of potentials class, that stores parameters of potentials, calculates them and their grad...
Definition: experiment.h:385
uint64_t wall_actions_total_
Total number of wall-crossings for current timestep.
Definition: experiment.h:665
uint64_t total_hypersurface_crossing_actions_
Total number of particles removed from the evolution in hypersurface crossing actions.
Definition: experiment.h:683
uint64_t interactions_total_
Total number of interactions for current timestep.
Definition: experiment.h:653
std::unique_ptr< DensityLattice > jmu_I3_lat_
Isospin projection density on the lattice.
Definition: experiment.h:437
uint64_t total_pauli_blocked_
Total number of Pauli-blockings for current timestep.
Definition: experiment.h:677
std::vector< bool > projectile_target_interact_
Whether the projectile and the target collided.
Definition: experiment.h:409
int64_t seed_
random seed for the next event.
Definition: experiment.h:705
uint64_t previous_interactions_total_
Total number of interactions for previous timestep.
Definition: experiment.h:659
uint64_t discarded_interactions_total_
Total number of discarded interactions, because they were invalidated before they could be performed.
Definition: experiment.h:689
bool kinematic_cuts_for_IC_output_
This indicates whether kinematic cuts are enabled for the IC output.
Definition: experiment.h:702
double total_energy_violated_by_Pythia_
Total energy violation introduced by Pythia.
Definition: experiment.h:699
void init_pythia_hadron_rndm()
Set PYTHIA random seeds to be desired values.
@ Frozen
Use fermi motion without potentials.
@ Fixed
Use fixed time step.
@ None
Don't use time steps; propagate from action to action.
Engine::result_type advance()
Advance the engine's state and return the generated value.
Definition: random.h:78
void set_seed(T &&seed)
Sets the seed of the random number engine.
Definition: random.h:71
EventInfo fill_event_info(const std::vector< Particles > &ensembles, double E_mean_field, double modus_impact_parameter, const ExperimentParameters &parameters, bool projectile_target_interact, bool kinematic_cut_for_SMASH_IC)
Generate the EventInfo object which is passed to outputs_.
Definition: experiment.cc:604
std::string format_measurements(const std::vector< Particles > &ensembles, uint64_t scatterings_this_interval, const QuantumNumbers &conserved_initial, SystemTimePoint time_start, double time, double E_mean_field, double E_mean_field_initial)
Generate a string which will be printed to the screen when SMASH is running.
Definition: experiment.cc:336
void update_lattice(RectangularLattice< T > *lat, const LatticeUpdate update, const DensityType dens_type, const DensityParameters &par, const std::vector< Particles > &ensembles, const bool compute_gradient)
Updates the contents on the lattice.
Definition: density.h:542
double calculate_mean_field_energy(const Potentials &potentials, RectangularLattice< smash::DensityOnLattice > &jmu_B_lat, RectangularLattice< std::pair< ThreeVector, ThreeVector >> *em_lattice, const ExperimentParameters &parameters)
Calculate the total mean field energy of the system; this will be printed to the screen when SMASH is...
Definition: experiment.cc:384
static constexpr int LExperiment
const std::string hline(113, '-')
String representing a horizontal line.
int n_ensembles
Number of parallel ensembles.
std::unique_ptr< Clock > outputclock
Output clock to keep track of the next output time.
std::unique_ptr< Clock > labclock
System clock (for simulation time keeping in the computational frame)

◆ run_time_evolution()

template<typename Modus >
void smash::Experiment< Modus >::run_time_evolution ( const double  t_end,
ParticleList &&  add_plist = {},
ParticleList &&  remove_plist = {} 
)

Runs the time evolution of an event with fixed-size time steps or without timesteps, from action to actions.

Within one timestep (fixed) evolution from action to action is invoked.

Parameters
[in]t_endTime until run_time_evolution is run, in SMASH this is the configured end_time, but it might differ if SMASH is used as an external library
[in]add_plistA by-default empty particle list which is added to the current particle content of the system
[in]remove_plistA by-default empty particle list which is removed from the current particle content of the system
Note
This function is meant to take over ownership of the to-be-added/removed particle lists and that's why these are passed by rvalue reference.

Definition at line 2224 of file experiment.h.

2226  {
2227  if (!add_plist.empty() || !remove_plist.empty()) {
2228  if (ensembles_.size() > 1) {
2229  throw std::runtime_error(
2230  "Adding or removing particles from SMASH is only possible when one "
2231  "ensemble is used.");
2232  }
2233  const double action_time = parameters_.labclock->current_time();
2234  if (!add_plist.empty()) {
2236  }
2237  if (!add_plist.empty()) {
2238  // Create and perform action to add particle(s)
2239  auto action_add_particles = std::make_unique<FreeforallAction>(
2240  ParticleList{}, add_plist, action_time);
2241  perform_action(*action_add_particles, 0);
2242  }
2243  if (!remove_plist.empty()) {
2244  validate_and_adjust_particle_list(remove_plist);
2245  const auto number_of_particles_to_be_removed = remove_plist.size();
2246  remove_plist.erase(
2247  std::remove_if(
2248  remove_plist.begin(), remove_plist.end(),
2249  [this, &action_time](const ParticleData &particle_to_remove) {
2250  return std::find_if(
2251  ensembles_[0].begin(), ensembles_[0].end(),
2252  [&particle_to_remove,
2253  &action_time](const ParticleData &p) {
2254  return are_particles_identical_at_given_time(
2255  particle_to_remove, p, action_time);
2256  }) == ensembles_[0].end();
2257  }),
2258  remove_plist.end());
2259  if (auto delta = number_of_particles_to_be_removed - remove_plist.size();
2260  delta > 0) {
2261  logg[LExperiment].warn(
2262  "When trying to remove particle(s) at the beginning ",
2263  "of the system evolution,\n", delta,
2264  " particle(s) could not be found and will be ignored.");
2265  }
2266  }
2267  if (!remove_plist.empty()) {
2268  // Create and perform action to remove particles
2269  auto action_remove_particles = std::make_unique<FreeforallAction>(
2270  remove_plist, ParticleList{}, action_time);
2271  perform_action(*action_remove_particles, 0);
2272  }
2273  }
2274 
2275  while (parameters_.labclock->current_time() < t_end) {
2276  const double dt = parameters_.labclock->timestep_duration();
2277  logg[LExperiment].debug("Timestepless propagation for next ", dt, " fm.");
2278 
2279  // Perform forced thermalization if required
2280  if (thermalizer_ &&
2281  thermalizer_->is_time_to_thermalize(parameters_.labclock)) {
2282  const bool ignore_cells_under_treshold = true;
2283  // Thermodynamics in thermalizer is computed from all ensembles,
2284  // but thermalization actions act on each ensemble independently
2285  thermalizer_->update_thermalizer_lattice(ensembles_, density_param_,
2286  ignore_cells_under_treshold);
2287  const double current_t = parameters_.labclock->current_time();
2288  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2289  thermalizer_->thermalize(ensembles_[i_ens], current_t,
2291  ThermalizationAction th_act(*thermalizer_, current_t);
2292  if (th_act.any_particles_thermalized()) {
2293  perform_action(th_act, i_ens);
2294  }
2295  }
2296  }
2297 
2298  std::vector<Actions> actions(parameters_.n_ensembles);
2299  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2300  actions[i_ens].clear();
2301  if (ensembles_[i_ens].size() > 0 && action_finders_.size() > 0) {
2302  /* (1.a) Create grid. */
2303  const double min_cell_length = compute_min_cell_length(dt);
2304  logg[LExperiment].debug("Creating grid with minimal cell length ",
2305  min_cell_length);
2306  /* For the hyper-surface-crossing actions also unformed particles are
2307  * searched and therefore needed on the grid. */
2308  const bool include_unformed_particles = IC_output_switch_;
2309  const auto &grid =
2310  use_grid_ ? modus_.create_grid(ensembles_[i_ens], min_cell_length,
2311  dt, parameters_.coll_crit,
2312  include_unformed_particles)
2313  : modus_.create_grid(ensembles_[i_ens], min_cell_length,
2314  dt, parameters_.coll_crit,
2315  include_unformed_particles,
2317 
2318  const double gcell_vol = grid.cell_volume();
2319  /* (1.b) Iterate over cells and find actions. */
2320  grid.iterate_cells(
2321  [&](const ParticleList &search_list) {
2322  for (const auto &finder : action_finders_) {
2323  actions[i_ens].insert(finder->find_actions_in_cell(
2324  search_list, dt, gcell_vol, beam_momentum_));
2325  }
2326  },
2327  [&](const ParticleList &search_list,
2328  const ParticleList &neighbors_list) {
2329  for (const auto &finder : action_finders_) {
2330  actions[i_ens].insert(finder->find_actions_with_neighbors(
2331  search_list, neighbors_list, dt, beam_momentum_));
2332  }
2333  });
2334  }
2335  }
2336 
2337  /* \todo (optimizations) Adapt timestep size here */
2338 
2339  /* (2) Propagate from action to action until next output or timestep end */
2340  const double end_timestep_time = parameters_.labclock->next_time();
2341  while (next_output_time() < end_timestep_time) {
2342  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2343  run_time_evolution_timestepless(actions[i_ens], i_ens,
2344  next_output_time());
2345  }
2346  ++(*parameters_.outputclock);
2347 
2349  }
2350  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2351  run_time_evolution_timestepless(actions[i_ens], i_ens, end_timestep_time);
2352  }
2353 
2354  /* (3) Update potentials (if computed on the lattice) and
2355  * compute new momenta according to equations of motion */
2356  if (potentials_) {
2358  update_momenta(ensembles_, parameters_.labclock->timestep_duration(),
2359  *potentials_, FB_lat_.get(), FI3_lat_.get(),
2360  EM_lat_.get());
2361  }
2362 
2363  /* (4) Expand universe if non-minkowskian metric; updates
2364  * positions and momenta according to the selected expansion */
2366  for (Particles &particles : ensembles_) {
2367  expand_space_time(&particles, parameters_, metric_);
2368  }
2369  }
2370 
2371  ++(*parameters_.labclock);
2372 
2373  /* (5) Check conservation laws.
2374  *
2375  * Check conservation of conserved quantities if potentials and string
2376  * fragmentation are off. If potentials are on then momentum is conserved
2377  * only in average. If string fragmentation is on, then energy and
2378  * momentum are only very roughly conserved in high-energy collisions. */
2381  std::string err_msg = conserved_initial_.report_deviations(ensembles_);
2382  if (!err_msg.empty()) {
2383  logg[LExperiment].error() << err_msg;
2384  throw std::runtime_error("Violation of conserved quantities!");
2385  }
2386  }
2387  }
2388 
2389  if (pauli_blocker_) {
2390  logg[LExperiment].info(
2391  "Interactions: Pauli-blocked/performed = ", total_pauli_blocked_, "/",
2393  }
2394 }
const ExpansionProperties metric_
This struct contains information on the metric to be used.
Definition: experiment.h:605
std::vector< std::unique_ptr< ActionFinderInterface > > action_finders_
The Action finder objects.
Definition: experiment.h:419
double next_output_time() const
Shortcut for next output time.
Definition: experiment.h:345
std::unique_ptr< GrandCanThermalizer > thermalizer_
Instance of class used for forced thermalization.
Definition: experiment.h:531
void run_time_evolution_timestepless(Actions &actions, int i_ensemble, const double end_time_propagation)
Performs all the propagations and actions during a certain time interval neglecting the influence of ...
Definition: experiment.h:2423
void intermediate_output()
Intermediate output during an event.
Definition: experiment.h:2487
std::unique_ptr< RectangularLattice< std::pair< ThreeVector, ThreeVector > > > FI3_lat_
Lattices for the electric and magnetic component of the symmetry force.
Definition: experiment.h:478
double compute_min_cell_length(double dt) const
Calculate the minimal size for the grid cells such that the ScatterActionsFinder will find all collis...
Definition: experiment.h:337
std::unique_ptr< RectangularLattice< std::pair< ThreeVector, ThreeVector > > > FB_lat_
Lattices for the electric and magnetic components of the Skyrme or VDF force.
Definition: experiment.h:474
std::unique_ptr< PauliBlocker > pauli_blocker_
An instance of PauliBlocker class that stores parameters needed for Pauli blocking calculations and c...
Definition: experiment.h:391
bool perform_action(Action &action, int i_ensemble, bool include_pauli_blocking=true)
Perform the given action.
const bool use_grid_
This indicates whether to use the grid.
Definition: experiment.h:602
void update_potentials()
Recompute potentials on lattices if necessary.
Definition: experiment.h:2633
const bool IC_output_switch_
This indicates whether the IC output is enabled.
Definition: experiment.h:617
std::string report_deviations(const std::vector< Particles > &ensembles) const
Checks if the current particle list has still the same values and reports about differences.
void update_momenta(std::vector< Particles > &particles, double dt, const Potentials &pot, RectangularLattice< std::pair< ThreeVector, ThreeVector >> *FB_lat, RectangularLattice< std::pair< ThreeVector, ThreeVector >> *FI3_lat, RectangularLattice< std::pair< ThreeVector, ThreeVector >> *EM_lat)
Updates the momenta of all particles at the current time step according to the equations of motion:
Definition: propagation.cc:111
void expand_space_time(Particles *particles, const ExperimentParameters &parameters, const ExpansionProperties &metric)
Modifies positions and momentum of all particles to account for space-time deformation.
Definition: propagation.cc:86
void validate_and_adjust_particle_list(ParticleList &particle_list)
Validate a particle list adjusting each particle to be a valid SMASH particle.
Definition: experiment.cc:626
@ Largest
Make cells as large as possible.
ExpansionMode mode_
Type of metric used.
Definition: propagation.h:28
const CollisionCriterion coll_crit
Employed collision criterion.
bool strings_switch
This indicates whether string fragmentation is switched on.
int testparticles
Number of test-particles.

◆ do_final_decays()

template<typename Modus >
void smash::Experiment< Modus >::do_final_decays

Performs the final decays of an event.

Exceptions
runtime_errorif found actions cannot be performed

Definition at line 2733 of file experiment.h.

2733  {
2734  /* At end of time evolution: Force all resonances to decay. In order to handle
2735  * decay chains, we need to loop until no further actions occur. */
2736  bool actions_performed, decays_found;
2737  uint64_t interactions_old;
2738  do {
2739  decays_found = false;
2740  interactions_old = interactions_total_;
2741  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2742  Actions actions;
2743 
2744  // Dileptons: shining of remaining resonances
2745  if (dilepton_finder_ != nullptr) {
2746  for (const auto &output : outputs_) {
2747  dilepton_finder_->shine_final(ensembles_[i_ens], output.get(), true);
2748  }
2749  }
2750  // Find actions.
2751  for (const auto &finder : action_finders_) {
2752  auto found_actions = finder->find_final_actions(ensembles_[i_ens]);
2753  if (!found_actions.empty()) {
2754  actions.insert(std::move(found_actions));
2755  decays_found = true;
2756  }
2757  }
2758  // Perform actions.
2759  while (!actions.is_empty()) {
2760  perform_action(*actions.pop(), i_ens, false);
2761  }
2762  }
2763  actions_performed = interactions_total_ > interactions_old;
2764  // Throw an error if actions were found but not performed
2765  if (decays_found && !actions_performed) {
2766  throw std::runtime_error("Final decays were found but not performed.");
2767  }
2768  // loop until no more decays occur
2769  } while (actions_performed);
2770 
2771  // Dileptons: shining of stable particles at the end
2772  if (dilepton_finder_ != nullptr) {
2773  for (const auto &output : outputs_) {
2774  for (Particles &particles : ensembles_) {
2775  dilepton_finder_->shine_final(particles, output.get(), false);
2776  }
2777  }
2778  }
2779 }
std::unique_ptr< DecayActionsFinderDilepton > dilepton_finder_
The Dilepton Action Finder.
Definition: experiment.h:422

◆ final_output()

template<typename Modus >
void smash::Experiment< Modus >::final_output

Output at the end of an event.

Definition at line 2782 of file experiment.h.

2782  {
2783  /* make sure the experiment actually ran (note: we should compare this
2784  * to the start time, but we don't know that. Therefore, we check that
2785  * the time is positive, which should heuristically be the same). */
2786  double E_mean_field = 0.0;
2787  if (likely(parameters_.labclock > 0)) {
2788  const uint64_t wall_actions_this_interval =
2790  const uint64_t interactions_this_interval = interactions_total_ -
2792  wall_actions_this_interval;
2793  if (potentials_) {
2794  // using the lattice is necessary
2795  if ((jmu_B_lat_ != nullptr)) {
2797  EM_lat_.get(), parameters_);
2798  }
2799  }
2800  if (std::abs(parameters_.labclock->current_time() - end_time_) >
2801  really_small) {
2802  logg[LExperiment].warn()
2803  << "SMASH not propagated until configured end time. Current time = "
2804  << parameters_.labclock->current_time()
2805  << "fm. End time = " << end_time_ << "fm.";
2806  } else {
2808  ensembles_, interactions_this_interval, conserved_initial_,
2810  }
2811  int total_particles = 0;
2812  for (const Particles &particles : ensembles_) {
2813  total_particles += particles.size();
2814  }
2815  if (IC_output_switch_ && (total_particles == 0)) {
2816  const double initial_system_energy_plus_Pythia_violations =
2818  const double fraction_of_total_system_energy_removed =
2819  initial_system_energy_plus_Pythia_violations / total_energy_removed_;
2820  // Verify there is no more energy in the system if all particles were
2821  // removed when crossing the hypersurface
2822  if (std::fabs(fraction_of_total_system_energy_removed - 1.) >
2823  really_small) {
2824  throw std::runtime_error(
2825  "There is remaining energy in the system although all particles "
2826  "were removed.\n"
2827  "E_remain = " +
2828  std::to_string((initial_system_energy_plus_Pythia_violations -
2830  " [GeV]");
2831  } else {
2832  logg[LExperiment].info() << hline;
2833  logg[LExperiment].info()
2834  << "Time real: " << SystemClock::now() - time_start_;
2835  logg[LExperiment].info()
2836  << "Interactions before reaching hypersurface: "
2839  logg[LExperiment].info()
2840  << "Total number of particles removed on hypersurface: "
2842  }
2843  } else {
2844  const double precent_discarded =
2846  ? static_cast<double>(discarded_interactions_total_) * 100.0 /
2848  : 0.0;
2849  std::stringstream msg_discarded;
2850  msg_discarded
2851  << "Discarded interaction number: " << discarded_interactions_total_
2852  << " (" << precent_discarded
2853  << "% of the total interaction number including wall crossings)";
2854 
2855  logg[LExperiment].info() << hline;
2856  logg[LExperiment].info()
2857  << "Time real: " << SystemClock::now() - time_start_;
2858  logg[LExperiment].debug() << msg_discarded.str();
2859 
2861  precent_discarded > 1.0) {
2862  // The choosen threshold of 1% is a heuristical value
2863  logg[LExperiment].warn()
2864  << msg_discarded.str()
2865  << "\nThe number of discarded interactions is large, which means "
2866  "the assumption for the stochastic criterion of\n1 interaction "
2867  "per particle per timestep is probably violated. Consider "
2868  "reducing the timestep size.";
2869  }
2870 
2871  logg[LExperiment].info() << "Final interaction number: "
2873  }
2874 
2875  // Check if there are unformed particles
2876  int unformed_particles_count = 0;
2877  for (const Particles &particles : ensembles_) {
2878  for (const ParticleData &particle : particles) {
2879  if (particle.formation_time() > end_time_) {
2880  unformed_particles_count++;
2881  }
2882  }
2883  }
2884  if (unformed_particles_count > 0) {
2885  logg[LExperiment].warn(
2886  "End time might be too small. ", unformed_particles_count,
2887  " unformed particles were found at the end of the evolution.");
2888  }
2889  }
2890 
2891  // Keep track of how many ensembles had interactions
2893 
2894  for (const auto &output : outputs_) {
2895  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2896  auto event_info = fill_event_info(
2897  ensembles_, E_mean_field, modus_.impact_parameter(), parameters_,
2899  output->at_eventend(ensembles_[i_ens],
2900  // Pretend each ensemble is an independent event
2901  event_ * parameters_.n_ensembles + i_ens, event_info);
2902  }
2903  // For thermodynamic output
2904  output->at_eventend(ensembles_, event_);
2905 
2906  // For thermodynamic lattice output
2908  output->at_eventend(event_, ThermodynamicQuantity::EckartDensity,
2910  output->at_eventend(ThermodynamicQuantity::EckartDensity);
2911  }
2912  if (printout_tmn_) {
2913  output->at_eventend(event_, ThermodynamicQuantity::Tmn,
2915  output->at_eventend(ThermodynamicQuantity::Tmn);
2916  }
2917  if (printout_tmn_landau_) {
2918  output->at_eventend(event_, ThermodynamicQuantity::TmnLandau,
2920  output->at_eventend(ThermodynamicQuantity::TmnLandau);
2921  }
2922  if (printout_v_landau_) {
2923  output->at_eventend(event_, ThermodynamicQuantity::LandauVelocity,
2925  output->at_eventend(ThermodynamicQuantity::LandauVelocity);
2926  }
2927  if (printout_j_QBS_) {
2928  output->at_eventend(ThermodynamicQuantity::j_QBS);
2929  }
2930  }
2931 }
void count_nonempty_ensembles()
Counts the number of ensembles in wich interactions took place at the end of an event.
Definition: experiment.h:2934
double x0() const
Definition: fourvector.h:313
FourVector momentum() const
@ Stochastic
Stochastic Criteiron.
#define likely(x)
Tell the branch predictor that this expression is likely true.
Definition: macros.h:14
constexpr double really_small
Numerical error tolerance.
Definition: constants.h:37

◆ first_ensemble()

template<typename Modus >
Particles* smash::Experiment< Modus >::first_ensemble ( )
inline

Provides external access to SMASH particles.

This is helpful if SMASH is used as a 3rd-party library.

Definition at line 256 of file experiment.h.

256 { return &ensembles_[0]; }

◆ all_ensembles()

template<typename Modus >
std::vector<Particles>* smash::Experiment< Modus >::all_ensembles ( )
inline

Getter for all ensembles.

Definition at line 258 of file experiment.h.

258 { return &ensembles_; }

◆ modus()

template<typename Modus >
Modus* smash::Experiment< Modus >::modus ( )
inline

Provides external access to SMASH calculation modus.

This is helpful if SMASH is used as a 3rd-party library.

Definition at line 264 of file experiment.h.

264 { return &modus_; }

◆ increase_event_number()

template<typename Modus >
void smash::Experiment< Modus >::increase_event_number

Increases the event number by one.

This function is helpful if SMASH is used as a 3rd-party library.

Definition at line 2969 of file experiment.h.

2969  {
2970  event_++;
2971 }

◆ perform_action()

template<typename Modus >
bool smash::Experiment< Modus >::perform_action ( Action action,
int  i_ensemble,
bool  include_pauli_blocking = true 
)
private

Perform the given action.

Parameters
[in]actionThe action to perform
[in]i_ensembleindex of ensemble in which action is performed
[in]include_pauli_blockingwheter to take Pauli blocking into account. Skipping Pauli blocking is useful for example for final decays.
Returns
False if the action is rejected either due to invalidity or Pauli-blocking, or true if it's accepted and performed.

◆ create_output()

template<typename Modus >
void smash::Experiment< Modus >::create_output ( const std::string &  format,
const std::string &  content,
const std::filesystem::path &  output_path,
const OutputParameters par 
)
private

Create a list of output files.

Parameters
[in]formatFormat of the output file (e.g. Root, Oscar, Vtk)
[in]contentContent of the output (e.g. particles, collisions)
[in]output_pathPath of the output file
[in]parOutput options.(e.g. Extended)

Definition at line 724 of file experiment.h.

727  {
728  logg[LExperiment].info() << "Adding output " << content << " of format "
729  << format << std::endl;
730 
731  if (format == "VTK" && content == "Particles") {
732  outputs_.emplace_back(
733  std::make_unique<VtkOutput>(output_path, content, out_par));
734  } else if (format == "Root") {
735 #ifdef SMASH_USE_ROOT
736  if (content == "Initial_Conditions") {
737  outputs_.emplace_back(
738  std::make_unique<RootOutput>(output_path, "SMASH_IC", out_par));
739  } else {
740  outputs_.emplace_back(
741  std::make_unique<RootOutput>(output_path, content, out_par));
742  }
743 #else
744  logg[LExperiment].error(
745  "Root output requested, but Root support not compiled in");
746 #endif
747  } else if (format == "Binary") {
748  if (content == "Collisions" || content == "Dileptons" ||
749  content == "Photons") {
750  outputs_.emplace_back(std::make_unique<BinaryOutputCollisions>(
751  output_path, content, out_par));
752  } else if (content == "Particles") {
753  outputs_.emplace_back(std::make_unique<BinaryOutputParticles>(
754  output_path, content, out_par));
755  } else if (content == "Initial_Conditions") {
756  outputs_.emplace_back(std::make_unique<BinaryOutputInitialConditions>(
757  output_path, content, out_par));
758  }
759  } else if (format == "Oscar1999" || format == "Oscar2013") {
760  outputs_.emplace_back(
761  create_oscar_output(format, content, output_path, out_par));
762  } else if (content == "Thermodynamics" && format == "ASCII") {
763  outputs_.emplace_back(
764  std::make_unique<ThermodynamicOutput>(output_path, content, out_par));
765  } else if (content == "Thermodynamics" &&
766  (format == "Lattice_ASCII" || format == "Lattice_Binary")) {
768  if (format == "Lattice_ASCII") {
770  }
771  if (format == "Lattice_Binary") {
773  }
774  outputs_.emplace_back(std::make_unique<ThermodynamicLatticeOutput>(
775  output_path, content, out_par, printout_full_lattice_ascii_td_,
777  } else if (content == "Thermodynamics" && format == "VTK") {
778  printout_lattice_td_ = true;
779  outputs_.emplace_back(
780  std::make_unique<VtkOutput>(output_path, content, out_par));
781  } else if (content == "Initial_Conditions" && format == "ASCII") {
782  outputs_.emplace_back(
783  std::make_unique<ICOutput>(output_path, "SMASH_IC", out_par));
784  } else if ((format == "HepMC") || (format == "HepMC_asciiv3") ||
785  (format == "HepMC_treeroot")) {
786 #ifdef SMASH_USE_HEPMC
787  if (content == "Particles") {
788  if ((format == "HepMC") || (format == "HepMC_asciiv3")) {
789  outputs_.emplace_back(std::make_unique<HepMcOutput>(
790  output_path, "SMASH_HepMC_particles", false, "asciiv3"));
791  } else if (format == "HepMC_treeroot") {
792 #ifdef SMASH_USE_HEPMC_ROOTIO
793  outputs_.emplace_back(std::make_unique<HepMcOutput>(
794  output_path, "SMASH_HepMC_particles", false, "root"));
795 #else
796  logg[LExperiment].error(
797  "Requested HepMC_treeroot output not available, "
798  "ROOT or HepMC3 ROOTIO missing or not found by cmake.");
799 #endif
800  }
801  } else if (content == "Collisions") {
802  if ((format == "HepMC") || (format == "HepMC_asciiv3")) {
803  outputs_.emplace_back(std::make_unique<HepMcOutput>(
804  output_path, "SMASH_HepMC_collisions", true, "asciiv3"));
805  } else if (format == "HepMC_treeroot") {
806 #ifdef SMASH_USE_HEPMC_ROOTIO
807  outputs_.emplace_back(std::make_unique<HepMcOutput>(
808  output_path, "SMASH_HepMC_collisions", true, "root"));
809 #else
810  logg[LExperiment].error(
811  "Requested HepMC_treeroot output not available, "
812  "ROOT or HepMC3 ROOTIO missing or not found by cmake.");
813 #endif
814  }
815  } else {
816  logg[LExperiment].error(
817  "HepMC only available for Particles and "
818  "Collisions content. Requested for " +
819  content + ".");
820  }
821 #else
822  logg[LExperiment].error(
823  "HepMC output requested, but HepMC support not compiled in");
824 #endif
825  } else if (content == "Coulomb" && format == "VTK") {
826  outputs_.emplace_back(
827  std::make_unique<VtkOutput>(output_path, "Fields", out_par));
828  } else if (content == "Rivet") {
829 #ifdef SMASH_USE_RIVET
830  // flag to ensure that the Rivet format has not been already assigned
831  static bool rivet_format_already_selected = false;
832  // if the next check is true, then we are trying to assign the format twice
833  if (rivet_format_already_selected) {
834  logg[LExperiment].warn(
835  "Rivet output format can only be one, either YODA or YODA-full. "
836  "Only your first valid choice will be used.");
837  return;
838  }
839  if (format == "YODA") {
840  outputs_.emplace_back(std::make_unique<RivetOutput>(
841  output_path, "SMASH_Rivet", false, out_par.rivet_parameters));
842  rivet_format_already_selected = true;
843  } else if (format == "YODA-full") {
844  outputs_.emplace_back(std::make_unique<RivetOutput>(
845  output_path, "SMASH_Rivet_full", true, out_par.rivet_parameters));
846  rivet_format_already_selected = true;
847  } else {
848  logg[LExperiment].error("Rivet format " + format +
849  "not one of YODA or YODA-full");
850  }
851 #else
852  logg[LExperiment].error(
853  "Rivet output requested, but Rivet support not compiled in");
854 #endif
855  } else {
856  logg[LExperiment].error()
857  << "Unknown combination of format (" << format << ") and content ("
858  << content << "). Fix the config.";
859  }
860 }
bool printout_full_lattice_ascii_td_
Whether to print the thermodynamics quantities evaluated on the lattices, point by point,...
Definition: experiment.h:520
bool printout_lattice_td_
Whether to print the thermodynamics quantities evaluated on the lattices.
Definition: experiment.h:516
bool printout_full_lattice_binary_td_
Whether to print the thermodynamics quantities evaluated on the lattices, point by point,...
Definition: experiment.h:524
FormattingHelper< T > format(const T &value, const char *unit, int width=-1, int precision=-1)
Acts as a stream modifier for std::ostream to output an object with an optional suffix string and wit...
Definition: logging.h:214
std::unique_ptr< OutputInterface > create_oscar_output(const std::string &format, const std::string &content, const std::filesystem::path &path, const OutputParameters &out_par)
Definition: oscaroutput.cc:810

◆ propagate_and_shine()

template<typename Modus >
void smash::Experiment< Modus >::propagate_and_shine ( double  to_time,
Particles particles 
)
private

Propagate all particles until time to_time without any interactions and shine dileptons.

Parameters
[in]to_timeTime at the end of propagation [fm]
[in,out]particlesParticles to be propagated

Definition at line 2397 of file experiment.h.

2398  {
2399  const double dt =
2400  propagate_straight_line(&particles, to_time, beam_momentum_);
2401  if (dilepton_finder_ != nullptr) {
2402  for (const auto &output : outputs_) {
2403  dilepton_finder_->shine(particles, output.get(), dt);
2404  }
2405  }
2406 }
double propagate_straight_line(Particles *particles, double to_time, const std::vector< FourVector > &beam_momentum)
Propagates the positions of all particles on a straight line to a given moment.
Definition: propagation.cc:44

◆ run_time_evolution_timestepless()

template<typename Modus >
void smash::Experiment< Modus >::run_time_evolution_timestepless ( Actions actions,
int  i_ensemble,
const double  end_time_propagation 
)
private

Performs all the propagations and actions during a certain time interval neglecting the influence of the potentials.

This function is called in either the time stepless cases or the cases with time steps.

Parameters
[in,out]actionsActions occur during a certain time interval. They provide the ending times of the propagations and are updated during the time interval.
[in]i_ensembleindex of ensemble to be evolved
[in]end_time_propagationtime until propagation should be performed

Definition at line 2423 of file experiment.h.

2424  {
2425  Particles &particles = ensembles_[i_ensemble];
2426  logg[LExperiment].debug(
2427  "Timestepless propagation: ", "Actions size = ", actions.size(),
2428  ", end time = ", end_time_propagation);
2429 
2430  // iterate over all actions
2431  while (!actions.is_empty()) {
2432  if (actions.earliest_time() > end_time_propagation) {
2433  break;
2434  }
2435  // get next action
2436  ActionPtr act = actions.pop();
2437  if (!act->is_valid(particles)) {
2439  logg[LExperiment].debug(~einhard::DRed(), "✘ ", act,
2440  " (discarded: invalid)");
2441  continue;
2442  }
2443  logg[LExperiment].debug(~einhard::Green(), "✔ ", act,
2444  ", action time = ", act->time_of_execution());
2445 
2446  /* (1) Propagate to the next action. */
2447  propagate_and_shine(act->time_of_execution(), particles);
2448 
2449  /* (2) Perform action.
2450  *
2451  * Update the positions of the incoming particles, because the information
2452  * in the action object will be outdated as the particles have been
2453  * propagated since the construction of the action. */
2454  act->update_incoming(particles);
2455  const bool performed = perform_action(*act, i_ensemble);
2456 
2457  /* No need to update actions for outgoing particles
2458  * if the action is not performed. */
2459  if (!performed) {
2460  continue;
2461  }
2462 
2463  /* (3) Update actions for newly-produced particles. */
2464 
2465  const double end_time_timestep = parameters_.labclock->next_time();
2466  // New actions are always search until the end of the current timestep
2467  const double time_left = end_time_timestep - act->time_of_execution();
2468  const ParticleList &outgoing_particles = act->outgoing_particles();
2469  // Grid cell volume set to zero, since there is no grid
2470  const double gcell_vol = 0.0;
2471  for (const auto &finder : action_finders_) {
2472  // Outgoing particles can still decay, cross walls...
2473  actions.insert(finder->find_actions_in_cell(outgoing_particles, time_left,
2474  gcell_vol, beam_momentum_));
2475  // ... and collide with other particles.
2476  actions.insert(finder->find_actions_with_surrounding_particles(
2477  outgoing_particles, particles, time_left, beam_momentum_));
2478  }
2479 
2481  }
2482 
2483  propagate_and_shine(end_time_propagation, particles);
2484 }
A stream modifier that allows to colorize the log output.
Definition: einhard.hpp:143
void propagate_and_shine(double to_time, Particles &particles)
Propagate all particles until time to_time without any interactions and shine dileptons.
Definition: experiment.h:2397
void check_interactions_total(uint64_t interactions_total)
Make sure interactions_total can be represented as a 32-bit integer.
Definition: experiment.h:2415

◆ intermediate_output()

template<typename Modus >
void smash::Experiment< Modus >::intermediate_output
private

Intermediate output during an event.

Auxiliary variable to communicate the time in the computational frame at the functions printing the thermodynamics lattice output

Definition at line 2487 of file experiment.h.

2487  {
2488  const uint64_t wall_actions_this_interval =
2491  const uint64_t interactions_this_interval = interactions_total_ -
2493  wall_actions_this_interval;
2495  double E_mean_field = 0.0;
2498  double computational_frame_time = 0.0;
2499  if (potentials_) {
2500  // using the lattice is necessary
2501  if ((jmu_B_lat_ != nullptr)) {
2503  EM_lat_.get(), parameters_);
2504  /*
2505  * Mean field calculated in a box should remain approximately constant if
2506  * the system is in equilibrium, and so deviations from its original value
2507  * may signal a phase transition or other dynamical process. This
2508  * comparison only makes sense in the Box Modus, hence the condition.
2509  */
2510  if (modus_.is_box()) {
2511  double tmp = (E_mean_field - initial_mean_field_energy_) /
2512  (E_mean_field + initial_mean_field_energy_);
2513  /*
2514  * This is displayed when the system evolves away from its initial
2515  * configuration (which is when the total mean field energy in the box
2516  * deviates from its initial value).
2517  */
2518  if (std::abs(tmp) > 0.01) {
2519  logg[LExperiment].info()
2520  << "\n\n\n\t The mean field at t = "
2521  << parameters_.outputclock->current_time()
2522  << " [fm] differs from the mean field at t = 0:"
2523  << "\n\t\t initial_mean_field_energy_ = "
2524  << initial_mean_field_energy_ << " [GeV]"
2525  << "\n\t\t abs[(E_MF - E_MF(t=0))/(E_MF + E_MF(t=0))] = "
2526  << std::abs(tmp)
2527  << "\n\t\t E_MF/E_MF(t=0) = "
2528  << E_mean_field / initial_mean_field_energy_ << "\n\n";
2529  }
2530  }
2531  }
2532  }
2533 
2535  ensembles_, interactions_this_interval, conserved_initial_, time_start_,
2536  parameters_.outputclock->current_time(), E_mean_field,
2538  const LatticeUpdate lat_upd = LatticeUpdate::AtOutput;
2539 
2540  // save evolution data
2541  if (!(modus_.is_box() && parameters_.outputclock->current_time() <
2542  modus_.equilibration_time())) {
2543  for (const auto &output : outputs_) {
2544  if (output->is_dilepton_output() || output->is_photon_output() ||
2545  output->is_IC_output()) {
2546  continue;
2547  }
2548  for (int i_ens = 0; i_ens < parameters_.n_ensembles; i_ens++) {
2549  auto event_info = fill_event_info(
2550  ensembles_, E_mean_field, modus_.impact_parameter(), parameters_,
2552 
2553  output->at_intermediate_time(ensembles_[i_ens], parameters_.outputclock,
2554  density_param_, event_info);
2555  computational_frame_time = event_info.current_time;
2556  }
2557  // For thermodynamic output
2558  output->at_intermediate_time(ensembles_, parameters_.outputclock,
2559  density_param_);
2560 
2561  // Thermodynamic output on the lattice versus time
2562  switch (dens_type_lattice_printout_) {
2563  case DensityType::Baryon:
2565  density_param_, ensembles_, false);
2566  output->thermodynamics_output(ThermodynamicQuantity::EckartDensity,
2568  output->thermodynamics_lattice_output(*jmu_B_lat_,
2569  computational_frame_time);
2570  break;
2572  update_lattice(jmu_I3_lat_.get(), lat_upd,
2574  ensembles_, false);
2575  output->thermodynamics_output(ThermodynamicQuantity::EckartDensity,
2577  *jmu_I3_lat_);
2578  output->thermodynamics_lattice_output(*jmu_I3_lat_,
2579  computational_frame_time);
2580  break;
2581  case DensityType::None:
2582  break;
2583  default:
2584  update_lattice(jmu_custom_lat_.get(), lat_upd,
2586  ensembles_, false);
2587  output->thermodynamics_output(ThermodynamicQuantity::EckartDensity,
2589  *jmu_custom_lat_);
2590  output->thermodynamics_lattice_output(*jmu_custom_lat_,
2591  computational_frame_time);
2592  }
2595  density_param_, ensembles_, false);
2596  if (printout_tmn_) {
2597  output->thermodynamics_output(ThermodynamicQuantity::Tmn,
2599  output->thermodynamics_lattice_output(
2600  ThermodynamicQuantity::Tmn, *Tmn_, computational_frame_time);
2601  }
2602  if (printout_tmn_landau_) {
2603  output->thermodynamics_output(ThermodynamicQuantity::TmnLandau,
2605  output->thermodynamics_lattice_output(
2607  computational_frame_time);
2608  }
2609  if (printout_v_landau_) {
2610  output->thermodynamics_output(ThermodynamicQuantity::LandauVelocity,
2612  output->thermodynamics_lattice_output(
2614  computational_frame_time);
2615  }
2616  }
2617  if (EM_lat_) {
2618  output->fields_output("Efield", "Bfield", *EM_lat_);
2619  }
2620  if (printout_j_QBS_) {
2621  output->thermodynamics_lattice_output(
2622  *j_QBS_lat_, computational_frame_time, ensembles_, density_param_);
2623  }
2624 
2625  if (thermalizer_) {
2626  output->thermodynamics_output(*thermalizer_);
2627  }
2628  }
2629  }
2630 }
LatticeUpdate
Enumerator option for lattice updates.
Definition: lattice.h:36

◆ update_potentials()

template<typename Modus >
void smash::Experiment< Modus >::update_potentials
private

Recompute potentials on lattices if necessary.

Definition at line 2633 of file experiment.h.

2633  {
2634  if (potentials_) {
2635  if (potentials_->use_symmetry() && jmu_I3_lat_ != nullptr) {
2640  parameters_.labclock->timestep_duration(), true);
2641  }
2642  if ((potentials_->use_skyrme() || potentials_->use_symmetry()) &&
2643  jmu_B_lat_ != nullptr) {
2648  parameters_.labclock->timestep_duration(), true);
2649  const size_t UBlattice_size = UB_lat_->size();
2650  for (size_t i = 0; i < UBlattice_size; i++) {
2651  auto jB = (*jmu_B_lat_)[i];
2652  const FourVector flow_four_velocity_B =
2653  std::abs(jB.rho()) > very_small_double ? jB.jmu_net() / jB.rho()
2654  : FourVector();
2655  double baryon_density = jB.rho();
2656  ThreeVector baryon_grad_j0 = jB.grad_j0();
2657  ThreeVector baryon_dvecj_dt = jB.dvecj_dt();
2658  ThreeVector baryon_curl_vecj = jB.curl_vecj();
2659  if (potentials_->use_skyrme()) {
2660  (*UB_lat_)[i] =
2661  flow_four_velocity_B * potentials_->skyrme_pot(baryon_density);
2662  (*FB_lat_)[i] =
2663  potentials_->skyrme_force(baryon_density, baryon_grad_j0,
2664  baryon_dvecj_dt, baryon_curl_vecj);
2665  }
2666  if (potentials_->use_symmetry() && jmu_I3_lat_ != nullptr) {
2667  auto jI3 = (*jmu_I3_lat_)[i];
2668  const FourVector flow_four_velocity_I3 =
2669  std::abs(jI3.rho()) > very_small_double
2670  ? jI3.jmu_net() / jI3.rho()
2671  : FourVector();
2672  (*UI3_lat_)[i] = flow_four_velocity_I3 *
2673  potentials_->symmetry_pot(jI3.rho(), baryon_density);
2674  (*FI3_lat_)[i] = potentials_->symmetry_force(
2675  jI3.rho(), jI3.grad_j0(), jI3.dvecj_dt(), jI3.curl_vecj(),
2676  baryon_density, baryon_grad_j0, baryon_dvecj_dt,
2677  baryon_curl_vecj);
2678  }
2679  }
2680  }
2681  if (potentials_->use_coulomb()) {
2684  for (size_t i = 0; i < EM_lat_->size(); i++) {
2685  ThreeVector electric_field = {0., 0., 0.};
2686  ThreeVector position = jmu_el_lat_->cell_center(i);
2687  jmu_el_lat_->integrate_volume(electric_field,
2689  potentials_->coulomb_r_cut(), position);
2690  ThreeVector magnetic_field = {0., 0., 0.};
2691  jmu_el_lat_->integrate_volume(magnetic_field,
2693  potentials_->coulomb_r_cut(), position);
2694  (*EM_lat_)[i] = std::make_pair(electric_field, magnetic_field);
2695  }
2696  } // if ((potentials_->use_skyrme() || ...
2697  if (potentials_->use_vdf() && jmu_B_lat_ != nullptr) {
2702  parameters_.labclock->timestep_duration(), true);
2705  fields_lat_.get(), old_fields_auxiliary_.get(),
2708  parameters_.labclock->timestep_duration());
2709  }
2710  const size_t UBlattice_size = UB_lat_->size();
2711  for (size_t i = 0; i < UBlattice_size; i++) {
2712  auto jB = (*jmu_B_lat_)[i];
2713  (*UB_lat_)[i] = potentials_->vdf_pot(jB.rho(), jB.jmu_net());
2716  (*FB_lat_)[i] = potentials_->vdf_force(
2717  jB.rho(), jB.drho_dxnu().x0(), jB.drho_dxnu().threevec(),
2718  jB.grad_rho_cross_vecj(), jB.jmu_net().x0(), jB.grad_j0(),
2719  jB.jmu_net().threevec(), jB.dvecj_dt(), jB.curl_vecj());
2720  break;
2722  auto Amu = (*fields_lat_)[i];
2723  (*FB_lat_)[i] = potentials_->vdf_force(
2724  Amu.grad_A0(), Amu.dvecA_dt(), Amu.curl_vecA());
2725  break;
2726  }
2727  } // for (size_t i = 0; i < UBlattice_size; i++)
2728  } // if potentials_->use_vdf()
2729  }
2730 }
std::unique_ptr< RectangularLattice< FourVector > > new_fields_auxiliary_
Auxiliary lattice for values of Amu at a time step t0 + dt.
Definition: experiment.h:498
std::unique_ptr< RectangularLattice< std::array< FourVector, 4 > > > fields_four_gradient_auxiliary_
Auxiliary lattice for calculating the four-gradient of Amu.
Definition: experiment.h:501
std::unique_ptr< FieldsLattice > fields_lat_
Mean-field A^mu on the lattice.
Definition: experiment.h:443
std::unique_ptr< RectangularLattice< FourVector > > UB_lat_
Lattices for Skyrme or VDF potentials (evaluated in the local rest frame) times the baryon flow 4-vel...
Definition: experiment.h:461
std::unique_ptr< RectangularLattice< FourVector > > old_fields_auxiliary_
Auxiliary lattice for values of Amu at a time step t0.
Definition: experiment.h:496
std::unique_ptr< DensityLattice > jmu_el_lat_
Electric charge density on the lattice.
Definition: experiment.h:440
static ThreeVector B_field_integrand(ThreeVector pos, DensityOnLattice &charge_density, ThreeVector point)
Integrand for calculating the magnetic field using the Biot-Savart formula.
Definition: potentials.h:235
static ThreeVector E_field_integrand(ThreeVector pos, DensityOnLattice &charge_density, ThreeVector point)
Integrand for calculating the electric field.
Definition: potentials.h:217
constexpr double very_small_double
A very small double, used to avoid division by zero.
Definition: constants.h:40
void update_fields_lattice(RectangularLattice< FieldsOnLattice > *fields_lat, RectangularLattice< FourVector > *old_fields, RectangularLattice< FourVector > *new_fields, RectangularLattice< std::array< FourVector, 4 >> *fields_four_grad_lattice, DensityLattice *jmu_B_lat, const LatticeUpdate fields_lat_update, const Potentials &potentials, const double time_step)
Updates the contents on the lattice of FieldsOnLattice type.
Definition: fields.cc:14
FieldDerivativesMode field_derivatives_mode
mode of calculating field derivatives

◆ compute_min_cell_length()

template<typename Modus >
double smash::Experiment< Modus >::compute_min_cell_length ( double  dt) const
inlineprivate

Calculate the minimal size for the grid cells such that the ScatterActionsFinder will find all collisions within the maximal transverse distance (which is determined by the maximal cross section).

Parameters
[in]dtThe current time step size [fm]
Returns
The minimal required size of cells

Definition at line 337 of file experiment.h.

337  {
340  }
341  return std::sqrt(4 * dt * dt + max_transverse_distance_sqr_);
342  }
double fixed_min_cell_length
Fixed minimal grid cell length (in fm).

◆ next_output_time()

template<typename Modus >
double smash::Experiment< Modus >::next_output_time ( ) const
inlineprivate

Shortcut for next output time.

Definition at line 345 of file experiment.h.

345  {
346  return parameters_.outputclock->next_time();
347  }

◆ count_nonempty_ensembles()

template<typename Modus >
void smash::Experiment< Modus >::count_nonempty_ensembles
private

Counts the number of ensembles in wich interactions took place at the end of an event.

Definition at line 2934 of file experiment.h.

2934  {
2935  for (bool has_interaction : projectile_target_interact_) {
2936  if (has_interaction) {
2938  }
2939  }
2940 }
int nonempty_ensembles_
Number of ensembles containing an interaction.
Definition: experiment.h:577

◆ is_finished()

template<typename Modus >
bool smash::Experiment< Modus >::is_finished
private

Checks wether the desired number events have been calculated.

Returns
wether the experiment is is_finished

Definition at line 2943 of file experiment.h.

2943  {
2945  return event_ >= nevents_;
2946  }
2949  return true;
2950  }
2951  if (event_ >= max_events_) {
2952  logg[LExperiment].warn()
2953  << "Maximum number of events (" << max_events_
2954  << ") exceeded. Stopping calculation. "
2955  << "The fraction of empty ensembles is "
2956  << (1.0 - static_cast<double>(nonempty_ensembles_) /
2958  << ". If this fraction is expected, try increasing the "
2959  "Maximum_Ensembles_Run.";
2960  return true;
2961  }
2962  return false;
2963  }
2964  throw std::runtime_error("Event counting option is invalid");
2965  return false;
2966 }
int minimum_nonempty_ensembles_
The number of ensembles, in which interactions take place, to be calculated.
Definition: experiment.h:563
const int nevents_
Number of events.
Definition: experiment.h:553
int max_events_
Maximum number of events to be calculated in order obtain the desired number of non-empty events usin...
Definition: experiment.h:583
EventCounting event_counting_
The way in which the number of calculated events is specified.
Definition: experiment.h:571

Friends And Related Function Documentation

◆ ExperimentBase

template<typename Modus >
friend class ExperimentBase
friend

Definition at line 188 of file experiment.h.

Member Data Documentation

◆ parameters_

template<typename Modus >
ExperimentParameters smash::Experiment< Modus >::parameters_
private

Struct of several member variables.

These variables are combined into a struct for efficient input to functions outside of this class.

Definition at line 367 of file experiment.h.

◆ density_param_

template<typename Modus >
DensityParameters smash::Experiment< Modus >::density_param_
private

Structure to precalculate and hold parameters for density computations.

Definition at line 370 of file experiment.h.

◆ modus_

template<typename Modus >
Modus smash::Experiment< Modus >::modus_
private

Instance of the Modus template parameter.

May store modus-specific data and contains modus-specific function implementations.

Definition at line 376 of file experiment.h.

◆ ensembles_

template<typename Modus >
std::vector<Particles> smash::Experiment< Modus >::ensembles_
private

Complete particle list, all ensembles in one vector.

Definition at line 379 of file experiment.h.

◆ potentials_

template<typename Modus >
std::unique_ptr<Potentials> smash::Experiment< Modus >::potentials_
private

An instance of potentials class, that stores parameters of potentials, calculates them and their gradients.

Definition at line 385 of file experiment.h.

◆ pauli_blocker_

template<typename Modus >
std::unique_ptr<PauliBlocker> smash::Experiment< Modus >::pauli_blocker_
private

An instance of PauliBlocker class that stores parameters needed for Pauli blocking calculations and computes phase-space density.

Definition at line 391 of file experiment.h.

◆ outputs_

template<typename Modus >
OutputsList smash::Experiment< Modus >::outputs_
private

A list of output formaters.

They will be called to write the state of the particles to file.

Definition at line 397 of file experiment.h.

◆ dilepton_output_

template<typename Modus >
OutputPtr smash::Experiment< Modus >::dilepton_output_
private

The Dilepton output.

Definition at line 400 of file experiment.h.

◆ photon_output_

template<typename Modus >
OutputPtr smash::Experiment< Modus >::photon_output_
private

The Photon output.

Definition at line 403 of file experiment.h.

◆ projectile_target_interact_

template<typename Modus >
std::vector<bool> smash::Experiment< Modus >::projectile_target_interact_
private

Whether the projectile and the target collided.

One value for each ensemble.

Definition at line 409 of file experiment.h.

◆ beam_momentum_

template<typename Modus >
std::vector<FourVector> smash::Experiment< Modus >::beam_momentum_ = {}
private

The initial nucleons in the ColliderModus propagate with beam_momentum_, if Fermi motion is frozen.

It's only valid in the ColliderModus, so is set as an empty vector by default.

Definition at line 416 of file experiment.h.

◆ action_finders_

template<typename Modus >
std::vector<std::unique_ptr<ActionFinderInterface> > smash::Experiment< Modus >::action_finders_
private

The Action finder objects.

Definition at line 419 of file experiment.h.

◆ dilepton_finder_

template<typename Modus >
std::unique_ptr<DecayActionsFinderDilepton> smash::Experiment< Modus >::dilepton_finder_
private

The Dilepton Action Finder.

Definition at line 422 of file experiment.h.

◆ photon_finder_

template<typename Modus >
std::unique_ptr<ActionFinderInterface> smash::Experiment< Modus >::photon_finder_
private

The (Scatter) Actions Finder for Direct Photons.

Definition at line 425 of file experiment.h.

◆ n_fractional_photons_

template<typename Modus >
int smash::Experiment< Modus >::n_fractional_photons_
private

Number of fractional photons produced per single reaction.

Definition at line 428 of file experiment.h.

◆ j_QBS_lat_

template<typename Modus >
std::unique_ptr<DensityLattice> smash::Experiment< Modus >::j_QBS_lat_
private

4-current for j_QBS lattice output

Definition at line 431 of file experiment.h.

◆ jmu_B_lat_

template<typename Modus >
std::unique_ptr<DensityLattice> smash::Experiment< Modus >::jmu_B_lat_
private

Baryon density on the lattice.

Definition at line 434 of file experiment.h.

◆ jmu_I3_lat_

template<typename Modus >
std::unique_ptr<DensityLattice> smash::Experiment< Modus >::jmu_I3_lat_
private

Isospin projection density on the lattice.

Definition at line 437 of file experiment.h.

◆ jmu_el_lat_

template<typename Modus >
std::unique_ptr<DensityLattice> smash::Experiment< Modus >::jmu_el_lat_
private

Electric charge density on the lattice.

Definition at line 440 of file experiment.h.

◆ fields_lat_

template<typename Modus >
std::unique_ptr<FieldsLattice> smash::Experiment< Modus >::fields_lat_
private

Mean-field A^mu on the lattice.

Definition at line 443 of file experiment.h.

◆ jmu_custom_lat_

template<typename Modus >
std::unique_ptr<DensityLattice> smash::Experiment< Modus >::jmu_custom_lat_
private

Custom density on the lattices.

In the config user asks for some kind of density for printout. Baryon and isospin projection density are anyway needed for potentials. If user asks for some other density type for printout, it will be handled using jmu_custom variable.

Definition at line 452 of file experiment.h.

◆ dens_type_lattice_printout_

template<typename Modus >
DensityType smash::Experiment< Modus >::dens_type_lattice_printout_ = DensityType::None
private

Type of density for lattice printout.

Definition at line 455 of file experiment.h.

◆ UB_lat_

template<typename Modus >
std::unique_ptr<RectangularLattice<FourVector> > smash::Experiment< Modus >::UB_lat_ = nullptr
private

Lattices for Skyrme or VDF potentials (evaluated in the local rest frame) times the baryon flow 4-velocity.

Definition at line 461 of file experiment.h.

◆ UI3_lat_

template<typename Modus >
std::unique_ptr<RectangularLattice<FourVector> > smash::Experiment< Modus >::UI3_lat_ = nullptr
private

Lattices for symmetry potentials (evaluated in the local rest frame) times the isospin flow 4-velocity.

Definition at line 467 of file experiment.h.

◆ FB_lat_

template<typename Modus >
std::unique_ptr<RectangularLattice<std::pair<ThreeVector, ThreeVector> > > smash::Experiment< Modus >::FB_lat_
private

Lattices for the electric and magnetic components of the Skyrme or VDF force.

Definition at line 474 of file experiment.h.

◆ FI3_lat_

template<typename Modus >
std::unique_ptr<RectangularLattice<std::pair<ThreeVector, ThreeVector> > > smash::Experiment< Modus >::FI3_lat_
private

Lattices for the electric and magnetic component of the symmetry force.

Definition at line 478 of file experiment.h.

◆ EM_lat_

template<typename Modus >
std::unique_ptr<RectangularLattice<std::pair<ThreeVector, ThreeVector> > > smash::Experiment< Modus >::EM_lat_
private

Lattices for electric and magnetic field in fm^-2.

Definition at line 482 of file experiment.h.

◆ Tmn_

template<typename Modus >
std::unique_ptr<RectangularLattice<EnergyMomentumTensor> > smash::Experiment< Modus >::Tmn_
private

Lattices of energy-momentum tensors for printout.

Definition at line 485 of file experiment.h.

◆ old_jmu_auxiliary_

template<typename Modus >
std::unique_ptr<RectangularLattice<FourVector> > smash::Experiment< Modus >::old_jmu_auxiliary_
private

Auxiliary lattice for values of jmu at a time step t0.

Definition at line 488 of file experiment.h.

◆ new_jmu_auxiliary_

template<typename Modus >
std::unique_ptr<RectangularLattice<FourVector> > smash::Experiment< Modus >::new_jmu_auxiliary_
private

Auxiliary lattice for values of jmu at a time step t0 + dt.

Definition at line 490 of file experiment.h.

◆ four_gradient_auxiliary_

template<typename Modus >
std::unique_ptr<RectangularLattice<std::array<FourVector, 4> > > smash::Experiment< Modus >::four_gradient_auxiliary_
private

Auxiliary lattice for calculating the four-gradient of jmu.

Definition at line 493 of file experiment.h.

◆ old_fields_auxiliary_

template<typename Modus >
std::unique_ptr<RectangularLattice<FourVector> > smash::Experiment< Modus >::old_fields_auxiliary_
private

Auxiliary lattice for values of Amu at a time step t0.

Definition at line 496 of file experiment.h.

◆ new_fields_auxiliary_

template<typename Modus >
std::unique_ptr<RectangularLattice<FourVector> > smash::Experiment< Modus >::new_fields_auxiliary_
private

Auxiliary lattice for values of Amu at a time step t0 + dt.

Definition at line 498 of file experiment.h.

◆ fields_four_gradient_auxiliary_

template<typename Modus >
std::unique_ptr<RectangularLattice<std::array<FourVector, 4> > > smash::Experiment< Modus >::fields_four_gradient_auxiliary_
private

Auxiliary lattice for calculating the four-gradient of Amu.

Definition at line 501 of file experiment.h.

◆ printout_tmn_

template<typename Modus >
bool smash::Experiment< Modus >::printout_tmn_ = false
private

Whether to print the energy-momentum tensor.

Definition at line 504 of file experiment.h.

◆ printout_tmn_landau_

template<typename Modus >
bool smash::Experiment< Modus >::printout_tmn_landau_ = false
private

Whether to print the energy-momentum tensor in Landau frame.

Definition at line 507 of file experiment.h.

◆ printout_v_landau_

template<typename Modus >
bool smash::Experiment< Modus >::printout_v_landau_ = false
private

Whether to print the 4-velocity in Landau frame.

Definition at line 510 of file experiment.h.

◆ printout_j_QBS_

template<typename Modus >
bool smash::Experiment< Modus >::printout_j_QBS_ = false
private

Whether to print the Q, B, S 4-currents.

Definition at line 513 of file experiment.h.

◆ printout_lattice_td_

template<typename Modus >
bool smash::Experiment< Modus >::printout_lattice_td_ = false
private

Whether to print the thermodynamics quantities evaluated on the lattices.

Definition at line 516 of file experiment.h.

◆ printout_full_lattice_ascii_td_

template<typename Modus >
bool smash::Experiment< Modus >::printout_full_lattice_ascii_td_ = false
private

Whether to print the thermodynamics quantities evaluated on the lattices, point by point, in ASCII format.

Definition at line 520 of file experiment.h.

◆ printout_full_lattice_binary_td_

template<typename Modus >
bool smash::Experiment< Modus >::printout_full_lattice_binary_td_ = false
private

Whether to print the thermodynamics quantities evaluated on the lattices, point by point, in Binary format.

Definition at line 524 of file experiment.h.

◆ printout_full_lattice_any_td_

template<typename Modus >
bool smash::Experiment< Modus >::printout_full_lattice_any_td_ = false
private

Whether to print the thermodynamics quantities evaluated on the lattices, point by point, in any format.

Definition at line 528 of file experiment.h.

◆ thermalizer_

template<typename Modus >
std::unique_ptr<GrandCanThermalizer> smash::Experiment< Modus >::thermalizer_
private

Instance of class used for forced thermalization.

Definition at line 531 of file experiment.h.

◆ process_string_ptr_

template<typename Modus >
StringProcess* smash::Experiment< Modus >::process_string_ptr_
private

Pointer to the string process class object, which is used to set the random seed for PYTHIA objects in each event.

Definition at line 537 of file experiment.h.

◆ nevents_

template<typename Modus >
const int smash::Experiment< Modus >::nevents_ = 0
private

Number of events.

Event is a single simulation of a physical phenomenon: elementary particle or nucleus-nucleus collision. Result of a single SMASH event is random (by construction) as well as result of one collision in nature. To compare simulation with experiment one has to take ensemble averages, i.e. perform simulation and real experiment many times and compare average results.

nevents_ is number of times single phenomenon (particle or nucleus-nucleus collision) will be simulated.

Definition at line 553 of file experiment.h.

◆ minimum_nonempty_ensembles_

template<typename Modus >
int smash::Experiment< Modus >::minimum_nonempty_ensembles_ = 0
private

The number of ensembles, in which interactions take place, to be calculated.

Can be specified as an inout instead of the number of events. In this case events will be calculated until this number of ensembles is reached.

Definition at line 563 of file experiment.h.

◆ event_counting_

template<typename Modus >
EventCounting smash::Experiment< Modus >::event_counting_ = EventCounting::Invalid
private

The way in which the number of calculated events is specified.

Can be either a fixed number of simulated events or a minimum number of events that contain interactions.

Definition at line 571 of file experiment.h.

◆ event_

template<typename Modus >
int smash::Experiment< Modus >::event_ = 0
private

Current event.

Definition at line 574 of file experiment.h.

◆ nonempty_ensembles_

template<typename Modus >
int smash::Experiment< Modus >::nonempty_ensembles_ = 0
private

Number of ensembles containing an interaction.

Definition at line 577 of file experiment.h.

◆ max_events_

template<typename Modus >
int smash::Experiment< Modus >::max_events_ = 0
private

Maximum number of events to be calculated in order obtain the desired number of non-empty events using the MinimumNonemptyEnsembles option.

Definition at line 583 of file experiment.h.

◆ end_time_

template<typename Modus >
const double smash::Experiment< Modus >::end_time_
private

simulation time at which the evolution is stopped.

Definition at line 586 of file experiment.h.

◆ delta_time_startup_

template<typename Modus >
const double smash::Experiment< Modus >::delta_time_startup_
private

The clock's timestep size at start up.

Stored here so that the next event will remember this.

Definition at line 593 of file experiment.h.

◆ force_decays_

template<typename Modus >
const bool smash::Experiment< Modus >::force_decays_
private

This indicates whether we force all resonances to decay in the last timestep.

Definition at line 599 of file experiment.h.

◆ use_grid_

template<typename Modus >
const bool smash::Experiment< Modus >::use_grid_
private

This indicates whether to use the grid.

Definition at line 602 of file experiment.h.

◆ metric_

template<typename Modus >
const ExpansionProperties smash::Experiment< Modus >::metric_
private

This struct contains information on the metric to be used.

Definition at line 605 of file experiment.h.

◆ dileptons_switch_

template<typename Modus >
const bool smash::Experiment< Modus >::dileptons_switch_
private

This indicates whether dileptons are switched on.

Definition at line 608 of file experiment.h.

◆ photons_switch_

template<typename Modus >
const bool smash::Experiment< Modus >::photons_switch_
private

This indicates whether photons are switched on.

Definition at line 611 of file experiment.h.

◆ bremsstrahlung_switch_

template<typename Modus >
const bool smash::Experiment< Modus >::bremsstrahlung_switch_
private

This indicates whether bremsstrahlung is switched on.

Definition at line 614 of file experiment.h.

◆ IC_output_switch_

template<typename Modus >
const bool smash::Experiment< Modus >::IC_output_switch_
private

This indicates whether the IC output is enabled.

Definition at line 617 of file experiment.h.

◆ time_step_mode_

template<typename Modus >
const TimeStepMode smash::Experiment< Modus >::time_step_mode_
private

This indicates whether to use time steps.

Definition at line 620 of file experiment.h.

◆ max_transverse_distance_sqr_

template<typename Modus >
double smash::Experiment< Modus >::max_transverse_distance_sqr_ = std::numeric_limits<double>::max()
private

Maximal distance at which particles can interact in case of the geometric criterion, squared.

Definition at line 626 of file experiment.h.

◆ conserved_initial_

template<typename Modus >
QuantumNumbers smash::Experiment< Modus >::conserved_initial_
private

The conserved quantities of the system.

This struct carries the sums of the single particle's various quantities as measured at the beginning of the evolution and can be used to regularly check if they are still good.

Definition at line 635 of file experiment.h.

◆ initial_mean_field_energy_

template<typename Modus >
double smash::Experiment< Modus >::initial_mean_field_energy_
private

The initial total mean field energy in the system.

Note: will only be calculated if lattice is on.

Definition at line 641 of file experiment.h.

◆ time_start_

template<typename Modus >
SystemTimePoint smash::Experiment< Modus >::time_start_ = SystemClock::now()
private

system starting time of the simulation

Definition at line 644 of file experiment.h.

◆ dens_type_

template<typename Modus >
DensityType smash::Experiment< Modus >::dens_type_ = DensityType::None
private

Type of density to be written to collision headers.

Definition at line 647 of file experiment.h.

◆ interactions_total_

template<typename Modus >
uint64_t smash::Experiment< Modus >::interactions_total_ = 0
private

Total number of interactions for current timestep.

For timestepless mode the whole run time is considered as one timestep.

Definition at line 653 of file experiment.h.

◆ previous_interactions_total_

template<typename Modus >
uint64_t smash::Experiment< Modus >::previous_interactions_total_ = 0
private

Total number of interactions for previous timestep.

For timestepless mode the whole run time is considered as one timestep.

Definition at line 659 of file experiment.h.

◆ wall_actions_total_

template<typename Modus >
uint64_t smash::Experiment< Modus >::wall_actions_total_ = 0
private

Total number of wall-crossings for current timestep.

For timestepless mode the whole run time is considered as one timestep.

Definition at line 665 of file experiment.h.

◆ previous_wall_actions_total_

template<typename Modus >
uint64_t smash::Experiment< Modus >::previous_wall_actions_total_ = 0
private

Total number of wall-crossings for previous timestep.

For timestepless mode the whole run time is considered as one timestep.

Definition at line 671 of file experiment.h.

◆ total_pauli_blocked_

template<typename Modus >
uint64_t smash::Experiment< Modus >::total_pauli_blocked_ = 0
private

Total number of Pauli-blockings for current timestep.

For timestepless mode the whole run time is considered as one timestep.

Definition at line 677 of file experiment.h.

◆ total_hypersurface_crossing_actions_

template<typename Modus >
uint64_t smash::Experiment< Modus >::total_hypersurface_crossing_actions_ = 0
private

Total number of particles removed from the evolution in hypersurface crossing actions.

Definition at line 683 of file experiment.h.

◆ discarded_interactions_total_

template<typename Modus >
uint64_t smash::Experiment< Modus >::discarded_interactions_total_ = 0
private

Total number of discarded interactions, because they were invalidated before they could be performed.

Definition at line 689 of file experiment.h.

◆ total_energy_removed_

template<typename Modus >
double smash::Experiment< Modus >::total_energy_removed_ = 0.0
private

Total energy removed from the system in hypersurface crossing actions.

Definition at line 694 of file experiment.h.

◆ total_energy_violated_by_Pythia_

template<typename Modus >
double smash::Experiment< Modus >::total_energy_violated_by_Pythia_ = 0.0
private

Total energy violation introduced by Pythia.

Definition at line 699 of file experiment.h.

◆ kinematic_cuts_for_IC_output_

template<typename Modus >
bool smash::Experiment< Modus >::kinematic_cuts_for_IC_output_ = false
private

This indicates whether kinematic cuts are enabled for the IC output.

Definition at line 702 of file experiment.h.

◆ seed_

template<typename Modus >
int64_t smash::Experiment< Modus >::seed_ = -1
private

random seed for the next event.

Definition at line 705 of file experiment.h.


The documentation for this class was generated from the following file: