Cyclotron.h 8.96 KB
Newer Older
gsell's avatar
gsell committed
1
//
2
// Class Cyclotron
3
//   Defines the abstract interface for a cyclotron.
gsell's avatar
gsell committed
4
//
5 6 7 8 9 10 11 12 13 14
// Copyright (c) 2007 - 2012, Jianjun Yang and Andreas Adelmann, Paul Scherrer Institut, Villigen PSI, Switzerland
// Copyright (c) 2013 - 2020, Paul Scherrer Institut, Villigen PSI, Switzerland
// All rights reserved
//
// Implemented as part of the PhD thesis
// "Beam dynamics in high intensity cyclotrons including neighboring bunch effects"
// and the paper
// "Beam dynamics in high intensity cyclotrons including neighboring bunch effects:
// Model, implementation, and application"
// (https://journals.aps.org/prab/pdf/10.1103/PhysRevSTAB.13.064201)
gsell's avatar
gsell committed
15
//
16
// This file is part of OPAL.
gsell's avatar
gsell committed
17
//
18 19 20 21 22 23 24 25 26 27
// OPAL is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// You should have received a copy of the GNU General Public License
// along with OPAL. If not, see <https://www.gnu.org/licenses/>.
//
#ifndef CLASSIC_Cyclotron_HH
#define CLASSIC_Cyclotron_HH
gsell's avatar
gsell committed
28 29

#include "AbsBeamline/Component.h"
snuverink_j's avatar
snuverink_j committed
30 31 32

#include <string>
#include <vector>
gsell's avatar
gsell committed
33 34

class Fieldmap;
35
class LossDataSink;
36
class TrimCoil;
gsell's avatar
gsell committed
37 38 39 40 41


struct BfieldData {
    std::string filename;
    // known from file: field and three theta derivatives
42 43 44 45
    std::vector<double> bfld;   //Bz
    std::vector<double> dbt;    //dBz/dtheta
    std::vector<double> dbtt;   //d2Bz/dtheta2
    std::vector<double> dbttt;  //d3Bz/dtheta3
gsell's avatar
gsell committed
46 47

    // to be calculated in getdiffs: all other derivatives:
48 49 50
    std::vector<double> dbr;    // dBz/dr
    std::vector<double> dbrr;   // ...
    std::vector<double> dbrrr;
gsell's avatar
gsell committed
51

52 53 54
    std::vector<double> dbrt;
    std::vector<double> dbrrt;
    std::vector<double> dbrtt;
gsell's avatar
gsell committed
55 56

    // used to get (Br,Btheta,Bz) at any off-plane point
57 58 59
    std::vector<double> f2;  // for Bz
    std::vector<double> f3;  // for Br
    std::vector<double> g3;  // for Btheta
gsell's avatar
gsell committed
60 61

    // Grid-Size
ext-calvo_p's avatar
ext-calvo_p committed
62 63 64
    int nrad, ntet; //need to be read from inputfile.
    int ntetS;      // one more grid line is stored in azimuthal direction
    int ntot;       // total grid points number.
gsell's avatar
gsell committed
65 66 67 68 69 70

    // Mean and Maximas
    double bacc, dbtmx, dbttmx, dbtttmx;
};

struct BPositions {
71
    // these 4 parameters are need to be read from field file.
gsell's avatar
gsell committed
72 73 74 75
    double  rmin, delr;
    double  tetmin, dtet;

    // Radii and step width of initial Grid
76
    std::vector<double> rarr;
gsell's avatar
gsell committed
77 78 79 80 81 82 83

    //  int     ThetaPeriodicity; // Periodicity of Magnetic field
    double  Bfact;      // MULTIPLICATION FACTOR FOR MAGNETIC FIELD
};


class Cyclotron: public Component {
ext-calvo_p's avatar
ext-calvo_p committed
84

gsell's avatar
gsell committed
85 86
public:

ext-calvo_p's avatar
ext-calvo_p committed
87 88 89 90 91 92 93 94 95 96
    enum class BFieldType {
        PSIBF,
        CARBONBF,
        ANSYSBF,
        AVFEQBF,
        FFABF,
        BANDRF,
        SYNCHRO
    };

gsell's avatar
gsell committed
97
    /// Constructor with given name.
ext-calvo_p's avatar
ext-calvo_p committed
98
    explicit Cyclotron(const std::string& name);
gsell's avatar
gsell committed
99 100 101

    Cyclotron();
    Cyclotron(const Cyclotron &);
102

gsell's avatar
gsell committed
103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119
    virtual ~Cyclotron();

    /// Apply visitor to Cyclotron.
    virtual void accept(BeamlineVisitor &) const;

