Massive particle. A force can be applied. More...

#include <mparticle.h>

Inheritance diagram for Heed::mparticle:

Public Member Functions
void	check_consistency () const
	Check consistency of kin_energy, gamma_1, speed, speed_of_light and mass.

virtual void	step (std::vector< gparticle * > &secondaries)

virtual void	curvature (int &fs_cf, vec &frelcen, vfloat &fmrange, vfloat prec)

virtual int	force (const point &pt, vec &f, vec &f_perp, vfloat &mrange)

void	new_speed ()
	Set new speed, direction and time for currpos.

	mparticle ()
	Default constructor.

	mparticle (manip_absvol *primvol, const point &pt, const vec &vel, vfloat time, double fmass)
	Constructor, $\gamma - 1$ calculated from the from velocity vector.

virtual void	print (std::ostream &file, int l) const

virtual mparticle *	copy () const

virtual	~mparticle ()
	Destructor.

Public Member Functions inherited from Heed::gparticle
	gparticle ()
	Default constructor.

	gparticle (manip_absvol *primvol, const point &pt, const vec &vel, vfloat time)
	Constructor.

virtual	~gparticle ()
	Destructor.

virtual void	step (std::vector< gparticle * > &secondaries)

virtual void	change_vol (void)

virtual void	curvature (int &fs_cf, vec &frelcen, vfloat &fmrange, vfloat prec)

virtual void	physics_after_new_speed (std::vector< gparticle * > &)

virtual void	physics (std::vector< gparticle * > &)

virtual void	physics_mrange (double &fmrange)

virtual stvpoint	calc_step_to_bord ()
	Produces nextpos.

stvpoint	switch_new_vol ()

virtual void	fly (std::vector< gparticle * > &secondaries)
	Transport the particle.

virtual void	print (std::ostream &file, int l) const

virtual gparticle *	copy () const

Public Member Functions inherited from Heed::RegPassivePtr
	RegPassivePtr (void)

	RegPassivePtr (char fs_ban_del, char fs_ban_sub, char fs_ban_cop=0)

	RegPassivePtr (const RegPassivePtr &f)

RegPassivePtr &	operator= (const RegPassivePtr &f)

CountPP_ns::CountPassivePtr *	book (void) const

void	clear_pointers (void) const

virtual RegPassivePtr *	copy () const

virtual	~RegPassivePtr ()

virtual void	print (std::ostream &file, int l=1) const

void	set_s_ban_del (char fs_ban_del)

char	get_s_ban_del (void) const

void	set_s_ban_sub (char fs_ban_sub)

char	get_s_ban_sub (void) const

void	set_s_ban_cop (char fs_ban_cop)

char	get_s_ban_cop (void) const

void	set_s_allow_del_at_zero_count (char fs_allow_del_at_zero_count)

char	get_s_allow_del_at_zero_count (void) const

long	get_total_number_of_references (void) const

Public Attributes
double	mass
	Mass (not mass * speed_of_light^2)

double	orig_kin_energy

double	orig_gamma_1

double	prev_kin_energy

double	prev_gamma_1

double	curr_kin_energy

double	curr_gamma_1

Public Attributes inherited from Heed::gparticle
bool	s_life

long	nstep
	Step number.

double	total_range_from_origin
	Range from origin to currpos.

long	n_zero_step
	Number of previous steps with zero range (including this step).

stvpoint	origin

stvpoint	prevpos

stvpoint	currpos

stvpoint	nextpos

vec	curr_relcen

Additional Inherited Members
Static Public Member Functions inherited from Heed::RegPassivePtr
static void	set_s_ban_del_ignore (char fs_ban_del_ignore)

static char	get_s_ban_del_ignore (void)

static void	set_s_print_adr_cpp (char fs_print_adr_cpp)

static char	get_s_print_adr_cpp (void)

Static Public Attributes inherited from Heed::gparticle
static long	max_q_zero_step = 100

Detailed Description

Massive particle. A force can be applied.

Definition at line 23 of file mparticle.h.

Constructor & Destructor Documentation

◆ mparticle() [1/2]

Heed::mparticle::mparticle ( )

inline

Default constructor.

Definition at line 69 of file mparticle.h.

69: gparticle(), mass(0.0) {}

Heed::gparticle::gparticle

gparticle()

Default constructor.

Definition: gparticle.h:159

