Kinematics in Magnetic field

Leave a comment

Generally, the best starting point of  any kinematics calculation for any reaction is from the CM frame. In the exit channel, the four-momentum in the CM frame is

P_{cm} = \begin{pmatrix} \sqrt{m^2 + p^2} = E_{cm} \\  \vec{p} \end{pmatrix}

the boost from CM frame to Lab frame is \vec{\beta} . Lets define the perpendicular vector on the plane of \vec{p} be \hat{n} , thus, the four-momentum in the Lab frame is

P = \begin{pmatrix} E \\ \vec{k} \end{pmatrix} =\begin{pmatrix} \gamma E_{cm} + \gamma \beta \hat{\beta}\cdot \vec{p} \\ \left(  \gamma \beta E_{cm} + \gamma (\hat{\beta} \cdot \vec{p}) \right) \hat{\beta} + (\hat{n} \cdot \vec{p} ) \hat{n} \end{pmatrix}

Suppose the magnetic field is parallel to z-axis

\vec{B} = B \hat{z}

The rotation radius is

\displaystyle \rho = \frac{\vec{k} \cdot \hat{xy} }{cZB}

where Z is the charge state of the particle, and \hat{xy} is the perpendicular unit vector on xy-plane.


The time for 1 cycle is

\displaystyle T_{cyc} = \frac{2\pi \rho}{v_{xy}} = \frac{2\pi}{cZB} \frac{\vec{k}\cdot\hat{xy}}{v_{xy}} = \frac{2\pi}{c^2 ZB} E

And the length for 1 cycle is

\displaystyle z_{cyc} = v_z T_{cyc} = \frac{2\pi}{cZB} \vec{k}\cdot \hat{xy} \frac{v_z}{v_{xy}} = \frac{2\pi}{cZB} \vec{k}\cdot \hat{z}

When the Lorentz boost is parallel to the B-field direction,

\vec{\beta} \cdot \hat{z} = 1, \hat{n} \cdot \hat{z} = 0

\displaystyle z_{cyc} = \frac{2\pi}{cZB} \left(\gamma \beta E_{cm} + \gamma (\hat{\beta} \cdot \vec{p})\right)

and the total energy is

\displaystyle E = \gamma E_{cm} + \gamma \beta (\hat{\beta}\cdot \vec{p})

By eliminating the \hat{\beta} \cdot \vec{p} ,

\displaystyle E = \frac{1}{\gamma}E_{cm} + \frac{cZB}{2\pi} \beta z_{cyc}

We can see that, the energy and z_{cyc} have a linear relation. That simplified a alot of thing. For instance, the intercept is related to E_{cm}, which is further related to the Q-value of the reaction. If we can plot the E-z plot and determine the intercept, we can extract the Q-value, which is related to the excitation energy or some sort.

Without the magnetic field, to extract the Q-value, we really need to know the scattering angle (the azimuthal angle can be skipped due to symmetry in some cases). And usually, the relation is not linear.


Magnetic Rigidity

Leave a comment

In ion optics, there is a concept called Magnetic Rigidity. When an charged ion passing through a magnetic field, it will bend and rotate. The magnetic rigidity is the B\rho, where B is the magnetic field strength in Tesla, and \rho is the rotation radius.

using basics physics, we have the centripetal force equals to the Lorentz force,

\displaystyle \frac{mv^2}{\rho} = Q \vec{v}\times \vec{B}


B\rho = mv/Q

The above is non-relativistic, and I always assume the mv can be replaced by relativistic momentum p = m \gamma \beta . Now I give a prove.

The relativistic centripetal force is

\displaystyle \vec {F} = m \gamma \vec{a_\perp} = - m \gamma\frac{ v^2}{\rho} \hat{\rho}


\displaystyle m \gamma \frac{v^2}{\rho} = QvB

\displaystyle m \gamma \beta = p = cQB \rho

where c is speed of light. In the last equation, the unit of p is MeV/c, c = 299.792458 mm/ns, and \rho in meter.

Tunneling through the earth core


I do something difference today, ha.

The problem is simple, How long does it take to tunneling through the center of the earth to the other side?

Assuming uniform density of the earth (\rho ), the acceleration inside the earth at radius r is

\displaystyle  - G \frac{4\pi}{3} \rho r = -k r

k = 1.5413\times 10^{-6} in 1/sec^2.

because only the mass within the radius matter. Thus, the equation of motion is

\displaystyle \frac{d^2 r}{dt^2} = - k r 

the solution is

