Commit 0cb54aac authored by kraus's avatar kraus Committed by kraus
Browse files

fix reference for case when all particles are in degrader material; refactor...

fix reference for case when all particles are in degrader material; refactor CollimatorPhysics and tried to clean up mix of units and convertions
parent 496b7ddf
......@@ -873,7 +873,7 @@ void ParallelTTracker::computeParticleMatterInteraction(IndexMap::value_t elemen
itsBunch_m->setdT(timeDifference / numSteps);
BorisPusher pusher(itsReference);
for (unsigned int i = 0; i < numSteps; ++ i) {
updateReferenceParticle(pusher);
updateRefToLabCSTrafo(pusher);
}
itsBunch_m->setdT(origdT);
}
......
......@@ -24,7 +24,9 @@
#include "Fields/Fieldmap.h"
#include "Structure/LossDataSink.h"
#include "Utilities/Options.h"
#include "Solvers/ParticleMatterInteractionHandler.hh"
#include "Physics/Physics.h"
#include "Utilities/Util.h"
#include <memory>
extern Inform *gmsg;
......@@ -122,6 +124,22 @@ bool Degrader::apply(const size_t &i, const double &t, Vector_t &E, Vector_t &B)
return false;
}
bool Degrader::applyToReferenceParticle(const Vector_t &R,
const Vector_t &P,
const double &t,
Vector_t &E,
Vector_t &B) {
if (!isInMaterial(R(2))) return false;
const double eV2GeV = 1e-9;
double Ekin = Util::getEnergy(P, RefPartBunch_m->getM() * eV2GeV);
bool isDead = getParticleMatterInteraction()->computeEnergyLoss(Ekin, RefPartBunch_m->getdT(), false);
double deltaP = Util::getP(Ekin, RefPartBunch_m->getM() * eV2GeV) - euclidean_norm(P);
E(2) += deltaP * RefPartBunch_m->getM() / (RefPartBunch_m->getdT() * RefPartBunch_m->getQ() * Physics::c);
return isDead;
}
void Degrader::initialise(PartBunchBase<double, 3> *bunch, double &startField, double &endField) {
RefPartBunch_m = bunch;
endField = startField + getElementLength();
......@@ -195,4 +213,4 @@ ElementBase::ElementType Degrader::getType() const {
string Degrader::getDegraderShape() {
return "DEGRADER";
}
}
\ No newline at end of file
......@@ -62,6 +62,12 @@ public:
virtual bool apply(const size_t &i, const double &t, Vector_t &E, Vector_t &B);
virtual bool applyToReferenceParticle(const Vector_t &R,
const Vector_t &P,
const double &t,
Vector_t &E,
Vector_t &B);
virtual void initialise(PartBunchBase<double, 3> *bunch, double &startField, double &endField);
virtual void initialise(PartBunchBase<double, 3> *bunch);
......
......@@ -56,7 +56,7 @@ public:
std::string getName() {return element_ref_m->getName();}
size_t getParticlesInMat() {return locPartsInMat_m;}
unsigned getRediffused() {return rediffusedStat_m;}
unsigned int getNumEntered() {return bunchToMatStat_m;}
inline void doPhysics(PartBunchBase<double, 3> *bunch);
void setCyclotron(Cyclotron* cycl) { cycl_m = cycl; };
......@@ -78,6 +78,10 @@ private:
void SecondaryParticles(PartBunchBase<double, 3> *bunch, size_t &i);
bool computeEnergyLoss(double &Eng,
const double deltat,
bool includeFluctuations = true) const { return false;}
Cyclotron *cycl_m;
BeamStripping *bstp_m;
......@@ -118,4 +122,4 @@ private:
static const double CSCoefSingleCapt_H[3][9];
};
#endif //BEAMSTRIPPINGPHYSICS_HH
#endif //BEAMSTRIPPINGPHYSICS_HH
\ No newline at end of file
This diff is collapsed.
......@@ -67,6 +67,14 @@ typedef struct { // struct for description of particle in material
} PART;
#endif
struct InsideTester {
virtual ~InsideTester()
{ }
virtual
bool checkHit(const Vector_t &R, const Vector_t &P, double dt) = 0;
};
class CollimatorPhysics: public ParticleMatterInteractionHandler {
public:
......@@ -83,30 +91,42 @@ public:
bool stillActive();
bool stillAlive(PartBunchBase<double, 3> *bunch);
inline double getTime() {return T_m;}
std::string getName() { return FN_m;}
size_t getParticlesInMat() { return locPartsInMat_m;}
unsigned getRediffused() { return rediffusedStat_m;}
unsigned int getNumEntered() {return bunchToMatStat_m;}
inline void doPhysics(PartBunchBase<double, 3> *bunch);
double getTime();
std::string getName();
size_t getParticlesInMat();
unsigned getRediffused();
unsigned int getNumEntered();
void doPhysics(PartBunchBase<double, 3> *bunch);
virtual bool computeEnergyLoss(double &Eng, const double deltat, bool includeFluctuations = true) const;
private:
void Material();
void CoulombScat(Vector_t &R, Vector_t &P, const double &deltat);
void EnergyLoss(double &Eng, bool &pdead, const double &deltat);
bool EnergyLoss(double &Eng, const double &deltat);
void configureMaterialParameters();
void computeCoulombScattering(Vector_t &R, Vector_t &P, const double &deltat);
void Rot(double &px, double &pz, double &x, double &z, double xplane, double Norm_P,
double thetacou, double deltas, int coord);
void applyRotation(double &px,
double &pz,
double &x,
double &z,
double xplane,
double Norm_P,
double thetacou,
double deltas,
int coord);
void copyFromBunch(PartBunchBase<double, 3> *bunch,
const std::pair<Vector_t, double> &boundingSphere);
void addBackToBunch(PartBunchBase<double, 3> *bunch, unsigned i);
void applyNonDKS(PartBunchBase<double, 3> *bunch,
const std::pair<Vector_t, double> &boundingSphere,
size_t numParticlesInSimulation);
#ifdef OPAL_DKS
void copyFromBunchDKS(PartBunchBase<double, 3> *bunch,
void applyDKS(PartBunchBase<double, 3> *bunch,
const std::pair<Vector_t, double> &boundingSphere,
size_t numParticlesInSimulation);
void copyFromBunchDKS(PartBunchBase<double, 3> *bunch,
const std::pair<Vector_t, double> &boundingSphere);
void addBackToBunchDKS(PartBunchBase<double, 3> *bunch, unsigned i);
......@@ -118,52 +138,47 @@ private:
void deleteParticleFromLocalVectorDKS();
#endif
void deleteParticleFromLocalVector();
void calcStat(double Eng);
double T_m; // own time, maybe larger than in the bunch object
double dT_m; // dt from bunch
double T_m; // own time, maybe larger than in the bunch object
double dT_m; // dt from bunch
gsl_rng *rGen_m; // random number generator
gsl_rng *rGen_m; // random number generator
std::string material_m; // type of material e.g. aluminum
std::string FN_m; // name of element
ElementBase::ElementType collshape_m; // the type of element (DEGRADER, CCOLLIMATOR or FLEXIBLECOLLIMATOR)
std::string collshapeStr_m; // the type of element as string
std::string material_m; // type of material e.g. aluminum
std::unique_ptr<InsideTester> hitTester_m; // tests whether particles are inside material
ElementBase::ElementType collshape_m; // the type of element (DEGRADER, CCOLLIMATOR or FLEXIBLECOLLIMATOR)
std::string collshapeStr_m; // the type of element as string
// material parameters
double Z_m; // the atomic number
double A_m; // the atomic mass
double rho_m; // the volumetric mass density
double X0_m; // the radiation length
double I_m; // the mean excitation energy
double Z_m; // the atomic number [1]
double A_m; // the atomic mass [u]
double rho_m; // the volumetric mass density in [g cm^-3]
double X0_m; // the radiation length in [m]
double I_m; // the mean excitation energy [eV]
double A2_c; // coefficients to fit model to measurement data
double A3_c; // see e.g. page 16-20 in H. H. Andersen, J. F. Ziegler,
double A4_c; // "Hydrogen Stopping Powers and Ranges in All Elements",
double A5_c; // Pergamon Press, 1977
double A2_c; // coefficients to fit model to measurement data
double A3_c; // see e.g. page 16-20 in H. H. Andersen, J. F. Ziegler,
double A4_c; // "Hydrogen Stopping Powers and Ranges in All Elements",
double A5_c; // Pergamon Press, 1977
// number of particles that enter the material in current step (count for single step)
unsigned bunchToMatStat_m;
unsigned int bunchToMatStat_m;
// number of particles that are stopped by the material in current step
unsigned stoppedPartStat_m;
unsigned int stoppedPartStat_m;
// number of particles that leave the material in current step
unsigned rediffusedStat_m;
unsigned int rediffusedStat_m;
// number of local particles that are in the material
size_t locPartsInMat_m;
// some statistics
double Eavg_m; // average kinetic energy
double Emax_m; // maximum kinetic energy
double Emin_m; // minimum kinetic energy
double Eavg_m; // average kinetic energy
double Emax_m; // maximum kinetic energy
double Emin_m; // minimum kinetic energy
std::vector<PART> locParts_m; // local particles that are in material
std::vector<PART> locParts_m; // local particles that are in material
std::unique_ptr<LossDataSink> lossDs_m;
......@@ -185,8 +200,6 @@ private:
IpplTimings::TimerRef DegraderApplyTimer_m;
IpplTimings::TimerRef DegraderLoopTimer_m;
// :FIXME: remove unused declaration
// IpplTimings::TimerRef DegraderInitTimer_m;
IpplTimings::TimerRef DegraderDestroyTimer_m;
};
......@@ -200,8 +213,28 @@ void CollimatorPhysics::calcStat(double Eng) {
}
inline
void CollimatorPhysics::EnergyLoss(double &Eng, bool &pdead, const double &deltat) {
pdead = EnergyLoss(Eng, deltat);
double CollimatorPhysics::getTime() {
return T_m;
}
inline
size_t CollimatorPhysics::getParticlesInMat() {
return locPartsInMat_m;
}
inline
unsigned int CollimatorPhysics::getRediffused() {
return rediffusedStat_m;
}
inline
unsigned int CollimatorPhysics::getNumEntered() {
return bunchToMatStat_m;
}
inline
const std::string CollimatorPhysics::getType() const {
return "CollimatorPhysics";
}
#endif //COLLIMATORPHYSICS_HH
......
......@@ -31,6 +31,7 @@ public:
void updateElement(ElementBase *newref);
ElementBase* getElement();
virtual bool computeEnergyLoss(double &Eng, const double deltat, bool includeFluctuations = true) const = 0;
protected:
ElementBase *element_ref_m;
bool allParticleInMat_m; ///< if all particles are in matter stay inside the particle matter interaction
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment