68 bool revBeam,
bool revTrack):
78 for(
int i = 0; i < 6; ++i) {
79 for(
int j = 0; j <= 6; ++j) {
97 for(
int i = 0; i < 6; ++i) {
98 for(
int j = 0; j <= 6; ++j) {
102 for(
int j = 7; j < 28; ++j) {
344 double Ex = scale * sep.
getEx();
345 double Ey = scale * sep.
getEy();
373 TptFun pz =
sqrt(pt * pt - px * px - py * py);
380 TptFun xt = C * itsMap[
X] + S * itsMap[
Y];
381 TptFun yt = C * itsMap[
Y] - S * itsMap[
X];
385 itsMap[
X] = C * xt + (S / ks) * pxt;
386 itsMap[
Y] = C * yt + (S / ks) * pyt;
387 itsMap[
PX] = C * pxt - (S * ks) * xt;
388 itsMap[
PY] = C * pyt - (S * ks) * yt;
389 itsMap[
T] += pt * (refTime / E - length / pz);
435 static void makeFocus
436 (
double k,
double L,
double &
c,
double &s,
double &d,
double &f) {
437 double t = k * L * L;
440 s = L * (1.0 - t / 6.0);
441 d = L * L * (0.5 - t / 24.0);
442 f = L * L * L * (1.0 / 6.0 - t / 120.0);
450 double r =
sqrt(- k);
466 TptFun pz =
sqrt(pt * pt - px * px - py * py);
471 itsMap[
T] += pt * (refTime / E - length / pz);
479 double hx = scale * field.
normal(1);
480 double ex = hx *
tan(angle);
483 itsMap[
PY] -= ey * itsMap[
Y];
491 double hx = scale * field.
normal(1);
492 double ex = hx *
tan(angle);
495 itsMap[
PY] -= ey * itsMap[
Y];
500 (
double length,
double refLength,
double h,
520 double ks = (Fx[2] - Fy[1]) / 2.0;
522 double hx = h - Fx[0];
526 double kin = itsReference.getM() / itsReference.getE();
527 double refTime = refLength * kin * kin;
532 double cx, sx, dx, fx, cy, sy, dy, fy;
533 makeFocus(kx, length, cx, sx, dx, fx);
534 makeFocus(ky, length, cy, sy, dy, fy);
535 double wx = - kx * sx;
536 double wy = - ky * sy;
539 itsMap[
X] = cx * x + sx * px + dx * hx + h * dx * pt + x0;
540 itsMap[
PX] = wx * x + cx * px + sx * hx + h * sx * pt;
541 itsMap[
Y] = cy * y + sy * py + dy * hy + y0;
542 itsMap[
PY] = wy * y + cy * py + sy * hy;
543 itsMap[
T] += refTime * pt + h * (sx * x + dx * px + fx * hx + length * x0);
546 double s1 = (kx + ky) / 2.0;
547 double d1 = (kx - ky) / 2.0;
548 double root =
sqrt(d1 * d1 + ks * ks);
549 double c2 = d1 / root;
550 double s2 = ks / root;
553 double cu, su, du, fu, cv, sv, dv, fv;
554 double ku = s1 + (d1 * d1 - ks * ks) / root;
555 double kv = s1 - (d1 * d1 - ks * ks) / root;
556 makeFocus(ku, length, cu, su, du, fu);
557 makeFocus(kv, length, cv, sv, dv, fv);
558 double wu = - ku * su;
559 double wv = - kv * sv;
562 TptFun u = c2 * x - s2 * y;
563 TptFun v = c2 * y + s2 * x;
564 TptFun pu = c2 * px - s2 * py;
565 TptFun pv = c2 * py + s2 * px;
566 double hu = c2 * (h + hx) - s2 * hy;
567 double hv = c2 * hy + s2 * (h + hx);
573 ((cu + cv) * x + (cu - cv) * u + (su + sv) * px + (su - sv) * pu +
574 (du + dv) * hx + (du - dv) * hu + ((du + dv) * h + (du - dv) * bu) * pt) / 2.0;
576 ((wu + wv) * x + (wu - wv) * u + (cu + cv) * px + (cu - cv) * pu +
577 (su + sv) * hx + (su - sv) * hu + ((su + sv) * h + (su - sv) * bu) * pt) / 2.0;
579 ((cu + cv) * y - (cu - cv) * v + (su + sv) * py - (su - sv) * pv +
580 (du + dv) * hy - (du - dv) * hv - (du - dv) * bv * pt) / 2.0;
582 ((wu + wv) * y - (wu - wv) * v + (cu + cv) * py - (cu - cv) * pv +
583 (su + sv) * hy - (su - sv) * hv - (su - sv) * bv * pt) / 2.0;
584 itsMap[
T] += refTime * pt;
595 int order = field.
order();
605 Fx = + field.
normal(order);
606 Fy = - field.
skew(order);
612 Fx = Fxt + field.
normal(order);
613 Fy = Fyt - field.
skew(order);
641 int order = field.
order();
650 TptFun kxt = x * kx - y * ky;
651 TptFun kyt = x * ky + y * kx;
652 kx = kxt + field.
normal(order);
653 ky = kyt - field.
skew(order);
674 itsMap[
PY] -= Fy[0] + Fy[1] * itsMap[
X] + Fy[2] * itsMap[
Y];
684 TptFun pz1 =
sqrt(pt * pt - px1 * px1 - py1 * py1);
686 itsMap[
PX] = euclid.
M(0, 0) * px1 + euclid.
M(1, 0) * py1 + euclid.
M(2, 0) * pz1;
687 itsMap[
PY] = euclid.
M(0, 1) * px1 + euclid.
M(1, 1) * py1 + euclid.
M(2, 1) * pz1;
688 TptFun pz2 = euclid.
M(0, 2) * px1 + euclid.
M(1, 2) * py1 + euclid.
M(2, 2) * pz1;
693 euclid.
M(0, 0) * x + euclid.
M(1, 0) * y - euclid.
M(2, 0) * euclid.
getZ();
695 euclid.
M(0, 1) * x + euclid.
M(1, 1) * y - euclid.
M(2, 1) * euclid.
getZ();
697 euclid.
M(0, 2) * x + euclid.
M(1, 2) * y - euclid.
M(2, 2) * euclid.
getZ();
704 itsMap[
Y] = y2 - sByPz * itsMap[
PY];
705 itsMap[
T] += pt * (refTime / E + sByPz);
712 int order = field.
order();
729 Ae = Ae * x + field.
normal(i);
730 Ao = Ao * x - field.
skew(i);
733 Ae = + (Ae * (1.0 + h * x)).
integral(0);
734 Ao = - (Ao * (1.0 + h * x));
742 Series2 factor = h / (1.0 + h * x);
749 if(++k > order)
break;
753 Ao = Ao.derivative(0);
756 if(++k > order)
break;
virtual double getExitFaceRotation() const =0
Get exit pole face rotation.
virtual void visitSeptum(const Septum &)
Apply the algorithm to a Septum.
virtual ElementBase * getElement() const
Return the contained element.
void applyThinSBend(const BMultipoleField &field, double scale, double h)
Thin SBend kick.
virtual double getArcLength() const
Get arc length.
Euclid3D & offset() const
Return the offset.
virtual BMultipoleField & getField() override=0
Get multipole expansion of field.
PETE_TUTree< FnAbs, typename T::PETE_Expr_t > abs(const PETE_Expr< T > &l)
virtual void visitDrift(const Drift &)
Apply the algorithm to a Drift.
virtual void visitSBend(const SBend &)
Apply the algorithm to a SBend.
double getX() const
Get displacement.
virtual BDipoleField & getField()=0
Return the corrector field.
constexpr double e
The value of .
virtual void visitMonitor(const Monitor &)
Apply the algorithm to a Monitor.
virtual void visitPatch(const Patch &pat)
Apply the algorithm to a patch.
virtual double getAmplitude() const =0
Get RF amplitude.
double normal(int) const
Get component.
The field of a magnetic dipole.
virtual BMultipoleField & getField() override=0
Get multipole expansion of field.
virtual void visitSeparator(const Separator &)
Apply the algorithm to a Separator.
Interface for septum magnet.
virtual void visitRFCavity(const RFCavity &)
Apply the algorithm to a RFCavity.
virtual BMultipoleField & getField() override=0
Get multipole field.
virtual void trackMap(FVps< double, 6 > &, const PartData &, bool revBeam, bool revTrack) const
Track a map.
Interface for electrostatic separator.
Euclid3D getExitPatch() const
Get patch.
virtual ~TransportMapper()
virtual void trackMap(FVps< double, 6 > &map, const PartData &, bool revBeam, bool revTrack) const
Track a map.
virtual double getArcLength() const
Get arc length.
virtual double getPhase() const =0
Get RF phase.
virtual void visitBeamStripping(const BeamStripping &)
Apply the algorithm to a BeamStripping.
Interface for beam position monitors.
Interface for RF Quadrupole.
A simple arc in the XZ plane.
Tps< T > sin(const Tps< T > &x)
Sine.
virtual void visitAlignWrapper(const AlignWrapper &)
Apply the algorithm to an offset beamline object wrapper.
virtual void visitBeamBeam(const BeamBeam &)
Apply the algorithm to a BeamBeam.
virtual double getBz() const =0
Get solenoid field Bz in Teslas.
Define the position of a misaligned element.
double getZ() const
Get displacement.
virtual double getElementLength() const override
Get design length.
virtual Euclid3D getEntranceTransform() const
Get entrance patch.
virtual void setMap(const LinearMap< double, 6 > &)
Reset the linear part of the accumulated map for restart.
virtual double getEx() const =0
Get horizontal component Ex of field in V/m.
Interface for general corrector.
virtual double getEy() const =0
Get vertical component Ey of field in V/m.
Interface for beam diagnostics.
Euclid3D getEntrancePatch() const
Get patch.
virtual double getBy() const
Get vertical component.
virtual PlanarArcGeometry & getGeometry() override=0
Get SBend geometry.
Tps< T > tan(const Tps< T > &x)
Tangent.
virtual double getEntryFaceRotation() const =0
Get pole entry face rotation.
Interface for drift space.
virtual void visitMarker(const Marker &)
Apply the algorithm to a Marker.
double getBeta() const
The relativistic beta per particle.
virtual void visitMapIntegrator(const MapIntegrator &)
Apply the algorithm to an integrator capable of mapping.
double M(int row, int col) const
Get component.
Interface for general multipole.
static Euclid3D YRotation(double angle)
Make rotation.
T deg(T x)
Convert radians to degrees.
Euclid3D Inverse(const Euclid3D &t)
Euclidean inverse.
virtual double getElementLength() const override
Get design length.
Tps< T > cosh(const Tps< T > &x)
Hyperbolic cosine.
virtual double getExitFaceRotation() const =0
Get exit pole face rotation.
double getBendAngle() const
Get angle.
void applyDrift(double length)
virtual Euclid3D getExitTransform() const
Get exit patch.
virtual void visitRBend(const RBend &)
Apply the algorithm to a RBend.
double skew(int) const
Get component.
Transport function in N variables of type T.
double getQ() const
The constant charge per particle.
virtual void getMap(LinearMap< double, 6 > &) const
Return the linear part of the accumulated map.
virtual void accept(BeamlineVisitor &visitor) const =0
Apply visitor.
virtual double getElementLength() const
Get design length.
virtual void visitCorrector(const Corrector &)
Apply the algorithm to a Corrector.
virtual double getFrequency() const =0
Get RF frequencey.
constexpr double c
The velocity of light in m/s.
virtual void visitCyclotronValley(const CyclotronValley &)
Apply the algorithm to a CyclotronValley.
virtual void visitCCollimator(const CCollimator &)
Apply the algorithm to a collimator.
Interface for cyclotron collimator.
virtual void visitFlexibleCollimator(const FlexibleCollimator &)
Apply the algorithm to a flexible collimator.
virtual void visitProbe(const Probe &)
Apply the algorithm to a probe.
Displacement and rotation in space.
Abstract beam-beam interaction.
void setCoefficient(int index, const T &value)
Set coefficient.
void applyEntranceFringe(double edge, const BMultipoleField &field, double scale)
bool isIdentity() const
Test for identity.
Interface for cyclotron valley.
virtual double getBx() const
Get horizontal component.
virtual double getElementLength() const
Get element length.
An abstract sequence of beam line components.
virtual double getElementLength() const
Get element length.
Tps< T > sqrt(const Tps< T > &x)
Square root.
virtual double getEntryFaceRotation() const =0
Get pole entry face rotation.
FTps derivative(int var) const
Partial derivative.
The geometry for a RBend element.
virtual void visitComponent(const Component &)
Apply the algorithm to an arbitrary component.
double getP() const
The constant reference momentum per particle.
virtual void visitDiagnostic(const Diagnostic &)
Apply the algorithm to a Diagnostic.
The magnetic field of a multipole.
virtual void visitLambertson(const Lambertson &)
Apply the algorithm to a Lambertson.
double getM() const
The constant mass per particle.
void applyTransform(const Euclid3D &, double refLength)
Apply transform.
void applyExitFringe(double edge, const BMultipoleField &field, double scale)
virtual void visitDegrader(const Degrader &)
Apply the algorithm to a Degrader.
Tps< T > cos(const Tps< T > &x)
Cosine.
virtual void visitMultipole(const Multipole &)
Apply the algorithm to a Multipole.
void applyTransportMap(double length, double refLength, double h, const FTps< double, 2 > &Fx, const FTps< double, 2 > &Fy)
static void setGlobalTruncOrder(int order)
Set the global truncation order.
FTps< double, 2 > Series2
void applySBendBody(double length, double refLength, double h, const BMultipoleField &field, double scale)
double getCurvature() const
Get curvature.
FTps integral(int var, int trunc=EXACT) const
Partial integral.
virtual RBendGeometry & getGeometry() override=0
Get RBend geometry.
virtual double getBendAngle() const
Get angle.
Truncated power series in N variables of type T.
TransportMap< double, 6 > itsMap
Interface for a geometric patch.
Interface for a single beam element.
void applyThinMultipole(const BMultipoleField &field, double factor)
Thin multipole kick.
FTps< double, 2 > buildSBendVectorPotential(const BMultipoleField &, double h)
Construct the vector potential for a SBend.
virtual const Euclid3D & getPatch() const =0
Get patch transform.
const PartData itsReference
The reference information.
double getY() const
Get displacement.
Tps< T > sinh(const Tps< T > &x)
Hyperbolic sine.
static FTps makeVariable(int var)
Make variable.
int order() const
Return order.
virtual void visitParallelPlate(const ParallelPlate &)
Apply the algorithm to a ParallelPlate.
virtual void visitRFQuadrupole(const RFQuadrupole &)
Apply the algorithm to a RFQuadrupole.
void applyMultipoleBody(double length, double refLength, const BMultipoleField &field, double scale)
TransportFun< double, 6 > TptFun
void setCoefficient(int index, const T &value)
Set coefficient.
Interface for a Lambertson septum.
virtual void visitSolenoid(const Solenoid &)
Apply the algorithm to a Solenoid.