\displaystyle r(t) = R \cos( \sqrt{k} t )

Thus, the time for a trip is T = \pi / \sqrt{k}   = 2530.5 sec or 42 min and 10.5 sec.

The maximum speed when passing the core. The speed is R \sqrt{k} = 7910 m/s

How about we use a realistic earth density?

The density, and acceleration can be found in the web, for example, here.

The travel time is 2291 sec, or 38 min 11 sec.


It is interesting that, in the uniform density calculation, the travel time is independent of the radius, but density. The peak velocity is

\displaystyle R \sqrt{k} < c \rightarrow \rho < \frac{c^2}{R^2}\frac{3}{4\pi} \frac{1}{G}

Let compare with Schwarzschild Radius:

\displaystyle R_S = \frac{2GM}{c^2} \rightarrow \rho_S = \frac{1}{2}\frac{c^2}{R^2}\frac{3}{4\pi} \frac{1}{G}

The Schwarzchild density is half of the maximum density by classical argument.

In case of point mass.

The equation of motion is

\displaystyle \frac{d^2r}{dt^2} = -\frac{GM}{r^2}

change of variable v(r) = \frac{dr}{dt}

\displaystyle \frac{d(v^2/2)}{dr} = -\frac{GM}{r^2}

\displaystyle v^2(r) = 2 \frac{GM}{r} - 2 \frac{GM}{R}

The above solution assume v(R) = 0 . We can see that, at r = 0 , the speed go to infinity. Something wrong…..

Update: that is not wrong at all. imagine two neutrino with head-on collision, the released gravitational energy will be infinity! But, because of uncertainly principle, their separation distance can never to be zero.

Taking account of the earth rotation with the centripetal force, the equation of motion becomes,

\displaystyle \frac{d^2 r}{dt^2} = - k r + \omega^2\cos(\theta) r

where \omega is the angular velocity of earth, and \theta is the polar angle from the north pole. The travel time would be

\displaystyle T = \frac{\pi} {\sqrt{k - \omega^2 \cos(\theta)}}

\displaystyle \omega = \frac{2\pi}{24 \times 60\times 60} = 7.3 \times 10^{-5} rad/s


\omega^2 = 5.3 \time 10^{-9} 1/sec^2.

which is a very small correction.


Simple model for 4He and NN-interaction

Leave a comment

Starting from deuteron, the binding energy, or the p-n interaction is 2.2 MeV.

From triton, 3H, the total binding energy is 8.5 MeV, in which, there are only 3 interactions, two p-n and one n-n. Assume the p-n interaction does not change, the n-n interaction is 4.1 MeV.

The total binding energy of 3He is 7.7 MeV. The p-p interaction is 3.3 MeV.

Notices that we neglected the 3-body force in 3H and 3He. And it is strange that the n-n and p-p interaction is stronger then p-n interaction.

In 4He, the total binding energy becomes 28.3 MeV. I try to decompose the energy in term of 2-body, 3-body, and 4-body interaction.

If we only assume 2-body interaction, the interaction strength from n-p, n-n, and p-p are insufficient. One way to look is the 1-particle separation energy.

The neutron separation energy is 20.6 MeV = 2(p-n) + (n-n).
The proton separation energy is 19.8 MeV = 2(p-n) + (p-p).
The total energy is 28.3 MeV = 4(p-n) + (n-n) + (p-p).

There is no solution for above 3 equations. Thus, only consider 2-body interaction is not enough.

The neutron separation energy is 20.6 MeV = 2(p-n) + (n-n) + 2(n-n-p) + (n-p-p)
The proton separation energy is 19.8 MeV = 2(p-n) +(p-p) + 2(n-p-p) + (n-n-p).
The total energy is 28.3 MeV = 4(p-n) + (n-n) + (p-p) + 2(n-n-p) + 2(n-p-p) + (n-n-p-p).

Assuming the 2-body terms are the same in 2H, 3H, and 3He, the (n-p-p) and (n-n-p) is 4.03 MeV, which is strange again, as the Coulomb repulsion should make the (n-p-p) interaction smaller then the (n-n-p) interaction. The (n-n-p-p) interaction is -4.03 MeV.

Lets also add 3-body force in 3H and 3He.

The neutron separation energy of 3H is 6.3 MeV = (p-n) + (n-n) + (n-n-p)
The toal energy of 3H is 8.5 MeV = 2(p-n) + (n-n) + (n-n-p)
The toal energy of 3He is 7.7 MeV = 2(p-n) + (p-p) + (n-p-p)
The neutron separation energy of 4He is 20.6 MeV = 2(p-n) + (n-n) + 2(n-n-p) + (n-p-p)
The proton separation energy is 4He 19.8 MeV = 2(p-n) +(p-p) + 2(n-p-p) + (n-n-p).
The total energy of 4He is 28.3 MeV = 4(p-n) + (n-n) + (p-p) + 2(n-n-p) + 2(n-p-p) + (n-n-p-p).

We have 6 equations, with 6 unknown [ (p-n) , (n-n) , (p-p) , (n-n-p), (n-p-p), and (n-n-p-p)]. Notice that the equation from proton separation energy of 3He is automatically satisfied. The solution is

(p-n) = 2.2 MeV
(p-p) = – 4.7 – (n-n) MeV
(n-n-p)  = 4.1 – (n-n) MeV
(n-p-p) = 8 + (n-n) MeV
(n-n-p-p) = 0 MeV

It is interesting that there is redundant equation. But still, the (p-n) interaction is 2.2 MeV, and 4-body (n-n-p-p) becomes 0 MeV. Also the (n-p-p) is more bound than (n-n-p) by 3.9 + 2(n-n) MeV. If (n-p-p) should be more unbound, than (n-n) must be negative and smaller than -1.95 MeV.

Since the interaction strength has to be on the s-orbit (mainly), by considering 2H, 3H, 3He, and 4He exhausted all possible equations (I think). We need other way to anchor either (n-n), (p-p), (n-p-p), and (n-n-p) interactions.

Use the Coulomb interaction, the Coulomb interaction should add -1.44 MeV on the NN pair (assuming the separation is a 1 fm). Lets assume the (n-n) – (p-p) = 1.44 MeV

(n-n) = -1.63 MeV
(p-p) = – 3.07 MeV
(n-n-p)  = 5.73 MeV
(n-p-p) = 6.37  MeV

The (p-p) is more unbound than (n-n) as expected, but the (n-p-p) is more bound than (n-n-p) by 0.64 MeV. This is surprising! We can also see that, the 3-body interaction play an important role in nuclear interaction.

According to this analysis, the main contribution of the binding energies of 3H and 3He are the 3-body force.

In 3H:  (n-n) + 2(n-p) + (n-n-p) = -1.6 + 4.4 + 5.7 = 8.5 MeV
In 3He: (p-p) + 2(n-p) + (n-p-p) = -3.1 + 4.4 + 6.4 = 7.7 MeV

Worked on the algebra, when ever the difference  (n-n) – (p-p)  > 0.8 MeV, the (n-p-p) will be more bound that (n-n-p). Thus, the average protons separation should be more than 1.8 fm. I plot the interactions energies with the change of Coulomb energy below.


The (n-n) and (p-p) are isoscalar pair, where tensor force is zero. While the (n-p) quasi-deuteron is isovector pair. Thus, the difference between (n-n) and (n-p) reflect the tensor force in s-orbit, which is 3.8 MeV. In s-orbit, there is no spin-orbital interaction, therefore, we can regard the tensor force is 3.8 MeV for (n-p) isovector pair.

Following this method, may be, we can explore the NN interaction in more complex system, say the p-shell nuclei. need an automatic method. I wonder the above analysis agreed present interaction theory or not. If not, why? 

Particle nano-Ampere (pnA)

Leave a comment

I am so confused between the unit of nA, pnA, enA, and pps.

pps = particle per second.

pnA  = particle nanoampere.
The electrical current in nanoamperes (10^{-9}A) that would be measured if all beam ions were singly charged. i.e. neglecting the actual charge state. 1 \textrm{pnA} =  6.25 \times 10^{9} \textrm{ions/second}

enA = pnA × Charge state.

In this sense, enA = nA, the actual current carried by the beam.

In General

\textrm{pnA} = \textrm{pps} \times e \times 10^9

\textrm{enA}  = \textrm{nA} = Z \times \textrm{pps} \times e \times 10^9

where e is electron charge.

For example, 12C with charge state of 6+ at 4 \times 10^{10} pps.

= 4 \times 10^{10} \times e \times 10^{9} = 6.41 \textrm{pnA}


= 6 \times 4 \times 10^{10} \times e \times 10^{9} = 38.5 \textrm{nA} = 38.5 \textrm{enA}

An ion beam of 5+ charge state at 1 pnA.

= 1 / e  10^{-9} = 6.24 \textrm{pps}


