24template<
unsigned int Dim>
29 template<
class T,
class FT,
class FT2>
31 grn = grnI[0] * hrsq[0] + grnI[1] * hrsq[1] + grnI[2] * hrsq[2];
32 grn = 1.0 /
sqrt(grn);
33 grn[0][0][0] = grn[0][0][1];
42 greensFunction_m(greensFunction),
63 for(i = 0; i < 3; i++) {
76 bool sineTransformDims[3];
77 for(
int d = 0; d < 3; ++d)
78 sineTransformDims[d] =
true;
84 for(i = 0; i < 3; ++i) {
111 mesh_m(&beam.getMesh()),
112 layout_m(&beam.getFieldLayout()),
115 greensFunction_m(greensFunction)
133 for(i = 0; i < 3; i++) {
142 bool sineTransformDims[3];
143 for(
int d = 0; d < 3; ++d)
144 sineTransformDims[d] =
true;
151 for(i = 0; i < 3; ++i) {
200 sine_m->transform(-1, rho);
218 sine_m->transform(+1, rho);
224 rho[I][J][K] = rho[I][J][
nr_m[2] - K - 1];
244 sine_m->transform(-1, rho);
261 sine_m->transform(+1, rho);
316 for(
int k = idx[2].first(); k <
std::min(
nr_m[2] + 2, idx[2].last() + 3); k++) {
317 for(
int j = idx[1].first(); j <
std::min(
nr_m[1] + 2, idx[1].last() + 3); j++) {
318 for(
int i = idx[0].first(); i <
std::min(
nr_m[0] + 2, idx[0].last() + 3); i++) {
325 r =
sqrt(vv(0) * vv(0) + vv(1) * vv(1) + vv(2) * vv(2));
326 tmpgrn = -vv(2) * vv(2) *
atan(vv(0) * vv(1) / (vv(2) * r)) / 2;
327 tmpgrn += -vv(1) * vv(1) *
atan(vv(0) * vv(2) / (vv(1) * r)) / 2;
328 tmpgrn += -vv(0) * vv(0) *
atan(vv(1) * vv(2) / (vv(0) * r)) / 2;
329 tmpgrn += vv(1) * vv(2) *
log(vv(0) + r);
330 tmpgrn += vv(0) * vv(2) *
log(vv(1) + r);
331 tmpgrn += vv(0) * vv(1) *
log(vv(2) + r);
400 for(
int k = idx[2].first(); k <
std::min(
nr_m[2] + 2, idx[2].last() + 3); k++) {
401 for(
int j = idx[1].first(); j <
std::min(
nr_m[1] + 2, idx[1].last() + 3); j++) {
402 for(
int i = idx[0].first(); i <
std::min(
nr_m[0] + 2, idx[0].last() + 3); i++) {
407 vv(2) = k *
hr_m[2] -
hr_m[2] / 2 + zshift;
409 r =
sqrt(vv(0) * vv(0) + vv(1) * vv(1) + vv(2) * vv(2));
410 tmpgrn = -vv(2) * vv(2) *
atan(vv(0) * vv(1) / (vv(2) * r)) / 2;
411 tmpgrn += -vv(1) * vv(1) *
atan(vv(0) * vv(2) / (vv(1) * r)) / 2;
412 tmpgrn += -vv(0) * vv(0) *
atan(vv(1) * vv(2) / (vv(0) * r)) / 2;
413 tmpgrn += vv(1) * vv(2) *
log(vv(0) + r);
414 tmpgrn += vv(0) * vv(2) *
log(vv(1) + r);
415 tmpgrn += vv(0) * vv(1) *
log(vv(2) + r);
426 for(
int k = idx[2].first(); k <
std::min(
nr_m[2] + 2, idx[2].last() + 3); k++) {
427 for(
int j = idx[1].first(); j <
std::min(
nr_m[1] + 2, idx[1].last() + 3); j++) {
428 for(
int i = idx[0].first(); i <
std::min(
nr_m[0] + 2, idx[0].last() + 3); i++) {
436 r =
sqrt(vv(0) * vv(0) + vv(1) * vv(1) + vv(2) * vv(2));
437 tmpgrn = -vv(2) * vv(2) *
atan(vv(0) * vv(1) / (vv(2) * r)) / 2;
438 tmpgrn += -vv(1) * vv(1) *
atan(vv(0) * vv(2) / (vv(1) * r)) / 2;
439 tmpgrn += -vv(0) * vv(0) *
atan(vv(1) * vv(2) / (vv(0) * r)) / 2;
440 tmpgrn += vv(1) * vv(2) *
log(vv(0) + r);
441 tmpgrn += vv(0) * vv(2) *
log(vv(1) + r);
442 tmpgrn += vv(0) * vv(1) *
log(vv(2) + r);
470 ggrn2[I][J][K] = grn2[I+1][J+1][K+1];
471 ggrn2[I][J][K] += grn2[I+1][J][K];
472 ggrn2[I][J][K] += grn2[I][J+1][K];
473 ggrn2[I][J][K] += grn2[I][J][K+1];
474 ggrn2[I][J][K] -= grn2[I+1][J+1][K];
475 ggrn2[I][J][K] -= grn2[I+1][J][K+1];
476 ggrn2[I][J][K] -= grn2[I][J+1][K+1];
477 ggrn2[I][J][K] -= grn2[I][J][K];
513 os <<
"* ************* F F T P o i s s o n S o l v e r ************************************ " <<
endl;
514 os <<
"* h " <<
hr_m <<
'\n';
515 os <<
"* ********************************************************************************** " <<
endl;
Tps< T > log(const Tps< T > &x)
Natural logarithm.
Tps< T > sqrt(const Tps< T > &x)
Square root.
T::PETE_Expr_t::PETE_Return_t min(const PETE_Expr< T > &expr, NDIndex< D > &loc)
PETE_TUTree< FnArcTan, typename T::PETE_Expr_t > atan(const PETE_Expr< T > &l)
PETE_TTTree< OpWhere, typename Cond_t::PETE_Expr_t, typename True_t::PETE_Expr_t, PETE_Scalar< Vektor< T, Dim > > > where(const PETE_Expr< Cond_t > &c, const PETE_Expr< True_t > &t, const Vektor< T, Dim > &f)
Inform & endl(Inform &inf)
bool lt(double x, double y)
static void calculate(Vektor< T, 3 > &hrsq, FT &grn, FT2 *grnI)
IpplTimings::TimerRef ShIntGreensFunctionTimer4_m
IpplTimings::TimerRef GreensFunctionTimer_m
void computePotential(Field_t &rho, Vector_t hr, double zshift)
Inform & print(Inform &os) const
IpplTimings::TimerRef ShIntGreensFunctionTimer2_m
std::string greensFunction_m
FFT< SineTransform, 3, double > SINE_t
IpplTimings::TimerRef IntGreensFunctionTimer4_m
FieldLayout_t * layout2_m
IpplTimings::TimerRef IntGreensFunctionTimer3_m
void shiftedIntGreensFunction(double zshift)
compute the shifted integrated Green function as described in Three-dimensional quasistatic model for...
void integratedGreensFunction()
compute the integrated Green function as described in Three-dimensional quasistatic model for high br...
IpplTimings::TimerRef GreensFunctionTimer4_m
IpplTimings::TimerRef IntGreensFunctionTimer2_m
FFTBoxPoissonSolver(PartBunch &bunch, std::string greensFuntion)
IpplTimings::TimerRef ShIntGreensFunctionTimer3_m
IpplTimings::TimerRef GreensFunctionTimer1_m
IpplTimings::TimerRef IntGreensFunctionTimer1_m
IpplTimings::TimerRef ShIntGreensFunctionTimer1_m
void initialize(Layout_t &)
const NDIndex< Dim > & getDomain() const
NDIndex< Dim > getLocalNDIndex()
const NDIndex< Dim > & getDomain() const
MFLOAT get_meshSpacing(unsigned d) const
static TimerRef getTimer(const char *nm)
static void stopTimer(TimerRef t)
static void startTimer(TimerRef t)
Vektor< double, 3 > Vector_t