    /// Get number of slices.
    //  Slices and stepsize used to determine integration step.
    virtual double getSlices() const = 0;

    /// Get stepsize.
    //  Slices and stepsize used to determine integration step.
    virtual double getStepsize() const = 0;

    void setFieldMapFN(std::string fmapfn);
    virtual std::string getFieldMapFN() const;

    void setRfFieldMapFN(std::vector<std::string> rffmapfn);
120
    void setRFFCoeffFN(std::vector<std::string> rff_coeff_fn);
121
    void setRFVCoeffFN(std::vector<std::string> rfv_coeff_fn);
122
    
ext-calvo_p's avatar
ext-calvo_p committed
123
    void setCyclotronType(std::string t);
124 125
    const std::string &getCyclotronType() const;
    virtual ElementBase::ElementType getType() const;
gsell's avatar
gsell committed
126

ext-calvo_p's avatar
ext-calvo_p committed
127
    virtual void getDimensions(double& zBegin, double& zEnd) const;
gsell's avatar
gsell committed
128

129 130
    unsigned int getNumberOfTrimcoils() const;

gsell's avatar
gsell committed
131 132 133 134
    void setCyclHarm(double h);
    virtual double getCyclHarm() const;

    void setRfPhi(std::vector<double> f);
135
    double getRfPhi(unsigned int i) const;
gsell's avatar
gsell committed
136 137

    void setRfFrequ(std::vector<double> f);
138
    double getRfFrequ(unsigned int i) const;
gsell's avatar
gsell committed
139 140 141 142 143 144 145 146 147 148 149 150 151

    void setSymmetry(double symmetry);
    virtual double getSymmetry() const;

    void   setRinit(double rinit);
    virtual double getRinit() const;

    void   setPRinit(double prinit);
    virtual double getPRinit() const;

    void   setPHIinit(double phiinit);
    virtual double getPHIinit() const;

152 153
    void   setZinit(double zinit);
    virtual double getZinit() const;
154

155 156 157
    void   setPZinit(double zinit);
    virtual double getPZinit() const;

gsell's avatar
gsell committed
158 159 160 161
    void   setBScale(double bs);
    virtual double getBScale() const;

    void   setEScale(std::vector<double> bs);
162
    virtual double getEScale(unsigned int i) const;
gsell's avatar
gsell committed
163

ext-calvo_p's avatar
ext-calvo_p committed
164
    void setTrimCoils(const std::vector<TrimCoil*>& trimcoils);
165

166
    void setSuperpose(std::vector<bool> flag);
167
    virtual bool getSuperpose(unsigned int i) const;
168

169 170 171 172 173 174 175 176 177 178
    void setMinR(double r);
    virtual double getMinR() const;
    void setMaxR(double r);
    virtual double getMaxR() const;

    void setMinZ(double z);
    virtual double getMinZ() const;
    void setMaxZ(double z);
    virtual double getMaxZ() const;

179 180 181 182
    void setFMLowE(double e);
    virtual double getFMLowE() const;
    void setFMHighE(double e);
    virtual double getFMHighE() const;
183

184 185 186
    void setTrimCoilThreshold(double);
    virtual double getTrimCoilThreshold() const;

187 188 189
    void setSpiralFlag(bool spiral_flag);
    virtual bool getSpiralFlag() const;

ext-calvo_p's avatar
ext-calvo_p committed
190
    virtual bool apply(const size_t& id, const double& t, Vector_t& E, Vector_t& B);
gsell's avatar
gsell committed
191

ext-calvo_p's avatar
ext-calvo_p committed
192
    virtual bool apply(const Vector_t& R, const Vector_t& P, const double& t, Vector_t& E, Vector_t& B);
gsell's avatar
gsell committed
193

194 195 196 197
    virtual void apply(const double& rad, const double& z,
                       const double& tet_rad, double& br,
                       double& bt, double& bz);

ext-calvo_p's avatar
ext-calvo_p committed
198
    virtual void initialise(PartBunchBase<double, 3>* bunch, double& startField, double& endField);
gsell's avatar
gsell committed
199

ext-calvo_p's avatar
ext-calvo_p committed
200
    virtual void initialise(PartBunchBase<double, 3>* bunch, const double& scaleFactor);
gsell's avatar
gsell committed
201 202 203 204 205 206 207 208

    virtual void finalise();

    virtual bool bends() const;

    virtual double getRmax() const;

    virtual double getRmin() const;
209

210 211 212 213 214 215
    bool interpolate(const double& rad,
                     const double& tet_rad,
                     double& br,
                     double& bt,
                     double& bz);
    
ext-calvo_p's avatar
ext-calvo_p committed
216 217 218 219 220 221
    void read(const double& scaleFactor);