Heed::mparticle::mass

double mass

Mass (not mass * speed_of_light^2)

Definition: mparticle.h:26

Referenced by copy().

◆ mparticle() [2/2]

Heed::mparticle::mparticle	(	manip_absvol *	primvol,
		const point &	pt,
		const vec &	vel,
		vfloat	time,
		double	fmass
	)

Constructor, $\gamma - 1$ calculated from the from velocity vector.

Definition at line 20 of file mparticle.cpp.

    : gparticle(primvol, pt, vel, time),
      mass(fmass) {
 
  mfunname("mparticle::mparticle(...)");
 
  orig_gamma_1 = lorgamma_1(origin.speed / c_light);
  orig_kin_energy = orig_gamma_1 * mass * c_squared;
  prev_gamma_1 = lorgamma_1(prevpos.speed / c_light);
  prev_kin_energy = prev_gamma_1 * mass * c_squared;
  curr_gamma_1 = lorgamma_1(currpos.speed / c_light);
  curr_kin_energy = curr_gamma_1 * mass * c_squared;
  check_consistency();
}

◆ ~mparticle()

virtual Heed::mparticle::~mparticle ( )

inlinevirtual

Destructor.

Definition at line 77 of file mparticle.h.

77{}

Member Function Documentation

◆ check_consistency()

void Heed::mparticle::check_consistency ( ) const

Check consistency of kin_energy, gamma_1, speed, speed_of_light and mass.

Definition at line 36 of file mparticle.cpp.

                                        {
  mfunname("void mparticle::check_consistency() const");
  check_econd11(vecerror, != 0, mcerr);
  double speed = c_light * lorbeta(orig_gamma_1);
  check_econd11a(fabs(speed - origin.speed) / (speed + origin.speed), > 1.0e-10,
                 (*this), mcerr);
  speed = c_light * lorbeta(prev_gamma_1);
  check_econd11a(fabs(speed - prevpos.speed) / (speed + prevpos.speed),
                 > 1.0e-10, (*this), mcerr);
  speed = c_light * lorbeta(curr_gamma_1);
  check_econd11a(fabs(speed - currpos.speed) / (speed + currpos.speed),
                 > 1.0e-10, (*this), mcerr);
  double kin_ener = orig_gamma_1 * mass * c_squared;
  if (kin_ener > 1000.0 * DBL_MIN) {
    check_econd11a(
        fabs(orig_kin_energy - kin_ener) / (orig_kin_energy + kin_ener),
        > 1.0e-9, "kin_ener=" << kin_ener << '\n' << (*this), mcerr);
  }
  kin_ener = prev_gamma_1 * mass * c_squared;
  if (kin_ener > 1000.0 * DBL_MIN) {
    check_econd11a(
        fabs(prev_kin_energy - kin_ener) / (prev_kin_energy + kin_ener),
        > 1.0e-9, "kin_ener=" << kin_ener << '\n' << (*this), mcerr);
  }
  kin_ener = curr_gamma_1 * mass * c_squared;
  if (kin_ener > 1000.0 * DBL_MIN) {
    check_econd11a(
        fabs(curr_kin_energy - kin_ener) / (curr_kin_energy + kin_ener),
        > 1.0e-9, "kin_ener=" << kin_ener << '\n' << (*this), mcerr);
  }
}

Referenced by mparticle(), and new_speed().

◆ copy()

virtual mparticle * Heed::mparticle::copy ( ) const

inlinevirtual

Reimplemented from Heed::gparticle.

Reimplemented in Heed::HeedDeltaElectron, Heed::HeedParticle, Heed::HeedParticle_BGM, and Heed::eparticle.

Definition at line 75 of file mparticle.h.

75{ return new mparticle(*this); }

Heed::mparticle::mparticle

mparticle()

Default constructor.

Definition: mparticle.h:69

◆ curvature()

void Heed::mparticle::curvature	(	int &	fs_cf,
		vec &	frelcen,
		vfloat &	fmrange,
		vfloat	prec
	)

virtual

Set curvature. Calls force(). If force is zero, returns fs_cf=0; frelcen=dv0; If force is zero, and currpos.dir==dv0, makes, in addition, fmrange=0; If currpos.dir==dv0, makes currpos.dir=unit_vec(f); If force is parallel or anti-parallel to dir, makes fs_cf=0; frelcen=dv0; If force is anti-parallel to dir, restricts range till exceeding kinetic energy.

