17 const std::string& mode,
18 const std::string& binning)
19 : onebunch_m(
OpalData::getInstance()->getInputBasename() +
"-onebunch.h5")
20 , numBunch_m(numBunch)
22 , coeffDBunches_m(para)
23 , radiusLastTurn_m(0.0)
24 , radiusThisTurn_m(0.0)
42 *gmsg <<
"***---------------------------- MULTI-BUNCHES MULTI-ENERGY-BINS MODE "
43 <<
"----------------------------*** " <<
endl;
51 const int BinCount = 1;
53 const int MaxBinNum = 1000;
56 size_t partInBin[MaxBinNum] = {0};
67 *gmsg <<
"In this restart job, the multi-bunches mode is forcely set to AUTO mode." <<
endl;
70 *gmsg <<
"In this restart job, the multi-bunches mode is forcely set to FORCE mode." << endl
71 <<
"If the existing bunch number is less than the specified number of TURN, "
72 <<
"readin the phase space of STEP#0 from h5 file consecutively" <<
endl;
84 *gmsg <<
"* Store beam to H5 file for multibunch simulation ... ";
107 std::map<std::string, double> additionalAttributes = {
108 std::make_pair(
"REFPR", 0.0),
109 std::make_pair(
"REFPT", 0.0),
110 std::make_pair(
"REFPZ", 0.0),
111 std::make_pair(
"REFR", 0.0),
112 std::make_pair(
"REFTHETA", 0.0),
113 std::make_pair(
"REFZ", 0.0),
114 std::make_pair(
"AZIMUTH", 0.0),
115 std::make_pair(
"ELEVATION", 0.0),
116 std::make_pair(
"B-ref_x", 0.0),
117 std::make_pair(
"B-ref_z", 0.0),
118 std::make_pair(
"B-ref_y", 0.0),
119 std::make_pair(
"E-ref_x", 0.0),
120 std::make_pair(
"E-ref_z", 0.0),
121 std::make_pair(
"E-ref_y", 0.0),
122 std::make_pair(
"B-head_x", 0.0),
123 std::make_pair(
"B-head_z", 0.0),
124 std::make_pair(
"B-head_y", 0.0),
125 std::make_pair(
"E-head_x", 0.0),
126 std::make_pair(
"E-head_z", 0.0),
127 std::make_pair(
"E-head_y", 0.0),
128 std::make_pair(
"B-tail_x", 0.0),
129 std::make_pair(
"B-tail_z", 0.0),
130 std::make_pair(
"B-tail_y", 0.0),
131 std::make_pair(
"E-tail_x", 0.0),
132 std::make_pair(
"E-tail_z", 0.0),
133 std::make_pair(
"E-tail_y", 0.0)
138 h5wrapper.
writeStep(beam, additionalAttributes);
141 *gmsg <<
"Done." <<
endl;
152 *gmsg <<
"* Read beam from H5 file for multibunch simulation ... ";
169 if ( numParticles == 0 ) {
176 size_t firstParticle = numParticlesPerNode *
Ippl::myNode();
177 size_t lastParticle = firstParticle + numParticlesPerNode - 1;
179 lastParticle = numParticles - 1;
183 numParticles = lastParticle - firstParticle + 1;
188 std::unique_ptr<PartBunchBase<double, 3> > tmpBunch =
nullptr;
191 if ( amrbunch_p !=
nullptr ) {
198 tmpBunch->
create(numParticles);
200 h5wrapper.
readStep(tmpBunch.get(), firstParticle, lastParticle);
203 beam->
create(numParticles);
205 for(
size_t ii = 0; ii < numParticles; ++ ii, ++ localNum) {
206 beam->
R[localNum] = tmpBunch->R[ii];
207 beam->
P[localNum] = tmpBunch->P[ii];
208 beam->
M[localNum] = tmpBunch->M[ii];
209 beam->
Q[localNum] = tmpBunch->Q[ii];
211 beam->
Bin[localNum] = bunchNum;
212 beam->
bunchNum[localNum] = bunchNum;
219 *gmsg <<
"Done." <<
endl;
228 bool& flagTransition)
243 *gmsg <<
"* MBM: Checking for automatically injecting new bunch ..." <<
endl;
255 double XYrms = std::hypot(Rrms[0], Rrms[1]);
262 flagTransition =
true;
267 *gmsg <<
"* MBM: XYrms = " << XYrms <<
" [m]" <<
endl;
283 *gmsg <<
"* MBM: Injecting a new bunch ..." <<
endl;
299 "We shouldn't be here in single bunch mode.");
305 <<
" injected, total particle number = "
334 if ( mbmode.compare(
"FORCE") == 0 ) {
335 *gmsg <<
"FORCE mode: The multi bunches will be injected consecutively" <<
endl
336 <<
" after each revolution, until get \"TURNS\" bunches." <<
endl;
338 }
else if ( mbmode.compare(
"AUTO") == 0 ) {
339 *gmsg <<
"AUTO mode: The multi bunches will be injected only when the" <<
endl
340 <<
" distance between two neighboring bunches is below" <<
endl
341 <<
" the limitation. The control parameter is set to "
346 "MBMODE name \"" + mbmode +
"\" unknown.");
354 if ( binning.compare(
"BUNCH") == 0 ) {
355 *gmsg <<
"Use 'BUNCH' injection for binnning." <<
endl;
357 }
else if ( binning.compare(
"GAMMA") == 0 ) {
358 *gmsg <<
"Use 'GAMMA' for binning." <<
endl;
362 "MB_BINNING name \"" + binning +
"\" unknown.");
375 *gmsg <<
"Radial position at restart position = ";
377 *gmsg <<
"Initial radial position = ";
391 const unsigned long localNum = beam->
getLocalNum();
393 long int bunchTotalNum = 0;
394 unsigned long bunchLocalNum = 0;
404 std::vector<double> local(7 * dim + 1);
408 for(
unsigned long k = 0; k < localNum; ++k) {
409 if ( beam->
bunchNum[k] != bunchNr ) {
419 for (
unsigned int i = 0; i < dim; ++i) {
421 double r = beam->
R[k](i);
422 double p = beam->
P[k](i);
428 local[i + dim + 1] += p;
432 local[i + 2 * dim + 1] += r2;
435 local[i + 3 * dim + 1] += p * p;
438 local[i + 4 * dim + 1] += r * p;
441 local[i + 5 * dim + 1] += r2 * r;
444 local[i + 6 * dim + 1] += r2 * r2;
449 allreduce(bunchTotalNum, 1, std::plus<long int>());
453 throw OpalException(
"MultiBunchHandler::calcBunchBeamParameters()",
454 "Bunch number " + std::to_string(bunchNr) +
455 " exceeds bunch index " + std::to_string(beam->
getNumBunch() - 1));
460 if ( bunchTotalNum == 0 )
464 const double m0 = beam->
getM() * 1.0e-6;
465 local[0] -= bunchLocalNum;
468 allreduce(local.data(), local.size(), std::plus<double>());
470 double invN = 1.0 / double(bunchTotalNum);
471 binfo.
ekin = local[0] * invN;
476 for (
unsigned int i = 0; i < dim; ++i) {
478 double w = local[i + 1] * invN;
479 double pw = local[i + dim + 1] * invN;
480 double w2 = local[i + 2 * dim + 1] * invN;
481 double pw2 = local[i + 3 * dim + 1] * invN;
482 double wpw = local[i + 4 * dim + 1] * invN;
483 double w3 = local[i + 5 * dim + 1] * invN;
484 double w4 = local[i + 6 * dim + 1] * invN;
490 binfo.
prms[i] = pw2 - pw * pw;
491 if ( binfo.
prms[i] < 0 ) {
499 binfo.
rrms[i] = w2 - w * w;
502 binfo.
emit[i] = (binfo.
rrms[i] * binfo.
prms[i] - wpw * wpw);
509 - 3.0 * w * w * w * w;
510 binfo.
halo[i] = tmp / ( binfo.
rrms[i] * binfo.
rrms[i] );
519 double denom = 1.0 / (binfo.
rrms[i] * binfo.
prms[i]);
520 binfo.
correlation[i] = (std::isfinite(denom)) ? wpw * denom : 0.0;
523 double tmp = 1.0 /
std::pow(binfo.
ekin / m0 + 1.0, 2.0);
540 binfo_m[b].pathlength += lpaths[b];
void setBinning(std::string binning)
ParticleAttrib< Vector_t > P
The global OPAL structure.
virtual void create(size_t)
void updateTime(const double &dt)
void setLocalNumPerBunch(size_t numpart, short n)
T::PETE_Expr_t::PETE_Return_t max(const PETE_Expr< T > &expr, NDIndex< D > &loc)
short injectBunch(PartBunchBase< double, 3 > *beam, const PartData &ref, bool &flagTransition)
short numBunch_m
The number of bunches specified in TURNS of RUN commond.
The base class for all OPAL exceptions.
ParticleAttrib< double > Q
Vector_t get_rrms() const
std::string toUpper(const std::string &str)
MultiBunchHandler(PartBunchBase< double, 3 > *beam, const int &numBunch, const double &eta, const double ¶, const std::string &mode, const std::string &binning)
ParticleAttrib< short > bunchNum
bool calcBunchBeamParameters(PartBunchBase< double, 3 > *beam, short bunchNr)
bool readBunch(PartBunchBase< double, 3 > *beam, const PartData &ref)
void saveBunch(PartBunchBase< double, 3 > *beam)
void setTotalNumPerBunch(size_t numpart, short n)
size_t getTotalNum() const
double dot(const Vector3D &lhs, const Vector3D &rhs)
Vector dot product.
virtual void writeHeader()
pbase_t * getAmrParticleBase()
void updatePathLength(const std::vector< double > &lpaths)
static OpalData * getInstance()
virtual void writeStep(PartBunchBase< double, 3 > *, const std::map< std::string, double > &additionalStepAttributes)
void calcBeamParameters()
ParticleAttrib< short > PType
virtual void readStep(PartBunchBase< double, 3 > *, h5_ssize_t firstParticle, h5_ssize_t lastParticle)
void allreduce(const T *input, T *output, int count, Op op)
static void startTimer(TimerRef t)
void setRadiusTurns(const double &radius)
long unsigned int nParticles
Tps< T > pow(const Tps< T > &x, int y)
Integer power.
size_t getLocalNum() const
Vector_t get_centroid() const
Tps< T > sqrt(const Tps< T > &x)
Square root.
std::vector< beaminfo_t > binfo_m
short getNumBunch() const
bool resetPartBinID2(const double eta)
reset Bin[] for each particle according to the method given in paper PAST-AB(064402) by G...
beaminfo_t & getBunchInfo(short bunchNr)
void setNumBunch(short n)
ParticleAttrib< int > Bin
ParticleAttrib< double > M
Timing::TimerRef TimerRef
static TimerRef getTimer(const char *nm)
static void stopTimer(TimerRef t)
size_t getNumParticles() const
static Communicate * Comm
void setMode(const std::string &mbmode)
set the working sub-mode for multi-bunch mode: "FORCE" or "AUTO"
void setPBins(PartBins *pbin)
Inform & endl(Inform &inf)
void updateParticleBins(PartBunchBase< double, 3 > *beam)