    void setBFieldType();

    BFieldType fieldType_m;

222
private:
223 224 225
    /// Apply trim coils (calculate field contributions) with smooth field transition
    void applyTrimCoil  (const double r, const double z, const double tet_rad, double& br, double& bz);
    /// Apply trim coils (calculate field contributions)
ext-calvo_p's avatar
ext-calvo_p committed
226 227
    void applyTrimCoil_m(const double r, const double z, const double tet_rad, double* br, double* bz);

228

frey_m's avatar
frey_m committed
229
protected:
ext-calvo_p's avatar
ext-calvo_p committed
230
   
frey_m's avatar
frey_m committed
231 232
    void   getdiffs();

ext-calvo_p's avatar
ext-calvo_p committed
233
    double gutdf5d(double* f, double dx, const int kor, const int krl, const int lpr);
frey_m's avatar
frey_m committed
234 235

    void   initR(double rmin, double dr, int nrad);
gsell's avatar
gsell committed
236

ext-calvo_p's avatar
ext-calvo_p committed
237 238 239 240 241 242 243
    void   getFieldFromFile_Ring(const double& scaleFactor);
    void   getFieldFromFile_Carbon(const double& scaleFactor);
    void   getFieldFromFile_CYCIAE(const double& scaleFactor);
    void   getFieldFromFile_AVFEQ(const double& scaleFactor);
    void   getFieldFromFile_FFA(const double& scaleFactor);
    void   getFieldFromFile_BandRF(const double& scaleFactor);
    void   getFieldFromFile_Synchrocyclotron(const double& scaleFactor);
frey_m's avatar
frey_m committed
244 245

    inline int idx(int irad, int ktet) {return (ktet + Bfield.ntetS * irad);}
246

ext-calvo_p's avatar
ext-calvo_p committed
247

gsell's avatar
gsell committed
248 249 250 251
private:

    std::string fmapfn_m; /* stores the filename of the fieldmap */
    std::vector<double> rffrequ_m;
252
    std::vector< std::vector<double> > rffc_m;
253 254
    std::vector<double> rfvrequ_m;
    std::vector< std::vector<double> > rfvc_m;
gsell's avatar
gsell committed
255
    std::vector<double> rfphi_m;
256
    std::vector<double> escale_m;  // a scale factor for the E-field
257
    std::vector<bool> superpose_m; // electric fields are superposed or not
gsell's avatar
gsell committed
258 259 260 261 262 263

    double symmetry_m;

    double rinit_m;
    double prinit_m;
    double phiinit_m;
264 265
    double zinit_m;
    double pzinit_m;
gsell's avatar
gsell committed
266

267
    bool spiral_flag_m;
268
    double trimCoilThreshold_m; ///< B-field threshold for applying trim coil
269

ext-calvo_p's avatar
ext-calvo_p committed
270
    std::string typeName_m; // name of the TYPE parameter in cyclotron
gsell's avatar
gsell committed
271 272 273 274
    double harm_m;

    double bscale_m; // a scale factor for the B-field

275 276
    /// Trim coils
    std::vector<TrimCoil*> trimcoils_m;
277

278 279 280 281 282 283
    double minr_m;
    double maxr_m;

    double minz_m;
    double maxz_m;

284 285 286
    double fmLowE_m;
    double fmHighE_m;

gsell's avatar
gsell committed
287
    // Not implemented.
288
    void operator=(const Cyclotron &) = delete;
gsell's avatar
gsell committed
289 290 291

    // RF field map handler
    //    Fieldmap *RFfield;
292
    std::vector<Fieldmap *> RFfields_m;
gsell's avatar
gsell committed
293
    std::vector<std::string> RFfilename_m;
294
    std::vector<std::string> RFFCoeff_fn_m;
295
    std::vector<std::string> RFVCoeff_fn_m;
296

297
    // handling for store the particle out of region
adelmann's avatar
adelmann committed
298
    std::unique_ptr<LossDataSink> lossDs_m;
299

300
    // Necessary for quick and dirty phase output -DW
301
    int waiting_for_gap = 1;
302

frey_m's avatar
frey_m committed
303
protected:
304
    // object of Matrices including magnetic field map and its derivates
frey_m's avatar
frey_m committed
305
    BfieldData Bfield;
gsell's avatar
gsell committed
306

frey_m's avatar
frey_m committed
307 308
    // object of parameters about the map grid
    BPositions BP;
gsell's avatar
gsell committed
309 310
};

311
#endif // CLASSIC_Cyclotron_HH