Reimplemented from Heed::gparticle.

Definition at line 96 of file mparticle.cpp.

                                       {
 
  pvecerror(
      "void mparticle::curvature(int& fs_cf, vec& frelcen, vfloat& fmrange)");
  vec f;
  vec f_perp_fl;
  int i = force(currpos.pt, f, f_perp_fl, fmrange);
  vec f_perp = currpos.speed * (currpos.dir || f_perp_fl);
  f += f_perp;
  if (i == 0 || f == dv0) {
    fs_cf = 0;
    frelcen = dv0;
    if (currpos.dir == dv0) fmrange = 0;  // to stay in the place
    return;
  }
  if (currpos.dir == dv0) {
    // starting to move in the direction of force
    currpos.dir = unit_vec(f);
  }
  const int j = check_par(currpos.dir, f, prec);
  if (j != 0) {
    fs_cf = 0;
    frelcen = dv0;
    if (j == -1) {
      // decelerate, search for stop point
      const double ran = curr_kin_energy / f.length();
      if (fmrange > ran) fmrange = ran;
    }
  } else {
    fs_cf = 1;
    vec fn = project_to_plane(f, currpos.dir);  // normal component
    frelcen = unit_vec(fn);
    double len = fn.length();
    vfloat rad =
        (currpos.speed * currpos.speed * (curr_gamma_1 + 1) * mass) / len;
    frelcen *= rad;
  }
  currpos.dirloc = currpos.dir;
  currpos.tid.up_absref(&currpos.dirloc);
}

◆ force()

int Heed::mparticle::force	(	const point &	pt,
		vec &	f,
		vec &	f_perp,
		vfloat &	mrange
	)

virtual

The force is considered to be split in two components. One component, namely f, can be in any direction and is capable of doing work. The other one is always normal to dir and cannot do work. The latter can represent the magnetic component of the Lorentz force. This splitting improve precision of calculation of kinetic energy. But the latter component is not the true force. To derive the force one should do vector multiplication of speed by f_perp, f_perp2 = currpos.speed * (currpos.dir || f_perp_fl2); Return 0 if there is no force, f is initialised to zero anyway. mrange is the distance at which the force should not change much. The dimension of f is [weight] * [lenght] / [time]^2 The dimension of f_perp is [weight] / [time];

Reimplemented in Heed::eparticle.

Definition at line 138 of file mparticle.cpp.

                                                                        {
  f = vec(0, 0, 0);
  f_perp = vec(0, 0, 0);
  mrange = max_vfloat;
  return 0;
}

Referenced by curvature(), and new_speed().

◆ new_speed()

void Heed::mparticle::new_speed ( )

Set new speed, direction and time for currpos.

