19 filters_m(filters.begin(), filters.end()),
36 std::pair<double, double> meshInfo;
38 const double &meshOrigin = meshInfo.first;
39 const double &meshSpacing = meshInfo.second;
42 if(
Ez_m.size() < numOfSlices) {
43 Ez_m.resize(numOfSlices, 0.0);
44 Chi_m.resize(numOfSlices, 0.0);
45 Grn_m.resize(numOfSlices, 0.0);
46 Psi_m.resize(numOfSlices, 0.0);
49 for(
unsigned int i = 0; i < numOfSlices; ++i) {
56 for(
unsigned int i = 1; i < numOfSlices; i++) {
57 double pathLengthOfSlice = minPathLength + i * meshSpacing;
69 for(
unsigned int i = 0; i < bunch->
getLocalNum(); ++i) {
71 unsigned int indexz = (
unsigned int)
floor((
R(2) - meshOrigin) / meshSpacing);
72 double leverz = (
R(2) - meshOrigin) / meshSpacing - indexz;
75 bunch->
Ef[i](2) += (1. - leverz) *
Ez_m[indexz] + leverz *
Ez_m[indexz + 1];
79 static std::string oldBendName;
80 static unsigned long counter = 0;
85 bool print_criterion = (counter + 1) % every == 0;
87 static unsigned int file_number = 0;
88 if(counter == 0) file_number = 0;
91 std::stringstream filename_str;
92 filename_str <<
"data/" <<
bendName_m <<
"-CSRWake" << std::setw(5) << std::setfill(
'0') << file_number <<
".txt";
93 std::ofstream csr(filename_str.str().c_str());
94 csr << spos <<
", " <<
FieldBegin_m <<
", " << smin(2) <<
", " << smax(2) <<
", " << meshSpacing*64 <<
std::endl;
96 csr << i *meshSpacing <<
"\t"
101 msg <<
"** wrote " << filename_str.str() <<
endl;
130 std::vector<Filter *>::const_iterator fit;
142 double xmu_const = 3.0 * gamma * gamma * gamma / (2.0 *
bendRadius_m);
143 double chi_const = 9.0 / 16.0 * (6.0 -
log(27.0 / 4.0));
145 for(
unsigned int i = 0; i < numOfSlices; ++i) {
147 double z = i * meshSpacing;
148 double xmu = xmu_const * z;
149 double b =
sqrt(xmu * xmu + 1.0) + xmu;
151 Chi_m[i] = chi_const + 0.5 *
pow(xmu, 2) - 7.0 / 54.0 *
pow(xmu, 4) + 140.0 / 2187.0 *
pow(xmu, 6);
153 Chi_m[i] = 9.0 / 16.0 * (3.0 * (-2.0 * xmu *
pow(b, 1.0/3.0) +
pow(b, 2.0/3.0) +
pow(b, 4.0/3.0)) +
154 log(
pow((1 -
pow(b, 2.0 / 3.0)) / xmu, 2) / (1 +
pow(b, 2.0 / 3.0) +
pow(b, 4.0 / 3.0))));
156 double grn_const = -16.0/(27.0 * gamma * gamma * meshSpacing);
158 Grn_m[numOfSlices - 1] = 0.0;
159 for(
unsigned int i = 1; i < numOfSlices - 1; ++i) {
167 int startSliceNum = 0;
168 for(
int j = sliceNumber; j >= startSliceNum; j--)
182 if(Ds_max2 / meshSpacing < sliceNumber) {
183 j = sliceNumber -
static_cast<int>(
floor(Ds_max2 / meshSpacing));
184 frac = Ds_max2 / meshSpacing - (sliceNumber - j);
190 if(Ds_max / meshSpacing < sliceNumber) {
191 j = sliceNumber -
static_cast<int>(
floor(Ds_max / meshSpacing));
192 frac = Ds_max / meshSpacing - (sliceNumber - j);
201 if(Ds_max / meshSpacing < sliceNumber) {
204 if(
Psi_m[k] > 0 &&
Psi_m[k] < totalBendAngle_m)
208 if(
Psi_m[0] > 0 &&
Psi_m[0] < totalBendAngle_m)
213 for(
unsigned int l = sliceNumber - k + 1; l < sliceNumber; ++ l) {
215 if(
Psi_m[l] > 0 &&
Psi_m[l] < totalBendAngle_m)
221 Psi_m[sliceNumber] =
calcPsi(
Psi_m[sliceNumber], angleOverlap, meshSpacing / 4.0);
222 if(
Psi_m[sliceNumber] > 0 &&
Psi_m[sliceNumber] < totalBendAngle_m)
226 Ez_m[sliceNumber] *= prefactor;
237 const int Nmax = 100;
238 const double eps = 1
e-10;
239 double residual = 0.0;
241 if(psiInitial != 0.0) psi = psiInitial;
243 for(
int i = 0; i < Nmax; ++i) {
244 residual =
bendRadius_m * psi * psi * psi * (psi + 4. * x) - 24. * Ds * psi - 24. * Ds * x;
254 ERRORMSG(
"In CSRWakeFunction::calcPsi(): exceed maximum number of iterations!" <<
endl);
258 return "CSRIGFWakeFunction";
void initialize(const ElementBase *ref)
std::vector< double > Grn_m
LineDensity lineDensity_m
PETE_TUTree< FnAbs, typename T::PETE_Expr_t > abs(const PETE_Expr< T > &l)
std::vector< Filter * > filters_m
constexpr double e
The value of .
virtual void getDimensions(double &sBegin, double &sEnd) const
ParticleAttrib< Vector_t > Ef
const unsigned int nBins_m
Interface for basic beam line object.
std::vector< double > Chi_m
LineDensity dlineDensitydz_m
void calculateContributionInside(size_t sliceNumber, double angleOfSlice, double meshSpacing)
double get_meanKineticEnergy() const
double getBendRadius() const
virtual const std::string & getName() const
Get element name.
Tps< T > log(const Tps< T > &x)
Natural logarithm.
virtual ElementType getType() const =0
Get element type std::string.
double getEffectiveLength() const
double getBendAngle() const
void calculateContributionAfter(size_t sliceNumber, double angleOfSlice, double meshSpacing)
virtual const std::string getType() const
double calcPsi(const double &psiInitial, const double &x, const double &Ds) const
constexpr double pi
The value of .
void apply(PartBunchBase< double, 3 > *bunch)
CSRIGFWakeFunction(const std::string &name, ElementBase *element, std::vector< Filter * > filters, const unsigned int &N)
Tps< T > pow(const Tps< T > &x, int y)
Integer power.
size_t getLocalNum() const
Tps< T > sqrt(const Tps< T > &x)
Square root.
std::shared_ptr< Filter > defaultFilter_m
double searchRoot(const double &tol)
bool hasPositiveRealRoots()
void calculateLineDensity(PartBunchBase< double, 3 > *bunch, std::pair< double, double > &meshInfo)
std::vector< double > Ez_m
constexpr double epsilon_0
The permittivity of vacuum in As/Vm.
double getEffectiveCenter() const
void calculateGreenFunction(PartBunchBase< double, 3 > *bunch, double meshSpacing)
void get_bounds(Vector_t &rmin, Vector_t &rmax)
virtual void calc_derivative(std::vector< double > &histogram, const double &h)=0
void calcLineDensity(unsigned int nBins, std::vector< double > &lineDensity, std::pair< double, double > &meshInfo)
calculates the 1d line density (not normalized) and append it to a file.
PETE_TUTree< FnFloor, typename T::PETE_Expr_t > floor(const PETE_Expr< T > &l)
Inform & endl(Inform &inf)
std::vector< double > Psi_m