OPAL (Object Oriented Parallel Accelerator Library)  2024.1
OPAL
DumpEMFields.cpp
Go to the documentation of this file.
1 //
2 // Class DumpEMFields
3 // DumpEMFields dumps the dynamically changing fields of a Ring in a user-
4 // defined grid.
5 //
6 // Copyright (c) 2017, Chris Rogers
7 // All rights reserved
8 //
9 // This file is part of OPAL.
10 //
11 // OPAL is free software: you can redistribute it and/or modify
12 // it under the terms of the GNU General Public License as published by
13 // the Free Software Foundation, either version 3 of the License, or
14 // (at your option) any later version.
15 //
16 // You should have received a copy of the GNU General Public License
17 // along with OPAL. If not, see <https://www.gnu.org/licenses/>.
18 //
20 
22 #include "AbsBeamline/Component.h"
23 #include "Attributes/Attributes.h"
25 #include "Physics/Units.h"
27 #include "Utilities/Util.h"
28 
29 #include <cmath>
30 #include <filesystem>
31 #include <fstream>
32 #include <map>
33 
34 extern Inform* gmsg;
35 
36 std::unordered_set<DumpEMFields*> DumpEMFields::dumpsSet_m;
37 
39  Action(SIZE, "DUMPEMFIELDS",
40  "The \"DUMPEMFIELDS\" statement dumps a field map to a user-defined "
41  "field file, for checking that fields are generated correctly. "
42  "The fields are written out on a grid in space and time."),
43  grid_m(nullptr),
44  filename_m("") {
45 
46  // would be nice if "steps" could be integer
48  ("FILE_NAME", "Name of the file to which field data is dumped");
49 
51  ("COORDINATE_SYSTEM", "Choose to use CARTESIAN or CYLINDRICAL coordinates", {"CARTESIAN", "CYLINDRICAL"}, "CARTESIAN");
52 
54  ("X_START", "(Cartesian) Start point in the grid in x [m]");
55 
56  itsAttr[DX] = Attributes::makeReal
57  ("DX", "(Cartesian) Grid step size in x [m]");
58 
59  itsAttr[X_STEPS] = Attributes::makeReal
60  ("X_STEPS", "(Cartesian) Number of steps in x");
61 
62  itsAttr[Y_START] = Attributes::makeReal
63  ("Y_START", "(Cartesian) Start point in the grid in y [m]");
64 
65  itsAttr[DY] = Attributes::makeReal
66  ("DY", "(Cartesian) Grid step size in y [m]");
67 
68  itsAttr[Y_STEPS] = Attributes::makeReal
69  ("Y_STEPS", "(Cartesian) Number of steps in y");
70 
71  itsAttr[Z_START] = Attributes::makeReal
72  ("Z_START", "Start point in the grid in z [m]");
73 
74  itsAttr[DZ] = Attributes::makeReal
75  ("DZ", "Grid step size in z [m]");
76 
77  itsAttr[Z_STEPS] = Attributes::makeReal
78  ("Z_STEPS", "Number of steps in z");
79 
80  itsAttr[T_START] = Attributes::makeReal
81  ("T_START", "Start point in the grid in time [ns]");
82 
83  itsAttr[DT] = Attributes::makeReal
84  ("DT", "Grid step size in time [ns]");
85 
86  itsAttr[T_STEPS] = Attributes::makeReal
87  ("T_STEPS", "Number of steps in time");
88 
89  itsAttr[R_START] = Attributes::makeReal
90  ("R_START", "(Cylindrical) Start point in the grid in radius [m]");
91 
92  itsAttr[DR] = Attributes::makeReal
93  ("DR", "(Cylindrical) Grid step size in radius [m]");
94 
95  itsAttr[R_STEPS] = Attributes::makeReal
96  ("R_STEPS", "(Cylindrical) Number of steps in radius");
97 
99  ("PHI_START", "(Cylindrical) Start point in the grid in phi [rad]");
100 
101  itsAttr[DPHI] = Attributes::makeReal
102  ("DPHI", "(Cylindrical) Grid step size in phi [rad]");
103 
105  ("PHI_STEPS", "(Cylindrical) Number of steps in phi");
106 
108 }
109 
110 DumpEMFields::DumpEMFields(const std::string& name, DumpEMFields* parent):
111  Action(name, parent), grid_m(nullptr)
112 {}
113 
115  delete grid_m;
116  dumpsSet_m.erase(this);
117 }
118 
119 DumpEMFields* DumpEMFields::clone(const std::string& name) {
120  DumpEMFields* dumper = new DumpEMFields(name, this);
121  if (grid_m != nullptr) {
122  dumper->grid_m = grid_m->clone();
123  }
124  dumper->filename_m = filename_m;
125  dumper->coordinates_m = coordinates_m;
126  if (dumpsSet_m.find(this) != dumpsSet_m.end()) {
127  dumpsSet_m.insert(dumper);
128  }
129  return dumper;
130 }
131 
133  static const std::map<std::string, CoordinateSystem> stringCoordinateSystem_s = {
134  {"CARTESIAN", CoordinateSystem::CARTESIAN},
135  {"CYLINDRICAL", CoordinateSystem::CYLINDRICAL}
136  };
137  coordinates_m = stringCoordinateSystem_s.at(Attributes::getString(itsAttr[COORDINATE_SYSTEM]));
138 }
139 
141  buildGrid();
142  // the routine for action (OpalParser/OpalParser) calls execute and then
143  // deletes 'this'; so we must build a copy that lasts until the field maps
144  // are constructed and we are ready for tracking (which is when the field
145  // maps are written). Hence the clone call below.
146  dumpsSet_m.insert(this->clone(""));
147 }
148 
150  std::vector<double> spacing(4);
151  std::vector<double> origin(4);
152  std::vector<int> gridSize(4);
154 
155  switch (coordinates_m) {
157  origin[0] = Attributes::getReal(itsAttr[X_START]);
158  spacing[0] = Attributes::getReal(itsAttr[DX]);
159  double nx = Attributes::getReal(itsAttr[X_STEPS]);
160  Util::checkInt(nx, "X_STEPS");
161  gridSize[0] = nx;
162 
163  origin[1] = Attributes::getReal(itsAttr[Y_START]);
164  spacing[1] = Attributes::getReal(itsAttr[DY]);
165  double ny = Attributes::getReal(itsAttr[Y_STEPS]);
166  Util::checkInt(ny, "Y_STEPS");
167  gridSize[1] = ny;
168 
169  break;
170  }
172  origin[0] = Attributes::getReal(itsAttr[R_START]);
173  spacing[0] = Attributes::getReal(itsAttr[DR]);
175  Util::checkInt(nr, "R_STEPS");
176  gridSize[0] = nr;
177 
178  origin[1] = Attributes::getReal(itsAttr[PHI_START]);
179  spacing[1] = Attributes::getReal(itsAttr[DPHI]);
180  double nphi = Attributes::getReal(itsAttr[PHI_STEPS]);
181  Util::checkInt(nphi, "PHI_STEPS");
182  gridSize[1] = nphi;
183 
184  break;
185  }
186  }
187 
188  origin[2] = Attributes::getReal(itsAttr[Z_START]);
189  spacing[2] = Attributes::getReal(itsAttr[DZ]);
190  double nz = Attributes::getReal(itsAttr[Z_STEPS]);
191  Util::checkInt(nz, "Z_STEPS");
192  gridSize[2] = nz;
193 
194  origin[3] = Attributes::getReal(itsAttr[T_START]);
195  spacing[3] = Attributes::getReal(itsAttr[DT]);
196  double nt = Attributes::getReal(itsAttr[T_STEPS]);
197  Util::checkInt(nt, "T_STEPS");
198  gridSize[3] = nt;
199 
200  if (grid_m != nullptr) {
201  delete grid_m;
202  }
203 
204  grid_m = new interpolation::NDGrid(4, &gridSize[0], &spacing[0], &origin[0]);
205 
207 }
208 
211  for (dump_iter it = dumpsSet_m.begin(); it != dumpsSet_m.end(); ++it) {
212  (*it)->writeFieldThis(field);
213  }
214 }
215 
216 void DumpEMFields::writeHeader(std::ofstream& fout) const {
217  fout << grid_m->end().toInteger() << "\n";
218  switch (coordinates_m) {
220  fout << 1 << " x [m]\n";
221  fout << 2 << " y [m]\n";
222  fout << 3 << " z [m]\n";
223  fout << 4 << " t [ns]\n";
224  fout << 5 << " Bx [kGauss]\n";
225  fout << 6 << " By [kGauss]\n";
226  fout << 7 << " Bz [kGauss]\n";
227  fout << 8 << " Ex [MV/m]\n";
228  fout << 9 << " Ey [MV/m]\n";
229  fout << 10 << " Ez [MV/m]\n";
230  break;
231  }
233  fout << 1 << " r [m]\n";
234  fout << 2 << " phi [deg]\n";
235  fout << 3 << " z [m]\n";
236  fout << 4 << " t [ns]\n";
237  fout << 5 << " Br [kGauss]\n";
238  fout << 6 << " Bphi [kGauss]\n";
239  fout << 7 << " Bz [kGauss]\n";
240  fout << 8 << " Er [MV/m]\n";
241  fout << 9 << " Ephi [MV/m]\n";
242  fout << 10 << " Ez [MV/m]\n";
243  break;
244  }
245  }
246  fout << 0 << std::endl;
247 }
248 
250  const Vector_t& pointIn,
251  const double& time,
252  std::ofstream& fout) const {
253  Vector_t centroid(0., 0., 0.);
254  Vector_t E(0., 0., 0.);
255  Vector_t B(0., 0., 0.);
256  Vector_t point = pointIn;
258  // pointIn is r, phi, z
259  point[0] = std::cos(pointIn[1])*pointIn[0];
260  point[1] = std::sin(pointIn[1])*pointIn[0];
261  }
262 
263  field->apply(point, centroid, time, E, B);
264  Vector_t Bout = B;
265  Vector_t Eout = E;
267  // pointIn is r, phi, z
268  Bout[0] = B[0]*std::cos(pointIn[1])+B[1]*std::sin(pointIn[1]);
269  Bout[1] = -B[0]*std::sin(pointIn[1])+B[1]*std::cos(pointIn[1]);
270  Eout[0] = E[0]*std::cos(pointIn[1])+E[1]*std::sin(pointIn[1]);
271  Eout[1] = -E[0]*std::sin(pointIn[1])+E[1]*std::cos(pointIn[1]);
272  fout << pointIn[0] << " " << pointIn[1]*Units::rad2deg << " " << pointIn[2] << " " << time << " ";
273  } else {
274  fout << pointIn[0] << " " << pointIn[1] << " " << pointIn[2] << " " << time << " ";
275  }
276 
277  fout << Bout[0] << " " << Bout[1] << " " << Bout[2] << " ";
278  fout << Eout[0] << " " << Eout[1] << " " << Eout[2] << "\n";
279 }
280 
282  if (grid_m == nullptr) {
283  throw OpalException("DumpEMFields::writeFieldThis",
284  "The grid was nullptr; there was a problem with the DumpEMFields initialisation.");
285  }
286  if (field == nullptr) {
287  throw OpalException("DumpEMFields::writeFieldThis",
288  "The field to be written was nullptr.");
289  }
290 
291  *gmsg << *this << endl;
292 
293  std::string fname;
294  if (std::filesystem::path(filename_m).is_absolute() == true) {
295  fname = filename_m;
296  } else {
297  fname = Util::combineFilePath({
299  filename_m
300  });
301  }
302 
303  std::vector<double> point_std(4);
304  Vector_t point(0., 0., 0.);
305  std::ofstream fout;
306  try {
307  fout.open(fname.c_str(), std::ofstream::out);
308  } catch (std::exception& exc) {
309  throw OpalException("DumpEMFields::writeFieldThis",
310  "Failed to open DumpEMFields file " + filename_m);
311  }
312  if (!fout.good()) {
313  throw OpalException("DumpEMFields::writeFieldThis",
314  "Failed to open DumpEMFields file " + filename_m);
315  }
316  // set precision
317  writeHeader(fout);
319  it < grid_m->end();
320  ++it) {
321  it.getPosition(&point_std[0]);
322  for (size_t i = 0; i < 3; ++i) {
323  point[i] = point_std[i];
324  }
325  double time = point_std[3];
326  writeFieldLine(field, point, time, fout);
327  }
328  if (!fout.good()) {
329  throw OpalException("DumpEMFields::writeFieldThis",
330  "Something went wrong during writing " + filename_m);
331  }
332  fout.close();
333 }
334 
335 void DumpEMFields::print(std::ostream& os) const {
336  os << "* ************* D U M P E M F I E L D S ****************************************** " << std::endl;
337  os << "* File name: '" << filename_m << "'\n";
339  os << "* Coordinate system: " << Attributes::getString(itsAttr[COORDINATE_SYSTEM]) << '\n'
340  << "* X_START = " << Attributes::getReal(itsAttr[X_START]) << " [m]\n"
341  << "* DX = " << Attributes::getReal(itsAttr[DX]) << " [m]\n"
342  << "* X_STEPS = " << Attributes::getReal(itsAttr[X_STEPS]) << '\n'
343  << "* Y_START = " << Attributes::getReal(itsAttr[Y_START]) << " [m]\n"
344  << "* DY = " << Attributes::getReal(itsAttr[DY]) << " [m]\n"
345  << "* Y_STEPS = " << Attributes::getReal(itsAttr[Y_STEPS]) << '\n';
347  os << "* Coordinate system: " << Attributes::getString(itsAttr[COORDINATE_SYSTEM]) << '\n'
348  << "* R_START = " << Attributes::getReal(itsAttr[R_START]) << " [m]\n"
349  << "* DR = " << Attributes::getReal(itsAttr[DR]) << " [m]\n"
350  << "* R_STEPS = " << Attributes::getReal(itsAttr[R_STEPS]) << '\n'
351  << "* PHI_START = " << Attributes::getReal(itsAttr[PHI_START]) << " [rad]\n"
352  << "* DPHI = " << Attributes::getReal(itsAttr[DPHI]) << " [rad]\n"
353  << "* PHI_STEPS = " << Attributes::getReal(itsAttr[PHI_STEPS]) << '\n';
354  }
355  os << "* Z_START = " << Attributes::getReal(itsAttr[Z_START]) << " [m]\n"
356  << "* DZ = " << Attributes::getReal(itsAttr[DZ]) << " [m]\n"
357  << "* Z_STEPS = " << Attributes::getReal(itsAttr[Z_STEPS]) << '\n'
358  << "* T_START = " << Attributes::getReal(itsAttr[T_START]) << " [ns]\n"
359  << "* DT = " << Attributes::getReal(itsAttr[DT]) << " [ns]\n"
360  << "* T_STEPS = " << Attributes::getReal(itsAttr[T_STEPS]) << '\n';
361  os << "* ********************************************************************************** " << std::endl;
362 }
Attribute makeReal(const std::string &name, const std::string &help)
Make real attribute.
Definition: Attributes.cpp:240
static OpalData * getInstance()
Definition: OpalData.cpp:196
void writeHeader(std::ofstream &fout) const
virtual void buildGrid()
and that you know you can do these things To protect your we need to make restrictions that forbid anyone to deny you these rights or to ask you to surrender the rights These restrictions translate to certain responsibilities for you if you distribute copies of the or if you modify it For if you distribute copies of such a whether gratis or for a you must give the recipients all the rights that you have You must make sure that receive or can get the source code And you must show them these terms so they know their rights We protect your rights with two distribute and or modify the software for each author s protection and we want to make certain that everyone understands that there is no warranty for this free software If the software is modified by someone else and passed we want its recipients to know that what they have is not the so that any problems introduced by others will not reflect on the original authors reputations any free program is threatened constantly by software patents We wish to avoid the danger that redistributors of a free program will individually obtain patent in effect making the program proprietary To prevent we have made it clear that any patent must be licensed for everyone s free use or not licensed at all The precise terms and conditions for distribution and modification follow GNU GENERAL PUBLIC LICENSE TERMS AND CONDITIONS FOR DISTRIBUTION AND MODIFICATION This License applies to any program or other work which contains a notice placed by the copyright holder saying it may be distributed under the terms of this General Public License The refers to any such program or and a work based on the Program means either the Program or any derivative work under copyright a work containing the Program or a portion of it
Definition: LICENSE:43
std::string getString(const Attribute &attr)
Get string value.
Definition: Attributes.cpp:343
void writeFieldLine(Component *field, const Vector_t &point, const double &time, std::ofstream &fout) const
virtual void writeFieldThis(Component *field)
Mesh::Iterator end() const
Definition: NDGrid.h:362
void print(std::ostream &os) const
virtual void execute()
Inform & endl(Inform &inf)
Definition: Inform.cpp:42
std::string::iterator iterator
Definition: MSLang.h:15
virtual bool apply(const size_t &i, const double &t, Vector_t &E, Vector_t &B)
Definition: Component.cpp:99
std::string filename_m
Definition: DumpEMFields.h:143
void registerOwnership(const AttributeHandler::OwnerType &itsClass) const
Definition: Object.cpp:191
The base class for all OPAL exceptions.
Definition: OpalException.h:28
std::string getAuxiliaryOutputDirectory() const
get the name of the the additional data directory
Definition: OpalData.cpp:666
Attribute makePredefinedString(const std::string &name, const std::string &help, const std::initializer_list< std::string > &predefinedStrings)
Make predefined string attribute.
Definition: Attributes.cpp:409
Definition: Inform.h:42
Attribute makeString(const std::string &name, const std::string &help)
Make string attribute.
Definition: Attributes.cpp:332
constexpr double rad2deg
Definition: Units.h:146
static void writeFields(Component *field)
std::vector< Attribute > itsAttr
The object attributes.
Definition: Object.h:216
Tps< T > cos(const Tps< T > &x)
Cosine.
Definition: TpsMath.h:129
const std::string name
c Accompany it with the information you received as to the offer to distribute corresponding source complete source code means all the source code for all modules it plus any associated interface definition plus the scripts used to control compilation and installation of the executable as a special exception
Definition: LICENSE:157
NDGrid * clone()
Definition: NDGrid.h:350
std::string combineFilePath(std::initializer_list< std::string > ilist)
Definition: Util.cpp:197
double getReal(const Attribute &attr)
Return real value.
Definition: Attributes.cpp:252
virtual ~DumpEMFields()
virtual DumpEMFields * clone(const std::string &name)
CoordinateSystem coordinates_m
Definition: DumpEMFields.h:145
Interface for a single beam element.
Definition: Component.h:50
void parseCoordinateSystem()
Mesh::Iterator begin() const
Definition: NDGrid.h:358
const int nr
Definition: ClassicRandom.h:24
Tps< T > sin(const Tps< T > &x)
Sine.
Definition: TpsMath.h:111
void checkInt(double real, std::string name, double tolerance)
Definition: Util.cpp:205
The base class for all OPAL actions.
Definition: Action.h:30
interpolation::NDGrid * grid_m
Definition: DumpEMFields.h:142
Inform * gmsg
Definition: Main.cpp:70
static std::unordered_set< DumpEMFields * > dumpsSet_m
Definition: DumpEMFields.h:147