Definition at line 145 of file mparticle.cpp.

                          {
  pvecerror("void mparticle::new_speed(void)");
  if (currpos.prange == 0.0) {
    check_consistency();
    return;
  }
  vec f1, f2, f_perp1, f_perp2, f_perp_fl1, f_perp_fl2;
  vec f_mean;
  vfloat r1, r2;  // ranges, do not need here
  int i = force(prevpos.pt, f1, f_perp_fl1, r1);
  int j = force(currpos.pt, f2, f_perp_fl2, r2);
  check_econd11a(vecerror, != 0, "position 1, after computing force\n", mcerr);
  f_perp1 = prevpos.speed * (prevpos.dir || f_perp_fl1);
  f_perp2 = currpos.speed * (currpos.dir || f_perp_fl2);
  // Later f_perp are ignored since they can not do the work;
  f_mean = 0.5 * (f1 + f2);
  check_econd11a(vecerror, != 0, "position 2, after computing f_perp\n", mcerr);
 
  if ((i == 0 && j == 0) || f_mean == dv0) {
    curr_kin_energy = prev_kin_energy;
    curr_gamma_1 = prev_gamma_1;
    // speed is preserved by gparticle
    currpos.speed = prevpos.speed;  
  } else {
    vec r = currpos.pt - prevpos.pt;
    double W = 0;  // force * range * cos() = work * cos() (may be negative)
    if (r != dv0) W = f_mean * r;
    // This is work which should lead to increase or decrease the speed
    if (W == 0) {
      curr_kin_energy = prev_kin_energy;
      curr_gamma_1 = prev_gamma_1;
      currpos.speed = prevpos.speed;
    } else {
      curr_kin_energy = prev_kin_energy + W;
      if (curr_kin_energy <= 0) {
        curr_kin_energy = 0;
        currpos.speed = 0;
        curr_gamma_1 = 0;
        currpos.dir = dv0;
      } else {
        double resten = mass * c_squared;
        curr_gamma_1 = curr_kin_energy / resten;
        currpos.speed = c_light * lorbeta(curr_gamma_1);
      }
    }
  }
  if (!(i == 0 && j == 0)) {
    vec fn1 = project_to_plane(f1, prevpos.dir);  // normal component
    vec fn2 = project_to_plane(f2, currpos.dir);  // normal component
    check_econd11a(vecerror, != 0, "position 3, after computing fn2\n", mcerr);
    // mean ortogonal component of working force
    vec mean_fn = 0.5 * (fn1 + fn2);  
    double mean_fn_len = mean_fn.length();
    vec fdir = prevpos.dir;
    if (mean_fn_len > 0.0) {
      vec relcen = unit_vec(mean_fn);
      double mean_speed = (prevpos.speed + currpos.speed) * 0.5;
      vfloat new_rad = (mean_speed * mean_speed *
                        ((prev_gamma_1 + curr_gamma_1) * 0.5 + 1) * mass) /
                       mean_fn_len;
      if (new_rad > 0.0) {
        vfloat ang = currpos.prange / new_rad;  // angle to turn
        fdir.turn(prevpos.dir || relcen, ang);  // direction at the end
      }
    }
    check_econd11a(vecerror, != 0, "position 4\n", mcerr);
    vec mean_f_perp_fl = 0.5 * (f_perp_fl1 + f_perp_fl2);
    double len_mean_f_perp_fl = mean_f_perp_fl.length();
    f_perp2 = currpos.speed * (currpos.dir || f_perp_fl2);
    double mean_f_perp = 0.5 * (f_perp1.length() + f_perp2.length());
    check_econd11a(vecerror, != 0, "position 5\n", mcerr);
    if (len_mean_f_perp_fl > 0.0) {
      vec fdir_proj = project_to_plane(prevpos.dir, mean_f_perp_fl);
      if (!apeq(fdir_proj.length(), 0.0)) {
        check_econd11a(vecerror, != 0, "position 6\n", mcerr);
        double length_proj = currpos.prange * cos2vec(prevpos.dir, fdir_proj);
        check_econd11a(vecerror, != 0, "position 7\n", mcerr);
        double acc =
            mean_f_perp / (((prev_gamma_1 + curr_gamma_1) * 0.5 + 1) * mass);
        double mean_speed = (prevpos.speed + currpos.speed) * 0.5;
        double new_rad = pow(mean_speed * fdir_proj.length(), 2) / acc;
        double ang = length_proj / new_rad;
        if (new_rad > 0 && ang > 0) {
          fdir.turn(mean_f_perp_fl, -ang);  // direction at the end
          check_econd11a(vecerror, != 0, "position 8\n", mcerr);
        }
      }
    }
    currpos.dir = fdir;
    check_econd11a(vecerror, != 0, "position 9, after turn\n", mcerr);
  }
  currpos.dirloc = currpos.dir;
  currpos.tid.up_absref(&currpos.dirloc);
  currpos.time =
      prevpos.time + currpos.prange / (0.5 * (prevpos.speed + currpos.speed));
  check_consistency();
}

Referenced by step().

◆ print()

void Heed::mparticle::print	(	std::ostream &	file,
		int	l
	)		const

virtual

Reimplemented from Heed::gparticle.

Reimplemented in Heed::HeedDeltaElectron, Heed::HeedParticle, Heed::HeedParticle_BGM, and Heed::eparticle.

Definition at line 242 of file mparticle.cpp.

                                                   {
  if (l < 0)
    Ifile << "mparticle: mass=" << mass << " (" << mass / CLHEP::kg << " kg, "
          << mass* c_squared / CLHEP::GeV << " GeV)\n";
  Ifile << "orig_kin_energy=" << orig_kin_energy << " ("
        << orig_kin_energy / CLHEP::GeV << " GeV)"
        << " orig_gamma_1=" << orig_gamma_1 << '\n';
  Ifile << "prev_kin_energy=" << prev_kin_energy << " ("
        << prev_kin_energy / CLHEP::GeV << " GeV)"
        << " prev_gamma_1=" << prev_gamma_1 << '\n';
  Ifile << "curr_kin_energy=" << curr_kin_energy << " ("
        << curr_kin_energy / CLHEP::GeV << " GeV)"
        << " curr_gamma_1=" << curr_gamma_1 << '\n';
  gparticle::print(file, l);
}