= 5 \textrm{enA}  = 5 \textrm{nA}

Sum of sine-square and n-root of unity

Leave a comment

Suddenly, I encounter this problem: find the sum

\displaystyle \sum_{n=1}^{N} \sin^2\left(x + \frac{2\pi}{N} n \right)

Here is my way of thinking,


\displaystyle \sin(a_n) = \frac{1}{2i} (e^{ia_n} - e^{-ia_n} )


a_n = x + \frac{2\pi}{N} n


\displaystyle \sin^2(a) = \frac{1}{4} (2- e^{2ia} - e^{-2ia} )

The sum break down to calculate

\displaystyle \sum_{n=1}^{N} \sin^2\left(x + \frac{2\pi}{N} n \right)  = \frac{N}{2} - \frac{1}{4} \sum_{n=1}^N (e^{2ia} + e^{-2ia})

The sum is related to n-root of unity as

\displaystyle \sum_{n=1}^N e^{2ia} = 0


\displaystyle z^N = 1 \rightarrow z^N-1 = 0 \rightarrow \Sigma_{n=1}^N (z - e^{i 2\pi n / N}) = 0

break down the product,

\displaystyle z^N + \sum_{n=1}^N e^{i 2\pi n/N} z^N + .... + \Sigma_{n=1}^N e^{i 2\pi n /N} = 0


\displaystyle  \sum_{n=1}^N e^{i 2\pi n/N}  = 0

And since \sqrt{1} = 1, thus, z^{N/2} = 1 and

\displaystyle  \sum_{n=1}^N e^{i 4\pi n/N}  = 0

The graphical imagination is that, the solutions of n-root of unity are the vectors,  isotropically  distributed from center of a unit circle to its circumference. Thus, the add up of all these vectors will form a N-size polygon, and the result of the sum is zero.

I found it is very interesting to transform one question into another question that can be solved graphically.

This problem is originated from 3-phase AC power that N = 3. I wonder, is it the total power be a constant for N >= 3, and the answer is yes.

The ultimate reason for that is, when generating AC power, the motor is rotating at constant speed and a constant power is inputting to the system, which can using 1-phase (N = 2), 2-phase (N = 4), 3-phase (N = 3), and so on, to generate the AC power. Thus, as long as the phase distributed equally, the sum of power must be a constant.

Stopping power and Bethe-Bloch Formula

Leave a comment

I am so surprised that this topic is not in this blog.  But I am not always posting what I did. Anyway…

The Stopping power is energy loss per length, in unit of MeV/cm.

\displaystyle S = -\frac{dE}{dx}

Since the stopping power is proportional to density, it is often normalized with the density in literature and the unit becomes MeV/(ug/cm^2).

The Bethe-Bloch formula based on classical argument.

In order to calculate the range, we can integrate the stopping power. Consider the particle moved by \Delta x distance, the energy loss is

E(x+\Delta x) = E(x) - S(E(x)) \Delta x

rearrange, gives

\displaystyle  \frac{dE(x)}{dx} = - S(E)

\displaystyle \int_{E_0}^{E} \frac{-1}{S(E)} dE = x(E ; E_0)

This is the relation between particle energy and range.

We can plot S(E(x)) to get the Bragg peak.

In here, I use the physical calculator in LISE++, SRIM, and Bethe-Bloch formula to calculate the stopping power for proton in CD2 target.

The atomic mass of deuteron is 2.014102 u. The density of CD2 target is 0.913 g/cm3.


For the Bethe curve, I adjusted the density to be 0.9 mg/cm3, and simple use the carbon charge. And the excitation potential to be 10^{-5} z, where z = 6 , so that the Bethe curve agree with the others.

It is well known that the Bethe curve fails at low energy.

Assume the proton is 5 MeV. The energy loss against distance isEvsX.PNG

and the stopping power against distance is


We can see the Bragg peak is very sharp and different models gives different stopping range. It is due to the rapid decrease of energy at small energy. In reality, at small energy, the microscopic effect becomes very important and statistical. thus, the curve will be smoothed and the Bragg peak will be broaden.

In very short range, the stopping power increase almost linearly, Because the incident energy is large, so that the stopping power is almost constant around that energy range.

For S(E_0) = h ,

\displaystyle  x = \frac{-1}{h}(E - E_0 )

\displaystyle E = E_0 - hx,  hx << E_0

\displaystyle S(E) = S(E_0 - hx) \approx S(E_0) - \frac{dS(E_0)}{dE}(hx)



Older Entries