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 135 of file experiment.h.

Inheritance diagram for smash::Experiment< Modus >:

Collaboration diagram for smash::Experiment< Modus >:

Public Member Functions
void	run () override
	Runs the experiment. More...
	Experiment (Configuration config, const bf::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 ()
	Runs the time evolution of an event with fixed-sized 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 (const int evt_num)
	Output at the end of an event. More...
Particles *	particles ()
	Provides external access to SMASH particles. More...
Modus *	modus ()
	Provides external access to SMASH calculation modus. 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
template<typename Container >
bool	perform_action (Action &action, const Container &particles_before_actions)
	Perform the given action. More...
void	create_output (const std::string &format, const std::string &content, const bf::path &output_path, const OutputParameters &par)
	Create a list of output files. More...
void	propagate_and_shine (double to_time)
	Propagate all particles until time to_time without any interactions and shine dileptons. More...
void	run_time_evolution_timestepless (Actions &actions)
	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...

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...
Particles	particles_
	Complete particle list. More...
std::unique_ptr< Potentials >	potentials_
	An instance of potentials class, that stores parameters of potentials, calculates them and their gradients. More...
std::unique_ptr< PauliBlocker >	pauli_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 >	nucleon_has_interacted_ = {}
	nucleon_has_interacted_ labels whether the particles in the nuclei have experienced any collisions or not. More...
bool	projectile_target_interact_ = false
	Whether the projectile and the target collided. More...
std::vector< FourVector >	beam_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< DecayActionsFinderDilepton >	dilepton_finder_
	The Dilepton Action Finder. More...
std::unique_ptr< ActionFinderInterface >	photon_finder_
	The (Scatter) Actions Finder for Direct Photons. More...
int	n_fractional_photons_ = 100
	Number of fractional photons produced per single reaction. More...
std::unique_ptr< DensityLattice >	jmu_B_lat_
	Baryon density on the lattices. More...
std::unique_ptr< DensityLattice >	jmu_I3_lat_
	Isospin projection density on the lattices. More...
std::unique_ptr< DensityLattice >	jmu_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 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 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< EnergyMomentumTensor > >	Tmn_
	Lattices of energy-momentum tensors for printout. 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 fram. More...
bool	printout_lattice_td_ = false
	Whether to print the thermodynamics quantities evaluated on the lattices. More...
std::unique_ptr< GrandCanThermalizer >	thermalizer_
	Instance of class used for forced thermalization. More...
StringProcess *	process_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_
	Number of events. 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	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, squared. More...
QuantumNumbers	conserved_initial_
	The conserved quantities of 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...
double	total_energy_removed_ = 0.0
	Total energy removed from the system in hypersurface crossing actions. 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< ExperimentBase >	create (Configuration config, const bf::path &output_path)
	Factory method that creates and initializes a new Experiment<Modus>. More...

Constructor & Destructor Documentation

template<typename Modus >

smash::Experiment< Modus >::Experiment ( Configuration  config, const bf::path &  output_path  )

explicit

Create a new Experiment.

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

Parameters

[in]	config	The 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_path	The directory where the output files are written.

Member Function Documentation

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 1983 of file experiment.h.

                            {
  const auto &mainlog = logger<LogArea::Main>();
  for (int j = 0; j < nevents_; j++) {
    mainlog.info() << "Event " << j;

    // Sample initial particles, start clock, some printout and book-keeping
    initialize_new_event();
    /* In the ColliderModus, if the first collisions within the same nucleus are
     * forbidden, 'nucleon_has_interacted_', which records whether a nucleon has
     * collided with another nucleon, is initialized equal to false. If allowed,
     * 'nucleon_has_interacted' is initialized equal to true, which means these
     * incoming particles have experienced some fake scatterings, they can
     * therefore collide with each other later on since these collisions are not
     * "first" to them. */
    if (modus_.is_collider()) {
      if (!modus_.cll_in_nucleus()) {
        nucleon_has_interacted_.assign(modus_.total_N_number(), false);
      } else {
        nucleon_has_interacted_.assign(modus_.total_N_number(), true);
      }
    }
    /* In the ColliderModus, if Fermi motion is frozen, assign the beam momenta
     * to the nucleons in both the projectile and the target. */
    if (modus_.is_collider() && modus_.fermi_motion() == FermiMotion::Frozen) {
      for (int i = 0; i < modus_.total_N_number(); i++) {
        const auto mass_beam = particles_.copy_to_vector()[i].effective_mass();
        const auto v_beam = i < modus_.proj_N_number()
                                ? modus_.velocity_projectile()
                                : modus_.velocity_target();
        const auto gamma = 1.0 / std::sqrt(1.0 - v_beam * v_beam);
        beam_momentum_.emplace_back(FourVector(gamma * mass_beam, 0.0, 0.0,
                                               gamma * v_beam * mass_beam));
      }
    }

    // Output at event start
    for (const auto &output : outputs_) {
      output->at_eventstart(particles_, j);
    }

    run_time_evolution();

    if (force_decays_) {
      do_final_decays();
    }

    // Output at event end
    final_output(j);
  }
}

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 1326 of file experiment.h.

                                             {
  const auto &log = logger<LogArea::Experiment>();

  random::set_seed(seed_);
  log.info() << "random number seed: " << seed_;
  /* Set seed for the next event. It has to be positive, so it can be entered
   * in the config.
   *
   * We have to be careful about the minimal integer, whose absolute value
   * cannot be represented. */
  int64_t r = random::advance();
  while (r == INT64_MIN) {
    r = random::advance();
  }
  seed_ = std::abs(r);
  /* Set the random seed used in PYTHIA hadronization
   * to be same with the SMASH one.
   * In this way we ensure that the results are reproducible
   * for every event if one knows SMASH random seed. */
  if (process_string_ptr_ != NULL) {
    process_string_ptr_->init_pythia_hadron_rndm();
  }

  particles_.reset();

  // Sample particles according to the initial conditions
  double start_time = modus_.initial_conditions(&particles_, parameters_);
  /* For box modus make sure that particles are in the box. In principle, after
   * a correct initialization they should be, so this is just playing it safe.
   */
  modus_.impose_boundary_conditions(&particles_, outputs_);
  // Reset the simulation clock
  double timestep = delta_time_startup_;

  switch (time_step_mode_) {
    case TimeStepMode::Fixed:
      break;
    case TimeStepMode::None:
      timestep = end_time_ - start_time;
      // Take care of the box modus + timestepless propagation
      const double max_dt = modus_.max_timestep(max_transverse_distance_sqr_);
      if (max_dt > 0. && max_dt < timestep) {
        timestep = max_dt;
      }
      break;
  }
  std::unique_ptr<UniformClock> clock_for_this_event;
  if (modus_.is_list() && (timestep < 0.0)) {
    throw std::runtime_error(
        "Timestep for the given event is negative. \n"
        "This might happen if the formation times of the input particles are "
        "larger than the specified end time of the simulation.");
  }
  clock_for_this_event = make_unique<UniformClock>(start_time, timestep);
  parameters_.labclock = std::move(clock_for_this_event);

  // Reset the output clock
  parameters_.outputclock->reset(start_time, true);
  // remove time before starting time in case of custom output times.
  parameters_.outputclock->remove_times_in_past(start_time);

  log.debug("Lab clock: t_start = ", parameters_.labclock->current_time(),
            ", dt = ", parameters_.labclock->timestep_duration());

  /* Save the initial conserved quantum numbers and total momentum in
   * the system for conservation checks */
  conserved_initial_ = QuantumNumbers(particles_);
  wall_actions_total_ = 0;
  previous_wall_actions_total_ = 0;
  interactions_total_ = 0;
  previous_interactions_total_ = 0;
  total_pauli_blocked_ = 0;
  projectile_target_interact_ = false;
  total_hypersurface_crossing_actions_ = 0;
  total_energy_removed_ = 0.0;
  // Print output headers
  log.info() << hline;
  log.info() << "Time [fm]   Ediff [GeV]    Scatt.|Decays   "
                "Particles         Timing";
  log.info() << hline;
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::run_time_evolution

(

)

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

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

Definition at line 1558 of file experiment.h.

                                           {
  Actions actions;

  const auto &log = logger<LogArea::Experiment>();

  log.info() << format_measurements(particles_, 0u, conserved_initial_,
                                    time_start_,
                                    parameters_.labclock->current_time());

  while (parameters_.labclock->current_time() < end_time_) {
    const double t = parameters_.labclock->current_time();
    const double dt =
        std::min(parameters_.labclock->timestep_duration(), end_time_ - t);
    log.debug("Timestepless propagation for next ", dt, " fm/c.");

    // Perform forced thermalization if required
    if (thermalizer_ &&
        thermalizer_->is_time_to_thermalize(parameters_.labclock)) {
      const bool ignore_cells_under_treshold = true;
      thermalizer_->update_thermalizer_lattice(particles_, density_param_,
                                               ignore_cells_under_treshold);
      const double current_t = parameters_.labclock->current_time();
      thermalizer_->thermalize(particles_, current_t,
                               parameters_.testparticles);
      ThermalizationAction th_act(*thermalizer_, current_t);
      if (th_act.any_particles_thermalized()) {
        perform_action(th_act, particles_);
      }
    }

    if (particles_.size() > 0 && action_finders_.size() > 0) {
      /* (1.a) Create grid. */
      double min_cell_length = compute_min_cell_length(dt);
      log.debug("Creating grid with minimal cell length ", min_cell_length);
      const auto &grid =
          use_grid_ ? modus_.create_grid(particles_, min_cell_length, dt)
                    : modus_.create_grid(particles_, min_cell_length, dt,
                                         CellSizeStrategy::Largest);

      const double cell_vol = grid.cell_volume();

      /* (1.b) Iterate over cells and find actions. */
      grid.iterate_cells(
          [&](const ParticleList &search_list) {
            for (const auto &finder : action_finders_) {
              actions.insert(finder->find_actions_in_cell(
                  search_list, dt, cell_vol, beam_momentum_));
            }
          },
          [&](const ParticleList &search_list,
              const ParticleList &neighbors_list) {
            for (const auto &finder : action_finders_) {
              actions.insert(finder->find_actions_with_neighbors(
                  search_list, neighbors_list, dt, beam_momentum_));
            }
          });
    }

    /* \todo (optimizations) Adapt timestep size here */

    /* (2) Propagation from action to action until the end of timestep */
    run_time_evolution_timestepless(actions);

    /* (3) Update potentials (if computed on the lattice) and
     *     compute new momenta according to equations of motion */
    if (potentials_) {
      update_potentials();
      update_momenta(&particles_, parameters_.labclock->timestep_duration(),
                     *potentials_, FB_lat_.get(), FI3_lat_.get());
    }

    /* (4) Expand universe if non-minkowskian metric; updates
     *     positions and momenta according to the selected expansion */
    if (metric_.mode_ != ExpansionMode::NoExpansion) {
      expand_space_time(&particles_, parameters_, metric_);
    }

    ++(*parameters_.labclock);

    /* (5) Check conservation laws.
     *
     * Check conservation of conserved quantities if potentials and string
     * fragmentation are off.  If potentials are on then momentum is conserved
     * only in average.  If string fragmentation is on, then energy and
     * momentum are only very roughly conserved in high-energy collisions. */
    if (!potentials_ && !parameters_.strings_switch &&
        metric_.mode_ == ExpansionMode::NoExpansion && !IC_output_switch_) {
      std::string err_msg = conserved_initial_.report_deviations(particles_);
      if (!err_msg.empty()) {
        log.error() << err_msg;
        throw std::runtime_error("Violation of conserved quantities!");
      }
    }
  }

  if (pauli_blocker_) {
    log.info("Interactions: Pauli-blocked/performed = ", total_pauli_blocked_,
             "/", interactions_total_ - wall_actions_total_);
  }
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::do_final_decays

(

)

Performs the final decays of an event.

Definition at line 1886 of file experiment.h.

                                        {
  /* At end of time evolution: Force all resonances to decay. In order to handle
   * decay chains, we need to loop until no further actions occur. */
  uint64_t interactions_old;
  const auto particles_before_actions = particles_.copy_to_vector();
  do {
    Actions actions;

    interactions_old = interactions_total_;

    // Dileptons: shining of remaining resonances
    if (dilepton_finder_ != nullptr) {
      for (const auto &output : outputs_) {
        dilepton_finder_->shine_final(particles_, output.get(), true);
      }
    }
    // Find actions.
    for (const auto &finder : action_finders_) {
      actions.insert(finder->find_final_actions(particles_));
    }
    // Perform actions.
    while (!actions.is_empty()) {
      perform_action(*actions.pop(), particles_before_actions);
    }
    // loop until no more decays occur
  } while (interactions_total_ > interactions_old);

  // Dileptons: shining of stable particles at the end
  if (dilepton_finder_ != nullptr) {
    for (const auto &output : outputs_) {
      dilepton_finder_->shine_final(particles_, output.get(), false);
    }
  }
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::final_output

(

const int

evt_num

)

Output at the end of an event.

Parameters

[in]

evt_num

Number of the event

Definition at line 1922 of file experiment.h.

                                                      {
  const auto &log = logger<LogArea::Experiment>();
  /* make sure the experiment actually ran (note: we should compare this
   * to the start time, but we don't know that. Therefore, we check that
   * the time is positive, which should heuristically be the same). */
  if (likely(parameters_.labclock > 0)) {
    const uint64_t wall_actions_this_interval =
        wall_actions_total_ - previous_wall_actions_total_;
    const uint64_t interactions_this_interval = interactions_total_ -
                                                previous_interactions_total_ -
                                                wall_actions_this_interval;
    log.info() << format_measurements(particles_, interactions_this_interval,
                                      conserved_initial_, time_start_,
                                      end_time_);
    if (IC_output_switch_ && (particles_.size() == 0)) {
      // Verify there is no more energy in the system if all particles were
      // removed when crossing the hypersurface
      const double remaining_energy =
          conserved_initial_.momentum().x0() - total_energy_removed_;
      if (remaining_energy > really_small) {
        throw std::runtime_error(
            "There is remaining energy in the system although all particles "
            "were removed.\n"
            "E_remain = " +
            std::to_string(remaining_energy) + " [GeV]");
      } else {
        log.info() << hline;
        log.info() << "Time real: " << SystemClock::now() - time_start_;
        log.info() << "Interactions before reaching hypersurface: "
                   << interactions_total_ - wall_actions_total_ -
                          total_hypersurface_crossing_actions_;
        log.info() << "Total number of particles removed on hypersurface: "
                   << total_hypersurface_crossing_actions_;
      }
    } else {
      log.info() << hline;
      log.info() << "Time real: " << SystemClock::now() - time_start_;
      log.info() << "Final interaction number: "
                 << interactions_total_ - wall_actions_total_;
    }

    // Check if there are unformed particles
    int unformed_particles_count = 0;
    for (const auto &particle : particles_) {
      if (particle.formation_time() > end_time_) {
        unformed_particles_count++;
      }
    }
    if (unformed_particles_count > 0) {
      log.warn("End time might be too small. ", unformed_particles_count,
               " unformed particles were found at the end of the evolution.");
    }
  }

  for (const auto &output : outputs_) {
    output->at_eventend(particles_, evt_num, modus_.impact_parameter(),
                        !projectile_target_interact_);
  }
}

Here is the call graph for this function:

template<typename Modus >

Particles* smash::Experiment< Modus >::particles ( )

inline

Provides external access to SMASH particles.

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

Definition at line 220 of file experiment.h.

{ return &particles_; }

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 226 of file experiment.h.

{ return &modus_; }

Here is the call graph for this function:

template<typename Modus >

template<typename Container >

bool smash::Experiment< Modus >::perform_action ( Action &  action, const Container &  particles_before_actions  )

private

Perform the given action.

Template Parameters

Container

type that holds the particles before the action.

Parameters

[in]	action	The action to perform. If it performs, it'll modify the private member particles_.
[in]	particles_before_actions	A container with the ParticleData from this time step before any actions were performed.

Returns

False if the action is rejected either due to invalidity or Pauli-blocking, or true if it's accepted and performed.

template<typename Modus >

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

private

Create a list of output files.

Parameters

[in]	format	Format of the output file (e.g. Root, Oscar, Vtk)
[in]	content	Content of the output (e.g. particles, collisions)
[in]	output_path	Path of the output file
[in]	par	Output options.(e.g. Extended)

Definition at line 564 of file experiment.h.

                                                                       {
  const auto &log = logger<LogArea::Experiment>();
  log.info() << "Adding output " << content << " of format " << format
             << std::endl;

  if (format == "VTK" && content == "Particles") {
    outputs_.emplace_back(
        make_unique<VtkOutput>(output_path, content, out_par));
  } else if (format == "Root") {
#ifdef SMASH_USE_ROOT
    if (content == "Initial_Conditions") {
      outputs_.emplace_back(
          make_unique<RootOutput>(output_path, "SMASH_IC", out_par));
    } else {
      outputs_.emplace_back(
          make_unique<RootOutput>(output_path, content, out_par));
    }
#else
    log.error("Root output requested, but Root support not compiled in");
#endif
  } else if (format == "Binary") {
    if (content == "Collisions" || content == "Dileptons" ||
        content == "Photons") {
      outputs_.emplace_back(
          make_unique<BinaryOutputCollisions>(output_path, content, out_par));
    } else if (content == "Particles") {
      outputs_.emplace_back(
          make_unique<BinaryOutputParticles>(output_path, content, out_par));
    } else if (content == "Initial_Conditions") {
      outputs_.emplace_back(make_unique<BinaryOutputInitialConditions>(
          output_path, content, out_par));
    }
  } else if (format == "Oscar1999" || format == "Oscar2013") {
    outputs_.emplace_back(
        create_oscar_output(format, content, output_path, out_par));
  } else if (content == "Thermodynamics" && format == "ASCII") {
    outputs_.emplace_back(
        make_unique<ThermodynamicOutput>(output_path, content, out_par));
  } else if (content == "Thermodynamics" && format == "VTK") {
    printout_lattice_td_ = true;
    outputs_.emplace_back(
        make_unique<VtkOutput>(output_path, content, out_par));
  } else if (content == "Initial_Conditions" && format == "ASCII") {
    outputs_.emplace_back(
        make_unique<ICOutput>(output_path, "SMASH_IC", out_par));
  } else {
    log.error() << "Unknown combination of format (" << format
                << ") and content (" << content << "). Fix the config.";
  }
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::propagate_and_shine ( double  to_time )

private

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

Parameters

[in]

to_time

Time at the end of propagation [fm/c]

Definition at line 1660 of file experiment.h.

                                                          {
  const double dt =
      propagate_straight_line(&particles_, to_time, beam_momentum_);
  if (dilepton_finder_ != nullptr) {
    for (const auto &output : outputs_) {
      dilepton_finder_->shine(particles_, output.get(), dt);
    }
  }
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::run_time_evolution_timestepless ( Actions &  actions )

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. In a time stepless case, the time interval should be equal to the whole evolution time, while in the case with time step, the intervals are given by the time steps.

Parameters

[in,out]

actions

Actions occur during a certain time interval. They provide the ending times of the propagations and are updated during the time interval.

Definition at line 1685 of file experiment.h.

                                                                        {
  const auto &log = logger<LogArea::Experiment>();

  const double start_time = parameters_.labclock->current_time();
  const double end_time =
      std::min(parameters_.labclock->next_time(), end_time_);
  double time_left = end_time - start_time;
  log.debug("Timestepless propagation: ", "Actions size = ", actions.size(),
            ", start time = ", start_time, ", end time = ", end_time);

  // iterate over all actions
  while (!actions.is_empty()) {
    // get next action
    ActionPtr act = actions.pop();
    if (!act->is_valid(particles_)) {
      log.debug(~einhard::DRed(), "✘ ", act, " (discarded: invalid)");
      continue;
    }
    if (act->time_of_execution() > end_time) {
      log.error(act, " scheduled later than end time: t_action[fm/c] = ",
                act->time_of_execution(), ", t_end[fm/c] = ", end_time);
    }
    log.debug(~einhard::Green(), "✔ ", act);

    while (next_output_time() <= act->time_of_execution()) {
      log.debug("Propagating until output time: ", next_output_time());
      propagate_and_shine(next_output_time());
      ++(*parameters_.outputclock);
      intermediate_output();
    }

    /* (1) Propagate to the next action. */
    log.debug("Propagating until next action ", act,
              ", action time = ", act->time_of_execution());
    propagate_and_shine(act->time_of_execution());

    /* (2) Perform action.
     *
     * Update the positions of the incoming particles, because the information
     * in the action object will be outdated as the particles have been
     * propagated since the construction of the action. */
    act->update_incoming(particles_);
    const bool performed = perform_action(*act, particles_);

    /* No need to update actions for outgoing particles
     * if the action is not performed. */
    if (!performed) {
      continue;
    }

    /* (3) Update actions for newly-produced particles. */

    time_left = end_time - act->time_of_execution();
    const ParticleList &outgoing_particles = act->outgoing_particles();
    // Cell volume set to zero, since there is no grid
    const double cell_vol = 0.0;
    for (const auto &finder : action_finders_) {
      // Outgoing particles can still decay, cross walls...
      actions.insert(finder->find_actions_in_cell(outgoing_particles, time_left,
                                                  cell_vol, beam_momentum_));
      // ... and collide with other particles.
      actions.insert(finder->find_actions_with_surrounding_particles(
          outgoing_particles, particles_, time_left, beam_momentum_));
    }

    check_interactions_total(interactions_total_);
  }

  while (next_output_time() <= end_time) {
    log.debug("Propagating until output time: ", next_output_time());
    propagate_and_shine(next_output_time());
    ++(*parameters_.outputclock);
    // Avoid duplicating printout at event end time
    if (parameters_.outputclock->current_time() < end_time_) {
      intermediate_output();
    }
  }
  log.debug("Propagating to time ", end_time);
  propagate_and_shine(end_time);
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::intermediate_output ( )

private

Intermediate output during an event.

Definition at line 1767 of file experiment.h.

                                            {
  const auto &log = logger<LogArea::Experiment>();
  const uint64_t wall_actions_this_interval =
      wall_actions_total_ - previous_wall_actions_total_;
  previous_wall_actions_total_ = wall_actions_total_;
  const uint64_t interactions_this_interval = interactions_total_ -
                                              previous_interactions_total_ -
                                              wall_actions_this_interval;
  previous_interactions_total_ = interactions_total_;
  log.info() << format_measurements(particles_, interactions_this_interval,
                                    conserved_initial_, time_start_,
                                    parameters_.outputclock->current_time());
  const LatticeUpdate lat_upd = LatticeUpdate::AtOutput;
  // save evolution data
  for (const auto &output : outputs_) {
    if (output->is_dilepton_output() || output->is_photon_output() ||
        output->is_IC_output()) {
      continue;
    }

    output->at_intermediate_time(particles_, parameters_.outputclock,
                                 density_param_);

    // Thermodynamic output on the lattice versus time
    switch (dens_type_lattice_printout_) {
      case DensityType::Baryon:
        update_lattice(jmu_B_lat_.get(), lat_upd, DensityType::Baryon,
                       density_param_, particles_, false);
        output->thermodynamics_output(ThermodynamicQuantity::EckartDensity,
                                      DensityType::Baryon, *jmu_B_lat_);
        break;
      case DensityType::BaryonicIsospin:
        update_lattice(jmu_I3_lat_.get(), lat_upd, DensityType::BaryonicIsospin,
                       density_param_, particles_, false);
        output->thermodynamics_output(ThermodynamicQuantity::EckartDensity,
                                      DensityType::BaryonicIsospin,
                                      *jmu_I3_lat_);
        break;
      case DensityType::None:
        break;
      default:
        update_lattice(jmu_custom_lat_.get(), lat_upd,
                       dens_type_lattice_printout_, density_param_, particles_,
                       false);
        output->thermodynamics_output(ThermodynamicQuantity::EckartDensity,
                                      dens_type_lattice_printout_,
                                      *jmu_custom_lat_);
    }
    if (printout_tmn_ || printout_tmn_landau_ || printout_v_landau_) {
      update_lattice(Tmn_.get(), lat_upd, dens_type_lattice_printout_,
                     density_param_, particles_);
      if (printout_tmn_) {
        output->thermodynamics_output(ThermodynamicQuantity::Tmn,
                                      dens_type_lattice_printout_, *Tmn_);
      }
      if (printout_tmn_landau_) {
        output->thermodynamics_output(ThermodynamicQuantity::TmnLandau,
                                      dens_type_lattice_printout_, *Tmn_);
      }
      if (printout_v_landau_) {
        output->thermodynamics_output(ThermodynamicQuantity::LandauVelocity,
                                      dens_type_lattice_printout_, *Tmn_);
      }
    }

    if (thermalizer_) {
      output->thermodynamics_output(*thermalizer_);
    }
  }
}

Here is the call graph for this function:

template<typename Modus >

void smash::Experiment< Modus >::update_potentials ( )

private

Recompute potentials on lattices if necessary.

Definition at line 1839 of file experiment.h.

                                          {
  if (potentials_) {
    if (potentials_->use_symmetry() && jmu_I3_lat_ != nullptr) {
      update_lattice(jmu_I3_lat_.get(), LatticeUpdate::EveryTimestep,
                     DensityType::BaryonicIsospin, density_param_, particles_,
                     true);
    }
    if ((potentials_->use_skyrme() || potentials_->use_symmetry()) &&
        jmu_B_lat_ != nullptr) {
      update_lattice(jmu_B_lat_.get(), LatticeUpdate::EveryTimestep,
                     DensityType::Baryon, density_param_, particles_, true);
      const size_t UBlattice_size = UB_lat_->size();
      assert(UBlattice_size == UI3_lat_->size());
      for (size_t i = 0; i < UBlattice_size; i++) {
        auto jB = (*jmu_B_lat_)[i];
        const FourVector flow_four_velocity_B =
            std::abs(jB.density()) > really_small ? jB.jmu_net() / jB.density()
                                                  : FourVector();
        double baryon_density = jB.density();
        ThreeVector baryon_grad_rho = jB.grad_rho();
        ThreeVector baryon_dj_dt = jB.dj_dt();
        ThreeVector baryon_rot_j = jB.rot_j();
        if (potentials_->use_skyrme()) {
          (*UB_lat_)[i] =
              flow_four_velocity_B * potentials_->skyrme_pot(baryon_density);
          (*FB_lat_)[i] = potentials_->skyrme_force(
              baryon_density, baryon_grad_rho, baryon_dj_dt, baryon_rot_j);
        }
        if (potentials_->use_symmetry() && jmu_I3_lat_ != nullptr) {
          auto jI3 = (*jmu_I3_lat_)[i];
          const FourVector flow_four_velocity_I3 =
              std::abs(jI3.density()) > really_small
                  ? jI3.jmu_net() / jI3.density()
                  : FourVector();
          (*UI3_lat_)[i] =
              flow_four_velocity_I3 *
              potentials_->symmetry_pot(jI3.density(), baryon_density);
          (*FI3_lat_)[i] = potentials_->symmetry_force(
              jI3.density(), jI3.grad_rho(), jI3.dj_dt(), jI3.rot_j(),
              baryon_density, baryon_grad_rho, baryon_dj_dt, baryon_rot_j);
        }
      }
    }
  }
}

Here is the call graph for this function:

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]

dt

The current time step size [fm/c]

Returns

The minimal required size of cells

Definition at line 290 of file experiment.h.

                                                  {
    return std::sqrt(4 * dt * dt + max_transverse_distance_sqr_);
  }

template<typename Modus >

double smash::Experiment< Modus >::next_output_time ( ) const

inlineprivate

Shortcut for next output time.

Definition at line 295 of file experiment.h.

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

Friends And Related Function Documentation

template<typename Modus >

friend class ExperimentBase

friend

Definition at line 165 of file experiment.h.

Member Data Documentation

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 304 of file experiment.h.

template<typename Modus >

DensityParameters smash::Experiment< Modus >::density_param_

private

Structure to precalculate and hold parameters for density computations.

Definition at line 307 of file experiment.h.

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 313 of file experiment.h.

template<typename Modus >

Particles smash::Experiment< Modus >::particles_

private

Complete particle list.

Definition at line 316 of file experiment.h.

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 322 of file experiment.h.

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 328 of file experiment.h.

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 334 of file experiment.h.

template<typename Modus >

OutputPtr smash::Experiment< Modus >::dilepton_output_

private

The Dilepton output.

Definition at line 337 of file experiment.h.

template<typename Modus >

OutputPtr smash::Experiment< Modus >::photon_output_

private

The Photon output.

Definition at line 340 of file experiment.h.

template<typename Modus >

std::vector<bool> smash::Experiment< Modus >::nucleon_has_interacted_ = {}

private

nucleon_has_interacted_ labels whether the particles in the nuclei have experienced any collisions or not.

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

Definition at line 347 of file experiment.h.

template<typename Modus >

bool smash::Experiment< Modus >::projectile_target_interact_ = false

private

Whether the projectile and the target collided.

Definition at line 351 of file experiment.h.

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 358 of file experiment.h.

template<typename Modus >

std::vector<std::unique_ptr<ActionFinderInterface> > smash::Experiment< Modus >::action_finders_

private

The Action finder objects.

Definition at line 361 of file experiment.h.

template<typename Modus >

std::unique_ptr<DecayActionsFinderDilepton> smash::Experiment< Modus >::dilepton_finder_

private

The Dilepton Action Finder.

Definition at line 364 of file experiment.h.

template<typename Modus >

std::unique_ptr<ActionFinderInterface> smash::Experiment< Modus >::photon_finder_

private

The (Scatter) Actions Finder for Direct Photons.

Definition at line 367 of file experiment.h.

template<typename Modus >

int smash::Experiment< Modus >::n_fractional_photons_ = 100

private

Number of fractional photons produced per single reaction.

Definition at line 370 of file experiment.h.

template<typename Modus >

std::unique_ptr<DensityLattice> smash::Experiment< Modus >::jmu_B_lat_

private

Baryon density on the lattices.

Definition at line 373 of file experiment.h.

template<typename Modus >

std::unique_ptr<DensityLattice> smash::Experiment< Modus >::jmu_I3_lat_

private

Isospin projection density on the lattices.

Definition at line 376 of file experiment.h.

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 385 of file experiment.h.

template<typename Modus >

DensityType smash::Experiment< Modus >::dens_type_lattice_printout_ = DensityType::None

private

Type of density for lattice printout.

Definition at line 388 of file experiment.h.

template<typename Modus >

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

private

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

Definition at line 394 of file experiment.h.

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 400 of file experiment.h.

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 force.

Definition at line 404 of file experiment.h.

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 408 of file experiment.h.

template<typename Modus >

std::unique_ptr<RectangularLattice<EnergyMomentumTensor> > smash::Experiment< Modus >::Tmn_

private

Lattices of energy-momentum tensors for printout.

Definition at line 411 of file experiment.h.

template<typename Modus >

bool smash::Experiment< Modus >::printout_tmn_ = false

private

Whether to print the energy-momentum tensor.

Definition at line 414 of file experiment.h.

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 417 of file experiment.h.

template<typename Modus >

bool smash::Experiment< Modus >::printout_v_landau_ = false

private

Whether to print the 4-velocity in Landau fram.

Definition at line 420 of file experiment.h.

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 423 of file experiment.h.

template<typename Modus >

std::unique_ptr<GrandCanThermalizer> smash::Experiment< Modus >::thermalizer_

private

Instance of class used for forced thermalization.

Definition at line 426 of file experiment.h.

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 432 of file experiment.h.

template<typename Modus >

const int smash::Experiment< Modus >::nevents_

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 448 of file experiment.h.

template<typename Modus >

const double smash::Experiment< Modus >::end_time_

private

simulation time at which the evolution is stopped.

Definition at line 451 of file experiment.h.

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 458 of file experiment.h.

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 464 of file experiment.h.

template<typename Modus >

const bool smash::Experiment< Modus >::use_grid_

private

This indicates whether to use the grid.

Definition at line 467 of file experiment.h.

template<typename Modus >

const ExpansionProperties smash::Experiment< Modus >::metric_

private

This struct contains information on the metric to be used.

Definition at line 470 of file experiment.h.

template<typename Modus >

const bool smash::Experiment< Modus >::dileptons_switch_

private

This indicates whether dileptons are switched on.

Definition at line 473 of file experiment.h.

template<typename Modus >

const bool smash::Experiment< Modus >::photons_switch_

private

This indicates whether photons are switched on.

Definition at line 476 of file experiment.h.

template<typename Modus >

const bool smash::Experiment< Modus >::IC_output_switch_

private

This indicates whether the IC output is enabled.

Definition at line 479 of file experiment.h.

template<typename Modus >

const TimeStepMode smash::Experiment< Modus >::time_step_mode_

private

This indicates whether to use time steps.

Definition at line 482 of file experiment.h.

template<typename Modus >

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

private

Maximal distance at which particles can interact, squared.

Definition at line 485 of file experiment.h.

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 494 of file experiment.h.

template<typename Modus >

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

private

system starting time of the simulation

Definition at line 497 of file experiment.h.

template<typename Modus >

DensityType smash::Experiment< Modus >::dens_type_ = DensityType::None

private

Type of density to be written to collision headers.

Definition at line 500 of file experiment.h.

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 506 of file experiment.h.

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 512 of file experiment.h.

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 518 of file experiment.h.

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 524 of file experiment.h.

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 530 of file experiment.h.

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 536 of file experiment.h.

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 542 of file experiment.h.

template<typename Modus >

int64_t smash::Experiment< Modus >::seed_ = -1

private

random seed for the next event.

Definition at line 545 of file experiment.h.

---

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

• src/include/smash/experiment.h