Referenced by Heed::HeedDeltaElectron::print(), Heed::HeedParticle::print(), Heed::HeedParticle_BGM::print(), Heed::eparticle::print(), and step().

◆ step()

void Heed::mparticle::step ( std::vector< gparticle * > & secondaries )

virtual

Assign prevpos = currpos and currpos = nextpos, calls change_vol if necessary and update nextpos =calc_step_to_bord().

Reimplemented from Heed::gparticle.

Definition at line 68 of file mparticle.cpp.

                                                       {
  // Make step to nextpos and calculate new step to border
  mfunname("void mparticle::step(...)");
  prevpos = currpos;
  prev_kin_energy = curr_kin_energy;
  prev_gamma_1 = curr_gamma_1;
  currpos = nextpos;
  curr_relcen = dv0;
  total_range_from_origin += currpos.prange;
  nstep++;
  if (currpos.prange == 0) {
    n_zero_step++;
    check_econd12a(n_zero_step, >, max_q_zero_step,
                   "too many zero steps, possible infinite loop\n";
                   print(mcout, 10);, mcerr);
  } else {
    n_zero_step = 0;
  }
  // Calculate new current speed, direction and time.
  new_speed();
  physics_after_new_speed(secondaries);
 
  if (s_life) {
    if (prevpos.tid != currpos.tid) change_vol();
    nextpos = calc_step_to_bord();
  }
}

Member Data Documentation

◆ curr_gamma_1

double Heed::mparticle::curr_gamma_1

Definition at line 33 of file mparticle.h.

Referenced by check_consistency(), curvature(), mparticle(), new_speed(), Heed::HeedParticle_BGM::physics(), Heed::HeedDeltaElectron::physics_after_new_speed(), print(), and step().

◆ curr_kin_energy

double Heed::mparticle::curr_kin_energy

Definition at line 32 of file mparticle.h.

Referenced by check_consistency(), curvature(), Garfield::TrackHeed::GetCluster(), mparticle(), new_speed(), Heed::HeedParticle::physics(), Heed::HeedParticle_BGM::physics(), Heed::HeedDeltaElectron::physics_after_new_speed(), Heed::HeedDeltaElectron::physics_mrange(), print(), step(), and Garfield::TrackHeed::TransportPhoton().

◆ mass

double Heed::mparticle::mass

Mass (not mass * speed_of_light^2)

Definition at line 26 of file mparticle.h.

Referenced by check_consistency(), curvature(), mparticle(), new_speed(), Heed::HeedParticle::physics(), Heed::HeedParticle_BGM::physics(), Heed::HeedDeltaElectron::physics_after_new_speed(), and print().

◆ orig_gamma_1

double Heed::mparticle::orig_gamma_1

Definition at line 29 of file mparticle.h.

Referenced by check_consistency(), mparticle(), and print().

◆ orig_kin_energy

double Heed::mparticle::orig_kin_energy

Definition at line 28 of file mparticle.h.

Referenced by check_consistency(), mparticle(), and print().

◆ prev_gamma_1

double Heed::mparticle::prev_gamma_1

Definition at line 31 of file mparticle.h.

Referenced by check_consistency(), mparticle(), new_speed(), print(), and step().

◆ prev_kin_energy

double Heed::mparticle::prev_kin_energy

Definition at line 30 of file mparticle.h.

Referenced by check_consistency(), mparticle(), new_speed(), Heed::HeedDeltaElectron::physics_after_new_speed(), print(), and step().

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

Public Member Functions

Public Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ mparticle() [1/2]

◆ mparticle() [2/2]

◆ ~mparticle()

Member Function Documentation

◆ check_consistency()

◆ copy()

◆ curvature()

◆ force()

◆ new_speed()

◆ print()

◆ step()

Member Data Documentation

◆ curr_gamma_1

◆ curr_kin_energy

◆ mass

◆ orig_gamma_1

◆ orig_kin_energy

◆ prev_gamma_1

◆ prev_kin_energy