Commit 2930269f authored by Daniel Winklehner's avatar Daniel Winklehner
Browse files

*) Fixed a problem with saving the lost particles if not all processors have...

*) Fixed a problem with saving the lost particles if not all processors have lost particles. \n*) Consolidated Tracker_RK4 and Tracker_LF in ParallelCyclotronTracker. Enclosed in #ifdef for now. Activate by uncommenting #define GENERICTRACKER in ParallelCyclotronTracker.h. \n*) Added the general OPAL OPTION boundpdestroyfreq for Cyclotron Tracker.
parent e8f7f39b
......@@ -281,6 +281,7 @@ bool Cyclotron::apply(const size_t &id, const double &t, double E[], double B[])
bool Cyclotron::apply(const size_t &id, const double &t, Vector_t &E, Vector_t &B) {
bool flagNeedUpdate = false;
const double rpos = sqrt(RefPartBunch_m->R[id](0) * RefPartBunch_m->R[id](0)
+ RefPartBunch_m->R[id](1) * RefPartBunch_m->R[id](1));
const double zpos = RefPartBunch_m->R[id](2);
......@@ -288,14 +289,14 @@ bool Cyclotron::apply(const size_t &id, const double &t, Vector_t &E, Vector_t &
if (zpos > maxz_m || zpos < minz_m || rpos > maxr_m || rpos < minr_m){
flagNeedUpdate = true;
Inform gmsgALL("OPAL ", INFORM_ALL_NODES);
gmsgALL<<getName() << ": particle "<< id <<" out of the global aperture of cyclotron!"<< endl;
gmsgALL << getName() << ": particle "<< id <<" out of the global aperture of cyclotron!"<< endl;
} else{
flagNeedUpdate = apply(RefPartBunch_m->R[id], Vector_t(0.0), t, E, B);
if(flagNeedUpdate){
Inform gmsgALL("OPAL ", INFORM_ALL_NODES);
gmsgALL<<getName() << ": particle "<< id <<" out of the field map boundary!"<< endl;
gmsgALL << getName() << ": particle "<< id <<" out of the field map boundary!"<< endl;
}
}
......@@ -303,6 +304,7 @@ bool Cyclotron::apply(const size_t &id, const double &t, Vector_t &E, Vector_t &
lossDs_m->addParticle(RefPartBunch_m->R[id], RefPartBunch_m->P[id],id);
RefPartBunch_m->Bin[id] = -1;
}
return flagNeedUpdate;
}
......@@ -499,9 +501,11 @@ bool Cyclotron::apply(const Vector_t &R, const Vector_t &centroid, const double
}
void Cyclotron::finalise() {
online_m = false;
lossDs_m->save();
*gmsg << "* Finalize cyclotron" << endl;
}
bool Cyclotron::bends() const {
......
......@@ -1117,8 +1117,7 @@ void PartBunch::computeSelfFields_cycl(double gamma) {
* Comments -DW:
* I have made some changes in here:
* -) Some refacturing to make more similar to computeSelfFields()
* -) Added meanR and quaternion to be handed to the function (TODO: fall back to meanR = 0 and unit quaternion
* if not specified) so that SAAMG solver knows how to rotate the boundary geometry correctly.
* -) Added meanR and quaternion to be handed to the function so that SAAMG solver knows how to rotate the boundary geometry correctly.
* -) Fixed an error where gamma was not taken into account correctly in direction of movement (y in cyclotron)
* -) Comment: There is no account for image charges in the cyclotron tracker (yet?)!
*/
......@@ -1186,7 +1185,7 @@ void PartBunch::computeSelfFields_cycl(double gamma,
//do the multiplication of the grid-cube volume coming
//from the discretization of the convolution integral.
//this is only necessary for the FFT solver
//FIXME: later move this scaling into FFTPoissonSolver
//TODO FIXME: later move this scaling into FFTPoissonSolver
if(fs_m->getFieldSolverType() == "FFT" || fs_m->getFieldSolverType() == "FFTBOX")
rho_m *= hr_scaled[0] * hr_scaled[1] * hr_scaled[2];
......
This diff is collapsed.
#ifndef OPAL_ParallelCyclotronTracker_HH
#define OPAL_ParallelCyclotronTracker_HH
// #define GENERICTRACKER
// ------------------------------------------------------------------------
// $RCSfile: ParallelCyclotronTracker.h,v $
// ------------------------------------------------------------------------
......@@ -252,6 +252,11 @@ private:
void Tracker_RK4();
void Tracker_MTS();
#ifdef GENERICTRACKER
void Tracker_Generic();
bool getFieldsAtPoint(const double &t, const size_t &Pindex, Vector_t &Efield, Vector_t &Bfield);
#endif
/*
Local Variables both used by the integration methods
*/
......@@ -418,6 +423,7 @@ private:
void initTrackOrbitFile();
void singleParticleDump();
void singleParticleDump_rk4(); // TODO: Consolidate the two -DW
void bunchDumpPhaseSpaceStatData();
......
......@@ -110,6 +110,10 @@ namespace Options {
// during each time step. If false then the time step during emission is set so that one
// energy bin of the beam is emitted during each time step.
bool fineEmission = true;
// Governs how often boundp_destroy is called to destroy lost particles
// Mainly used in the CyclotronTracker as of now -DW
int boundpDestroyFreq = 10;
}
......@@ -152,6 +156,7 @@ namespace {
FINEEMISSION,
ENABLEHDF5,
ASCIIDUMP,
BOUNDPDESTROYFREQ,
SIZE
};
}
......@@ -247,6 +252,9 @@ Option::Option():
itsAttr[ASCIIDUMP] = Attributes::makeBool
("ASCIIDUMP", "If true, some of the elements dump in ASCII instead of HDF5", false);
itsAttr[BOUNDPDESTROYFREQ] = Attributes::makeReal
("BOUNDPDESTROYFREQ", "The frequency to do boundp_destroy to delete lost particles. Default value is 10.");
FileStream::setEcho(echo);
rangen.init55(seed);
}
......@@ -289,6 +297,7 @@ Option::Option(const string &name, Option *parent):
Attributes::setReal(itsAttr[NLHS], nLHS);
Attributes::setBool(itsAttr[ENABLEHDF5], enableHDF5);
Attributes::setBool(itsAttr[ASCIIDUMP], asciidump);
Attributes::setReal(itsAttr[BOUNDPDESTROYFREQ], boundpDestroyFreq);
}
......@@ -411,6 +420,10 @@ void Option::execute() {
fineEmission = bool(Attributes::getBool(itsAttr[FINEEMISSION]));
}
if(itsAttr[BOUNDPDESTROYFREQ]) {
boundpDestroyFreq = int(Attributes::getReal(itsAttr[BOUNDPDESTROYFREQ]));
}
// Set message flags.
FileStream::setEcho(echo);
......
......@@ -141,11 +141,13 @@ void LossDataSink::addParticle(const Vector_t x, const Vector_t p, const size_t
void LossDataSink::save() {
if (element_m == std::string("")) return;
if (hasNoParticlesToDump()) return;
if (hasNoParticlesToDump()) return;
if (h5hut_mode_m) {
if (!Options::enableHDF5) return;
fn_m = element_m + std::string(".h5");
INFOMSG("Save " << fn_m << endl);
openH5();
......@@ -153,19 +155,23 @@ void LossDataSink::save() {
saveH5();
H5CloseFile(H5file_m);
Ippl::Comm->barrier();
}
else {
fn_m = element_m + std::string(".loss");
INFOMSG("Save " << fn_m << endl);
if(OpalData::getInstance()->inRestartRun())
append();
else
open();
writeHeaderASCII();
saveASCII();
close();
Ippl::Comm->barrier();
fn_m = element_m + std::string(".loss");
INFOMSG("Save " << fn_m << endl);
if(OpalData::getInstance()->inRestartRun())
append();
else
open();
writeHeaderASCII();
saveASCII();
close();
Ippl::Comm->barrier();
}
x_m.clear();
y_m.clear();
z_m.clear();
......@@ -177,13 +183,27 @@ void LossDataSink::save() {
time_m.clear();
}
// Note: This was changed to calculate the global number of dumped particles
// because there are two cases to be considered:
// 1. ALL nodes have 0 lost particles -> nothing to be done.
// 2. Some nodes have 0 lost particles, some not -> H5 can handle that but all
// nodes HAVE to participate, otherwise H5 waits endlessly for a response from
// the nodes that didn't enter the saveH5 function. -DW
bool LossDataSink::hasNoParticlesToDump() {
size_t nLoc = x_m.size();
reduce(nLoc, nLoc, OpAddAssign());
return nLoc == 0;
}
void LossDataSink::saveH5() {
h5_int64_t rc;
size_t nLoc = x_m.size();
INFOMSG("saveH5 n= " << nLoc << endl);
INFOMSG("saveH5 n = " << nLoc << endl);
std::unique_ptr<char[]> varray(new char[(nLoc)*sizeof(double)]);
double *farray = reinterpret_cast<double *>(varray.get());
......
......@@ -65,15 +65,18 @@ private:
}
}
bool hasNoParticlesToDump() {
return x_m.size() == 0;
}
//bool hasNoParticlesToDump() {
//return x_m.size() == 0;
//}
bool hasNoParticlesToDump();
void writeHeaderH5();
void openH5();
void saveH5();
void saveASCII();
private:
......
......@@ -39,7 +39,6 @@ namespace Options {
// true if in bet mode
extern bool bet;
/// Trace flag.
// If true, print CPU time before and after each command.
extern bool mtrace;
......@@ -81,6 +80,10 @@ namespace Options {
// if true, dump phase space after each turn
extern bool psDumpEachTurn;
// Governs how often boundp_destroy is called to destroy lost particles
// Mainly used in the CyclotronTracker as of now -DW
extern int boundpDestroyFreq;
/// flag to decide in which coordinate frame the phase space will be dumped for OPAL-cycl
// if true, in local Cartesian frame, otherwise in global Cartesian frame
extern bool psDumpLocalFrame;
......@@ -95,7 +98,6 @@ namespace Options {
// the particle will be deleted artifically to hold the accuracy of space charge calculation. The default setting of -1 stands for no deletion.
extern int remotePartDel;
/// this allows to repeat tracks starting always at the begining of the lattice and
/// generates a new distribution
......
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