Sensitivity of Woods-Saxon parameters in Glauber model

July 15, 2021

GoLuckyRyan Basic , Leave a comment

Experimentally, people measure the reaction cross section and compare with Glauber model. Since the inputs of the Glauber model are the nuclear density function and the NN cross section. The NN cross section is fixed by experiment. For an assumed density function, one or two parameters have to be fitted to the experimental cross section.

This time, we explore the sensitivity of the Woods-Saxon parameters in Glauber model using the method.

From the above figure, the sensitivity is largest for r_0 and a_0 is around 10 MeV/u. And the sensitivity drops for higher energy.

Note that a_0 = 0.60 match the experimental data even better.

We plot the ratio

And we can see that, all energy range has the same % sensitivity for r_0 , but different sensitivity for a_0 .

Next, we fixed the energy at 10 MeV/u and study the sensitivity of a_0

We can see that,

\displaystyle \frac{\Delta a_0 }{\Delta \sigma_R} \approx \frac{0.05 \textrm{fm}}{ 200 \textrm{mb}}

So, if the measured cross section is smaller than 40 mb, the uncertainty of a_0 is about 0.01 fm at 10 MeV/u for 12C+12C. The last plot at below is the sensitivity for r_0 .

\displaystyle \frac{\Delta r_0 }{\Delta \sigma_R} \approx \frac{0.05 \textrm{fm}}{ 50 \textrm{mb}}

By measuring the reaction cross section at two energies with uncertainty less than 50 mb, if the model explain the data, the Woods-Saxon parameters should be extracted with small uncertainty.

I did not consider the uncertainty from the density. The error propagation from the density to the Woods-Saxon parameters need to be studied.

Woods-Saxon density in Glauber model

July 14, 2021

GoLuckyRyan Basic , Leave a comment

Following the receipt of the Glauber model. I did the calculation on 12C + 12C using Woods-Saxon density, and compared with the experimental result.

I did it using Mathematica. In order to avoid high dimension integration, the T_{A}(s), T_{B}(s) , T_{AB}(b), T(b) are first tablelized and interpolated. This reduced the complexity and speeded things up.

The NN cross section is taken from the reference S. K. Charagi and S. K. Gupta, PRC 41 (1990) 1610. It gives a parametrized from for the NN cross section. But I did not cross check the value for < 10 MeV.

I also compared the calculation with and without Coulomb correction. Also, finite-range interaction was also calculated.

Here is the results. The WS radius is R = 1.25 A^{1/3}, a = 0.67 .

The T_{AB}(b) in this case is very close to a Gaussian, and can be approximated as T_{AB}(b) = 0.02124 \exp\left( - \frac{b^2}{3.884^2} \right)

In order to match the 12C+12C cross section, I have to scale the T_{AB}(b) by 80 times.

The above is the 12C+12C result from Glauber model using zero-range. The calculation is close to the experimental result from the reference. And the low energy is the Coulomb correction. I guess my Coulomb correction is not correct? But for low energy, k^2 = 2m(K)/\hbar^2 should be correct, I am more less confident for the high energy. Anyway, the calculation is not that far away. Below is a figure taken from the reference for the 12C+12C reaction

Taken From S. K. Charagi and S. K. Gupta, PRC 41 (1990) 1610

I take the Coulomb KE to be half of the incident KE, i.e. using Center of mass KE. We get a better result at low energy. It seems that the center of mass KE should be used, but why? [update 2021-07-18] it is because the target will be pushed by the Coulomb force from the projectile. Since the net force on the two-particle system is zero, transform the coordinate into the CM frame, and the CM frame KE is used to calculate the closed distance.

Next I compare the calculation using a Gaussian finite-range with finite range parameter of 0.9 fm, below is the T_{AB}(b) :

The finite-range interaction thickness function takes a similar shape but only different amplitude. This is understandable as the NN interaction not only when the nucleon-nucleon at the same position, but also can interaction at 0.9 fm apart. Since the transparency function is insensitive (I will discuss this at below) for the inner part of the thickness function. The increase of the amplitude would not cause much different, but the increase width will make a different.

At the end, I like to add some comments on the calculation.

Frist, in the calculation, I used a normalized Woods-Saxon density. Is it necessary? Since the probability of a NN-collision at impact parameter b is T_{AB}(b) \sigma_{NN}, and \sigma_{NN} for low energy is 100 ~ 1000 mb or 10 ~ 100 fm2. This term T_{AB} \sigma_{NN} < 1 is often not fulfilled.

Also, the transparency function, if we take the original form

\displaystyle T(b) = (1 - \sigma_{NN} T_{AB}(b) ) ^{AB}

This requires T_{AB}(b) \sigma_{NN} < 1 , but as I discussed above it is usually not the case. However, if we use the exp form

\displaystyle T(b) = \exp \left(- \sigma_{NN} T_{AB}(b) \right)

then there seems to be no limitation for the thickness function T_{AB}(b) . So, if the normalization doesn’t matter, then what is the meaning of the thickness function times NN cross section?? not probability anymore? [update July-16, 2021] refer to the end of this post. The density function is normalized to mass number.

Second, I tried to follow the step in the reference S.K. Charagi and S.K. Gupta, but I cannot reproduce the Gaussian parameters as they did. In the paper, they match the Gaussian density to the surface of the Woods-Saxon function. The problem is the normalization of the Woods-Saxon density. I can use the table to interpolate the density and radius, and the Gaussian density is not normalized. Even I use their table values for the density and radius, I still unable to reproduce their result, except the 12C+12C, I can reproduce the trend, but the curve is still not the same as their.

Third, in the Coulomb correction, the wave number is non-relativistic.

Fourth, at energy below the Coulomb barrier, it becomes mostly elastic scattering, and the effective impact parameter, is very large, the reaction cross section should becomes zero. The Coulomb corrected Glauber model should reflect this.

Fifth, the interplay between the NN cross section, finite-range, the thickness function, and the transparency function. The transparency function is

\displaystyle T(b) = \exp \left(- \sigma_{NN} T_{AB}(b) \right)

We use a Gaussian approximation for the T_{AB}(b) for illustration. And the transparency function becomes

\displaystyle T(b) = \exp \left(- \alpha \exp\left(-\frac{b^2}{\beta^2}\right) \right)

This is the exactly the form from using Gaussian density.

In below, alpha = 4, \beta = 2 ,

We can see that, when \sigma_{NN} T_{AB}(b) > 4 , T(b) \sim 0 and \sigma_{NN} T_{AB}(b) = \ln(2) \rightarrow T(b) = 0.5 . We can see that, the increase of the NN cross section, the increase of the WS density, and the finite-range effect, has similar effect and all contribute to increase the width of \sigma_{NN} T_{AB}(b) . When the width increased, the position for T(b) = 0.5 also increased, and the integration or the reaction cross section increased.

Here is the calculation for the 40Ca + 40Ca, the WS density also 80 times more.

Sixth, the 1- T(b) looks like a Woods-Saxon function, can we find a simple relation between the density parameters to the 1-T(b) ?

[update July-16, 2021]

In the above, we have to multiple the T_{AB}(b) by 80 to reproduce the experimental result. At the calculation, the T_{AB}(b) is calculated using Woods-Saxon function that is normalized by 1. Thus, this 80 is the multiplication of the mass number. But the multiplication of the mass number is 12 \times 12 = 144 . And 80 is only 55%. This 55% could be a correction to the density in the compensation for the radius and diffusiveness parameters.. For Woods-Saxon potential that normalized to the mass number, r_0 = 1.25, a_0 = 0.67

We can see that the density has 20% difference, and we know that the density should be almost constant as the short-range nuclear force saturated at the center of nucleus. And the experimental nuclear density is 0.16 fm-3.

I redo the calculation using r_0 = 1.0, a_0 = 0.6 , and the scaling is 144.

The result is closer to the experimental result.

Gaussian density for Glauber model

July 9, 2021

GoLuckyRyan Basic , Leave a comment

We are going to calculate the reaction cross section using Glauber model with Gaussian density. This already done in reference S. K. Charagi and S. K. Gupta, PRC 41 (1990) 1610, But I will provide more detail of the calculation. And also, in the reference, there is a factor 10, that is the conversion 1 \textrm{fm}^2 = 10 \textrm{mb}

The density function is

\displaystyle \rho(r) = \rho_0 e^{-r^2/R^2}, \rho_0 = \frac{1}{(\sqrt{\pi}R)^3}

The finite-range function takes the form

\displaystyle f(r) = \frac{1}{\pi R^2} e^{-r^2/R^2}

when R \rightarrow 0 , the function becomes 2-D delta function .

The goal is to compute

\displaystyle \sigma^{NN} = 2\pi \int bdb \left( 1 - \exp(-\sigma_{NN} T_{AB}(b) \right)

\displaystyle T_{AB}(b) = \int d z_A d z_B d^2s d^2s' \\ \rho_A\left(\sqrt{s^2+z_A^2}\right) \rho_B\left(\sqrt{|\vec{s'}-\vec{b}|^2 + z_B^2}\right) f(|\vec{s}-\vec{s'}|)

The T_{AB} is a 6-dimension integral, but thanks to the Gaussian density and Gaussian finite-range function, all 6-dimension can be integrated individually in Cartesian coordinate.

Let integrate the beam direction.

\displaystyle \int_{-\infty}^{-\infty} \rho\left(\sqrt{s^2+z^2}\right) dz = \rho_0 \sqrt{\pi} R \exp\left( - \frac{s^2}{R^2} \right)

Thus,

\displaystyle T_{AB}(b) = \rho_1 \rho_2 \frac{R_1 R_2}{R_0^2} \int d^2s d^2s' \exp \left( - \frac{s^2}{R_1^2} \right) \\ \exp \left( - \frac{|\vec{s'}-\vec{b}|^2}{R_2^2} \right) \exp \left( - \frac{|\vec{s}-\vec{s'}|^2}{R_0^2} \right)

In Cartesian coordinate,

\displaystyle T_{AB}(b) = \rho_1 \rho_2 \frac{R_1 R_2}{R_0^2} \int dx_1 dy_1 dx_2 dy_2 \exp \left( - \frac{x_1^2+y_1^2}{R_1^2} \right) \\ \exp \left( - \frac{(x_2 - x_b)^2+(y_2 - y_b)^2}{R_2^2} \right) \exp \left( - \frac{(x_1-x_2)^2 + (y_1-y_2)^2}{R_0^2} \right)

Lets integrate x_2

\displaystyle \int_{-\infty}^{\infty} \exp \left( - \frac{(x_2 -x_b)^2}{R_2^2} - \frac{(x_1-x_2^2)}{R_0^2} \right) dx_2 \\ = \frac{\sqrt{\pi}R_0 R_2}{\sqrt{R_0^2+R_2^2}} \exp\left(-\frac{(x_1-x_b)^2}{R_0^2+R_2^2} \right)

Similarly for the y_2 , and we have

\displaystyle T_{AB}(b) = \rho_1 \rho_2 \frac{R_1 R_2}{R_0^2} \frac{\pi R_0^2 R_2^2}{R_0^2+R_2^2} \int dx_1 dy_1 \\ \exp \left( - \frac{x_1^2+y_1^2}{R_1^2} \right) \exp\left(-\frac{(x_1-x_b)^2}{R_0^2+R_2^2} -\frac{(y_1-y_b)^2}{R_0^2+R_2^2} \right)

We can see a similar pattern in here, rearrange and simplify

\displaystyle T_{AB}(b) = \rho_1 \rho_2 \pi \frac{R_1 R_2^3}{R_0^2+R_2^2} \int dx_1 dy_1 \\ \exp \left( - \frac{x_1^2}{R_1^2} -\frac{(x_1-x_b)^2}{R_0^2+R_2^2}\right) \exp\left(- \frac{y_1^2}{R_1^2}-\frac{(y_1-y_b)^2}{R_0^2+R_2^2} \right)

Thus,

\displaystyle T_{AB}(b) = \rho_1 \rho_2 \pi \frac{R_1 R_2^3}{R_0^2+R_2^2} \frac{\pi R_1^2 (R_0^2+R_2^2)}{R_0^2+R_1^2+R_2^2} \exp \left( -\frac{x_b^2+ y_b^2}{R_0^2+R_1^2+R_2^2}\right)

\displaystyle T_{AB}(b) = \rho_1 \rho_2 \pi^2 \frac{R_1^3 R_2^3}{R_0^2+R_1^2+R_2^2} \exp \left( -\frac{b^2}{R_0^2+R_1^2+R_2^2}\right)

When R_0 \rightarrow 0 , it reduces to zero-range.

We simplify T_{AB}(b) :

\displaystyle T_{AB}(b) = \alpha \exp\left( - \frac{b^2}{\beta^2} \right) = \frac{1}{\pi \beta^2} \exp \left( - \frac{b^2}{\beta^2}\right)

where

\displaystyle \beta^2 = R_0^2 + R_1^2 + R_2^2, \alpha = \rho_1 \rho_2 \pi^2 \frac{R_1^3 R_2^3}{\beta^2} = \frac{1}{\pi \beta^2 }

We substitute \rho_i = 1/(\sqrt{\pi}R_i)^3 at the last step. We can check that

\displaystyle 2\pi\int_0^{\infty} T_{AB}(b) b db = 1

which is expected as T_{AB}(b) is the probability of a NN-collision at impact parameter b.

The reaction cross section take the form

\displaystyle \sigma^{AB} = 2\pi \int_{0}^{\infty} \left( 1- \exp(-\sigma_{NN} \alpha e^{-\frac{b^2}{\beta^2}}) \right)b db

This is exponential of exponential, and integrate. I never deal with this before. Lets get a feeling on the \exp(-a \exp(-x^2) . When x is large, the inside exp approaches to zero, so, the outside \exp(-0) gives 1. When x is small, the inside exp gives a finite number, and the outside \exp(-|y|) < 1 . In fact, if we plot using \beta = 1, \sigma_{NN} = 2

We can see that, the T_{AB} is the probability of a NN-collision, and the \exp( -a T_{AB}) is the transparency function.

Mathematica gives the integral as:

\displaystyle \sigma^{AB} =\pi \beta^2 \left( \gamma + \ln(\alpha \sigma_{NN}) - E_i(-\alpha \sigma_{NN}) \right) \\= \pi \beta^2\left( \gamma + \ln\left(\frac{\sigma_{NN}}{\pi \beta^2}\right) - E_i\left(-\frac{\sigma_{NN}}{\pi \beta^2}\right) \right)

where \gamma = 0.577216 is the Euler Gamma constant and E_i(z) is the exponential integral function.

\displaystyle E_i(z) = - \int_{-z}^{\infty} \frac{e^{-t}}{t} dt

The expression agree with the reference. We plot the \sigma^{AB} below,

For a fixed \beta, i.e. the size of the target and projectile, the reaction cross section change with \sigma_{NN} non-linearly. But for a fixed \sigma_{NN} the reaction cross section saturated to \sigma_{NN}.

In fact, mathematically,

\displaystyle \lim_{\beta \rightarrow \infty} \sigma^{AB} = \sigma_{NN}

For some I don’t know why reason, the radial distance R with

\displaystyle \exp\left( -\sigma_{NN} T_{AB}(R) \right) = 0.5

is called the strong absorption radius. In here, the density normalization is to the mass number, not 1.

\displaystyle \exp\left(\sigma_{NN} T_{AB}(R) \right) = 2

\displaystyle \frac{\sigma_{NN}}{\pi \beta^2} \exp\left( - \frac{R^2}{\beta^2} \right) = \chi_0 \exp\left( - \frac{R^2}{\beta^2} \right) =\ln(2)

\displaystyle - \frac{R^2}{\beta^2} =\ln(\ln(2)) - \ln(\chi_0)

\displaystyle R^2 =\beta^2 \left( \ln(\chi_0) - \ln(\ln(2)) \right)

R^2 > 0 requires \chi_0 > \ln(2) = 0.693147 or

\displaystyle \sigma_{NN} AB > \ln(2) \pi \beta^2 = \ln(2) \pi (R_1^2 +R_2^2 + R_0^2 )

When the energy is more than 100 MeV, the \sigma_{NN} \sim 40 \textrm{mb} = 4 \textrm{fm}^2 , this impiles \beta^2 < 1.8 A B ~ \textrm{fm}^2 , which is satisfied for most nuclei.

Then, a classical geometrical cross section is

\displaystyle \sigma^G = \pi R^2 = \pi \beta^2 \left( \ln\left( \frac{\sigma_{NN}}{\pi \beta^2} \right) - \ln(\ln(2)) \right)

Compare with

\displaystyle \sigma^{AB} = \pi \beta^2\left( \ln\left(\frac{\sigma_{NN}}{\pi \beta^2}\right) + \gamma - E_i\left(-\frac{\sigma_{NN}}{\pi \beta^2}\right) \right)

RMS radius of various nuclear density functions

July 8, 2021

GoLuckyRyan Basic , Leave a comment

In the Glauber model, one of the input is the density functions. People will adjust the parameters of the density function to fit the experimental cross section.

There are many choice of nuclear density functions, for example, Woods-Saxon, Harmonic Oscillator, Square well, and Gaussian etc. Among all of them, the root-mean-squared (RMS) radius of these functions usually represent the nuclear RMS radius. So, there are many RMS radius depends on the density function. In here, I will make some comparisons.

All density functions are normalized

\displaystyle 4\pi \int_{0}^{\infty} f(r) r^2 dr = 1

and the mean-squared radius is

\displaystyle \left<r^2\right> = \frac{\int_{0}^{\infty} f(r) r^4 dr}{\int_{0}^{\infty} f(r) r^2 dr } = 4\pi \int_{0}^{\infty} f(r) r^4 dr

Using the above definition, I have

NameDensity functionRMS radius
Square Well\frac{3}{4\pi R^3}, 0<r<R \sqrt{\frac{3}{5}}R
Woods-Saxon\frac{1}{-8 a^3 \pi \textrm{PolyLog}(3, -e^{R/a})} \frac{1}{1+e^{(r-R)/a}} a \sqrt{12\frac{\textrm{PolyLog}(5, -e^{R/a})}{\textrm{PolyLog}(3, -e^{R/a})}}
Gaussian\frac{1}{\sqrt{\pi^3}R^3} \exp\left(- \frac{r^2}{R^2}\right) \sqrt{\frac{3}{2}}R

The rms radius of Square well and Gaussian is simple ratio to the radius parameter.

As we shown before, the Woods-Saxon will approach to Square well, as it should be. But the Gaussian will be different. How to translate the rms radius between each of them? For example, If I say the rms radius is 4 fm +- 0.1 fm with Gaussian density function, what is the rms radius for Woods-Saxon ?

Here is other comparisons

I guess, to translate the rms radius between different density functions in Glauber model, we have to compare the T_A(s) . If the T_A(s) is similar. To define the similarity of functions, we using

\displaystyle 4\pi \int \sqrt{f(r) g(r)} r^2 dr

When f(x) = g(x) , the integral is equal to 1.

I compare the Woods-Saxon and Gaussian,

The 3rd order approximation for the peak (red) line is

\displaystyle 1.80144 + 0.212716 x + 0.0444465 x^2 - 0.00193038 x^3

Here is the plot for the Woods-Saxon with latex R=4 , a = 0.67$, and the Gaussian with R = 3. 24 ,

Finite-range Glauber Model and Coulomb correction

July 7, 2021

GoLuckyRyan Basic Leave a comment

In the past two post, the NN interaction/collision takes place when the two nucleons at the same position. In order to account the finite size of nucleon and also the finite range of the strong force, a finite-range correction is needed. Under the correction, the T_{AB}(b) becomes,

\displaystyle T_{AB}(b) = \int dz_A \int dz_B \int s ds \int d\theta \int ds' f\left(|\vec{s}-\vec{s'}|\right) \\ \rho_A\left(\sqrt{s^2 + z_A^2}\right) \rho_B\left( \sqrt{s'^2 + b^2 -2s'b \cos\theta + z_B^2} \right)

The finite-range function f\left(|\vec{s}-\vec{s'}|\right) is introduced with

\displaystyle \int f(s) ds = 1

We can see that, when the finite-range function becomes a delta function, \vec{s} = \vec{s'} and it restored to the “old” Glauber model.

The origin Glauber model is for high energy collision, that the nucleons motion are independent and linear, the Coulomb force can be neglected. But in low energy, Coulomb force will distorted the motion and the effective impact parameter before the collision will be different from the impact parameter at infinity.

Taken from S. K. Charagi and S. K. Gupta, PRC 41 (1990) 1610

In non-relativistic, the wave number

\displaystyle k = \frac{\sqrt{2mE}}{\hbar} = \frac{p}{\hbar} = \frac{mv}{\hbar} = \frac{2\pi}{\lambda}

The Coulomb force is a central force, it conserves the angular momentum. At infinity and at the closed point, the momentum must be perpendicular to the impact parameter.

\displaystyle p b = p'b'

And by conservation of energy

\displaystyle E = E' + \frac{Z_A Z_B e^2}{4\pi \epsilon b'}

Solve the above equations and using the Summerfeld parameter

\displaystyle \eta = \frac{Z_AZ_B e^2}{4\pi \epsilon \hbar v}

we have

\displaystyle kb' = \eta + \sqrt{\eta^2 + k^2 b^2 }

It is easy to verifly at b = 0, b' = 2\eta/k .

———————- [update 2021-07-18]

In the center of momentum frame, the total KE is

\displaystyle T_{cm} = \frac{m_1 m_2}{m_1+m_2} K

, where K is the MeV/u in Lab frame of the projectile. The expression of KE is for non-relativistic. The modified Summerfeld parameter is

\displaystyle \eta' = \frac{Z_1 Z_2 1.44}{ 2 T_{cm} }

And the Coulomb correction is

\displaystyle b' = \eta' + \sqrt{\eta'^2 + b^2}

———————————————-

The Coulomb correction is replacing b \rightarrow b' , the T_{AB}(b) becomes

\displaystyle T_{AB}(b') = T_{AB}\left( \frac{\eta + \sqrt{\eta^2 + k^2b^2}}{k} \right)

and the reaction cross section, using the transparency function is

\displaystyle \sigma^{AB} = 2\pi \int b db ( 1- T(b') )

where the transparency function is

\displaystyle T(b) = \exp\left( -T_{AB}(b) \sigma_{NN} \right)

Glauber Model – 2

July 6, 2021

GoLuckyRyan Basic , Leave a comment

I compared the calculations from hard sphere and Woods-Saxon type density functions, to study the sensitivity of the Glauber model. The Glauber model calculate, for spherical nuclei

\displaystyle T_{AB}(b)= \int d z_A \int d z_B \int d\theta \int sds \\\rho_A\left(\sqrt{s^2+z_A^2}\right) \rho_B\left(\sqrt{s^2 + b^2 - 2 sb \cos\theta + z_B^2}\right)

which is a 4-dimension integration. Once the T_{AB}(b) is calculated, all other things follow. The reaction cross section is

\displaystyle \sigma^{AB} = 2\pi \int b db \left( 1 - \left(1-T_{AB}(b) \sigma_{NN}\right)^{AB} \right)

When the T_{AB}(b) \sigma_{NN} is small and AB is large, make analog,

\displaystyle \lim_{n\rightarrow \infty} \left( 1+ \frac{x}{n} \right)^n = e^x

an approximation:

\displaystyle \sigma^{AB} = 2\pi \int b db \left(1 - e^{- T_{AB}(b) \sigma_{NN}} \right)

and call the exponential be transparency function T(b)

\displaystyle T(b) = e^{- T_{AB}(b) \sigma_{NN}}

Because the form of the reaction cross section is also like

\displaystyle \sigma_a = \frac{\pi}{k^2} \sum (2l+1) (1-|S_l|^2)

And the meaning of |S_l|^2 is also the transparency. Also, in semi-classical, the average angular momentum is l+1/2 = k b , Thus, people would link

\displaystyle T(b) \approx S_l .

The hard sphere is relatively less time consuming, and the density function is easy to play with. It is easy to know that

\displaystyle \rho(r) = \frac{3}{4\pi R^3} , 0<r<R,

the corresponding T_A is

\displaystyle T_A(s) = \int_{-R}^{R} \rho(\sqrt{s^2 +z_A^2}) dz_A = \frac{3}{2\pi R^3}\sqrt{R^2-s^2}, 0<s<R

But it is hard to get the analytical form of the T_A for Woods-Saxon, so I have to integrate the density numerically using 4-dimensional 25-points Gauss quadrature. Because it is 4-dimensional, the number of sum is 25^4 = 390625 , comparing to 2-dimensional 25-point GQ, which only need 25^2 = 625 .

The normalized Woods-Saxon function is

\displaystyle WS(r) = \frac{1}{4\pi} \frac{1}{-2a^2 PolyLog(3, -e^{R/a})} \frac{1}{1+e^{(r-R)/a}}

Here is the comparison for radius 2 and radius 1, and diffusiveness is 0.3.

The simple hard-sphere means that I use 2-dimension GQ, “hard-sphere” means the 4-D GQ was used. And that verified the calculations for 4-dimension GQ works well.

The Woods-Saxon gives a smooth curve because the function itself is smooth. In Hard-sphere, the function is not smooth and that could create problem when using GQ method.

If we interpolate the two functions and integrate, i.e.

\displaystyle \int T_{AB}(b) d db

The hard-sphere gives 0.160056 and WS gives 0.158891. So, I think I can conclude that the Glauber model is not sensitive in these two cases. i.e., If someone try to extract the nucleus radius by fitting the experimental reaction cross section using Glauber model, using hard-sphere and WS is basically the same.

[Update July-16, 2021]

The “approximation” of \sigma^{AB} mentioned above is not correct.

\displaystyle \sigma^{AB} = 2\pi \int b db \left( 1 - \left(1-T_{AB}(b) \sigma_{NN}\right)^{AB} \right)

For \sigma_{NN} T_{AB}(b) < 1 ,

\displaystyle \sigma^{AB} = 2\pi \int b db \left( 1 - \exp\left(1-AB T_{AB}(b) \sigma_{NN}\right) \right)

And we can absorb the AB into T_{AB}(b) and this lead to re-define of the density function.

\displaystyle \int \rho_A(r) d^3 r = A

That means, the density function is no longer normalized to 1, but it normalized to equal the number of nucleons.

If we assume the density function is normalized to 1 and takes a Woods-Saxon shape with parameters r_0 = 1.25, a_0 = 0.67 , the maximum \rho_z(0, A=1) = 0.082 \textrm{fm}^{-2}. And this limited the NN cross section be \sigma_{NN} < 12.2 \textrm{fm}^2 = 122 \textrm{mb} . And that translate to E_{Lab} >\approx 100 \textrm{MeV}/u . This estimation is only based on one nucleus. T_{AB}(0) = \int \rho_z(A) \rho_z(B) d^2s < \int \rho_z^2(\textrm{min}(A,B)) d^2s , We calculated the energy limit versus mass number using

\displaystyle T_{AB}(0) = \int \rho_z^2(A) d^2s

\displaystyle \sigma_{NN}(E) T_{AB}(0) = 1

We can see that for most of the nuclei, the lab energy can be as low as 4 MeV/u (green line). And that is close to the Coulomb barrier.

Glauber Model for reaction cross section

July 5, 2021

GoLuckyRyan Basic Leave a comment

Before we go to the Glauber model, lets review the cross section from microscopic point of view.

Suppose we have a target at the center, the target density can be non-uniform, the number of target in area A and thickness is

\displaystyle N_t = \int_A{\rho(r) dt rdrd\theta}

we have a beam, the beam is not uniform and has some distribution, the total number of beam particle per unit area is

\displaystyle \frac{dN_B}{dS} = f_B(r)

and the total number of beam and target in area A and S, or, the Luminosity [number of beam and target per unit area],

\displaystyle L(A) = \int_A {f_B(r) \rho(r) dt rdrd\theta}

The cross section for incident energy is \sigma(E) . With impact parameter b and the scattering angle \theta

\displaystyle \sigma(E) = \int \frac{d \sigma(E, b)}{db} db^2 = \int \frac{d\sigma(E,\theta)}{d\theta} d\theta

The yield for the reaction is

\displaystyle Y(E) = L(A) \sigma(E) = \int {\int_A {f_B(r) \rho(r) dt r dr d\theta}} \sigma( E)

This is the most generic expression of the reaction yield for non-uniform target density and beam profile. To be more general, we should include the impact parameter for each collision, but we assumed the number of target and beam are large enough that all impact parameters will be exhausted and the result is a simple cross section.

Usually, we assume a uniform target thickness, so, integration on t becomes multiplication. And we have a very big target, much bigger than the beam size, so the integration or area will be defined by the beam profile function. If not, the area is defined by the product of \rho(r) f_B(r) .

The Nobel prize (2005) winner, Roy J. Glauber, proposed a theoretical model to describe nuclear reaction in the 1950s.

The scattering problem involve the impact parameter, number of nucleons (both target and projectile), and the number of NN-collisions.

The Glauber model has two input the nuclear density and the inelastic nucleon-nucleon cross section. The nuclear density function usually take the Woods-Saxon shape.

\displaystyle \rho(r) = \rho_0 \frac{1}{1+\exp(\frac{r-R}{a})}

The Glauber model assumed

  • high incident energy that the nucleons are un-deflected.
  • the nucleons are independent

The following is taken from M. L. Miller et al., Annu. Rev. Nucl. Part. Sci. 2007. 57:205-43.

Suppose the target is nucleus A and the projectile is nucleus B. The impact parameter b.

Taken from M. L. Miller et al., Annu. Rev. Nucl. Part. Sci. 2007. 57:205-43

For a nucleon, the probability per unit perpendicular (to the beam direction) area can be found at location \vec{s} in target A is

\displaystyle T_A(\vec{s}) = \int \rho_A(\vec{s}, z_A)dz_A

the integration of the density \rho_A should give unity

\displaystyle \int \rho_A(\vec{s}, z_A) dz_A ds^2 = 1

Similarity, the probability of a nucleon appears at location \vec{s} - \vec{b} in projectile B is

\displaystyle T_B(\vec{s}-\vec{b}) = \int \rho_B(\vec{s}-\vec{b}, z_B)dz_B

Thus, the total probability of a nucleon at A and an other nucleon at B at location \vec{s} is

\displaystyle T_A(\vec{s}) T_B(\vec{s}-\vec{b}) ds^2

Integrate all possible position, the probability with impact parameter b is

\displaystyle T_{AB}(\vec{b}) = \int T_A(\vec{s}) T_B(\vec{s}-\vec{b}) ds^2

T_{AB} is also called as thickness with dimension of inverse area.

The probability of a NN collision is the thickness times nucleon-nucleon cross section \sigma_{NN} . Since the incident energy is large, the dominated cross section is from the inelastic collision.

The probability of n N-N interaction at impact parameter b, with normalized thickness, is

\displaystyle P(n,\vec{b}) = \begin{pmatrix} AB \\ n \end{pmatrix} (T_{AB}(b) \sigma_{NN})^n (1-T_{AB}(b) \sigma_{NN})^{AB-n}

The total probability of the reaction is

\displaystyle \frac{d^2\sigma^{AB}}{db^2} = p^{AB}(b) = \sum_{n=1}^{AB} P(n,b) = 1- (1-T_{AB}(\vec{b}) \sigma_{NN})^{AB}

For non-polarized beam and target, the total cross section is

\displaystyle \sigma^{AB} = 2\pi \int{ bdb (1-(1-T_{AB}(b)\sigma_{NN})^{AB})}

The average number of collisions is

\displaystyle N_{col}(b) = \sum_{n=1}^{AB} n P(n,b) = AB T_{AB}(b) \sigma_NN

The number of participants (or wounded nucleons) is

\displaystyle N_{part}(\vec{b}) = A \int T_A(\vec{s})\left( 1 - \left( 1- T_B(\vec{s} -\vec{b}) \right)^B \right) ds^2 + \\ B \int T_B(\vec{s} -\vec{b})\left( 1 - \left( 1- T_A(\vec{s}) \right)^A \right)

The above is the Glauber model

Now, we are going to do some calculation using the model. Lets use a hard-sphere

\displaystyle \rho_A(\vec{r}) = \rho_B(\vec{r}) = \frac{3}{4\pi R^3}, |r| < R

Using simple geometry,

\displaystyle T_A(\vec{s}) = T_B(\vec{s}) = 2 \sqrt{R^2 - |\vec{s}|^2}, |\vec{s}| < R

The thickness is

\displaystyle T_{AB}(\vec{b}) = \int T_A(\vec{s}) T_B(\vec{s}-\vec{b}) ds^2

without loss pf generality, set \vec{b} = (b, 0)

\displaystyle T_{AB}(b) = 4 \int_{0}^{2\pi} \int_{0}^{3R} \sqrt{R^2 - s^2} \sqrt{R^2-s^2 -b^2 + 2sb \cos\phi} sdsd\phi

The above integral is the overlap area of 2 circles.

Using 30-points Gauss Quadrature, for R=1, we have:

The orange line is a Gaussian fit with \sigma = 0.687717. Since it is 2 hard spheres with radius 1, when b > 2 the overlap should be zero.

In Mathematica, the calculation code is: 

(* Define Gauss Quadrature *)
GQ2[f_Symbol, n_, rangeX_, rangeY_] := Module[{a, W, x},
  a = x /. NSolve[LegendreP[n, x] == 0, x];
  W = Table[
    2/((1 - a[[i]]^2) Evaluate[D[LegendreP[n, x], x]]^2 /. {x -> 
       a[[i]]}), {i, 1, n}];
  (rangeY[[2]] - rangeY[[1]])/2 (rangeX[[2]] - rangeX[[1]])/
   2 Sum[W[[i]] W[[
      j]] f[(rangeX[[2]] - rangeX[[1]])/2 a[[i]] + (
       rangeX[[2]] + rangeX[[1]])/
       2, (rangeY[[2]] - rangeY[[1]])/2 a[[j]] + (
       rangeY[[2]] + rangeY[[1]])/2], {i, 1, n}, {j, 1, n}]
  ]

 (* define TA *)
\[Rho][r_, R_] := Piecewise[{{2 Sqrt[R^2 - r^2], 0 < r < R}}, 0]

(* calculate TAB for given b *)
Data = Monitor[Table[
    g[r_, \[Theta]_] := \[Rho][r, 1] \[Rho][Sqrt[
       r^2 + b^2 - 2 r b Cos[\[Theta]]], 1] r ;
    {b, GQ2[g, 30, {0, 15}, {0, 2 \[Pi]}]},
    {b, 0, 15, 0.1}], b];

OK, now I can calculate the T_{AB}(b), then the total probability can be calculated. And then integrate the impact parameter, we should have the total reaction cross section.

In summary, the input of the Glauber model is the nucleon density function \rho(r) , and the NN cross section \sigma_{NN}. The NN cross section and total reaction cross section can be measured experimentally, thus, we can adjust the parameters in the nucleon density function to fit the data.

If using Woods-Saxon function as the density function, how sensitive is the Glauber model?

Paper reading : One-neutron halo structure in 15C

April 7, 2021

GoLuckyRyan papers reading , , Leave a comment

The paper DOI is https://doi.org/10.1103/PhysRevC.69.034613

This is my personal understanding of the paper.

Reaction : 14,15C one- or two- neutron removal on a carbon target (~ 375 mg/cm2)

Energy : 51 or 83 MeV (at the middle of the target), the 14,15C was produced using 110 A MeV 22Ne on 2 mm thick Be target.

Measurement : Reaction cross section and longitudinal momentum of carbon fragments (15C,14,13C) and (14C,13C) in the beam rest frame.

Result: the widths of the longitudinal momentum from (15C,14C) and (15C,13C) are 71 ± 9 MeV/c and 223 ± 28 MeV/c respectively, the width for that of (14C,13C) is 195 ± 21 MeV/c.

  • Due to the small neutron separation energy (Sn = -1.218 MeV) and the ground state spin-parity of 1/2+, 15C is thought to be a one-neutron halo.
  • The ground state spin parity of 1/2+ suggests that the valence neutron is in the 1s1/2 orbital.
  • 14C(d,p) reactions measured the spectroscopic factor for the s-wave around around 0.76 to 1.03, depending on the DWBA setting.

The authors motive is follow:

  • momentum distributions of (15C,14C) at beam energy ~ 100 MeV/u are measured at NSCL, and the result widths are smaller than the Glodhaber’s model. That suggests a larger spatial distribution, which is a characteristic of halo.
  • Another characteristic of halo is enhanced reaction cross section
  • The reaction cross section on carbon target (945 mb) at ~740 A MeV for 15C showed no enhancement, but the cross section (853, 862, 880, 1036, 1056, 1104, 1231, 1187 mb) for 12,13,14,16,17,18,19,20C are measured at ~950 A MeV [A. Ozawa, NPA, 693 (2001) 32]
  • The reaction cross section measured at ~ 60 MeV/u show some enhancement. [D.Q.Fang, PRC, 61 (2000) 064311]
  • The charge exchange cross section at 930 MeV/u show enhancement.

The authors stated that, low energy reaction is more sensitive to the tail of the nuclear density distribution than the high energy. They decided to measure the reaction cross section and the momentum distribution simultaneously at 51 and 83 MeV/u to find a consistence picture of the halo structure in 15C.

My comment : the high energy probe cannot see the enhancement is probably because of the probe is sensitive to all nucleus, or as the author suggested not sensitive to the tail. The high energy probe has de Broglie wave length shorter than nucleon can also able to removal the 0s-state. The 1s neutron only 1/9 of the all neutrons, so the effect of it is overwhelm by the other neutrons.

Also, may be I am not fully understood, I feel that the cross section is somewhat beyond intuition, for example, in every text book, the elastic cross section is 4 times bigger than the geometry cross section of the actual hard sphere, how and why ? I mean, I understand the math, but how to understand it in physical sense? And for the thermal neutron energy, the absorption cross sections for various elements do not make any sense. (do they?) I think the measured reaction cross section, not only depended on the nuclei size and the projectile de Broglie wave length, but also the nuclear interaction, and that makes the interpretation difficult and ambiguous.

The 1st part of the experimental result is the longitudinal momentum distribution. The widths of the longitudinal momentum from (15C,14C) and (15C,13C) are 71 ± 9 MeV/c and 223 ± 28 MeV/c respectively, the width for that of (14C,13C) is 195 ± 21 MeV/c.

My comment : The result could be well understand as the 1s1/2 spatial distribution, which is probed by the (15C,14C) alone, is larger than that of the p-wave, which is probed by (15C,13C) and (14C,13C). The (15C,13C) reaction should probe the mixed s- and p-wave distribution.

From the (15C,14C) distribution, they gave the s- and d- wave component of 66% and 4%, respectively, by assuming a fixed ratio between the s- and p- wave components. And I guess the remaining 30% is from the p-wave. This indicate the d-wave component is small, which is agreed with gamma-ray spectroscopy, that gave the d-wave is 2%.

The 2nd part of the experimental result is the reaction cross section. The result show a clearly enhancement on 15C at 81 MeV/u.

The Interesting part of the paper is the study of the relation between beam energy and reaction cross section using the finite-range Glauber model. The authors doing in this manner:

  • A harmonic-oscillator (HO) type density distribution were used to reproduced the cross section at relativistic energy (~ 950 meV/u)
  • The reaction cross section for lower energy was calculated using the high energy fitted parameter.
  • adding a p1/2 neutron tail for 14C, then fit the model with all energies.
  • Also tried HO-type + Yukawa tail, where the Yukawa tail gives a longer tail.

The study show that the Yukawa tail gives a better agreement to the 14C reaction cross section, although the other two methods also give similar result.

Then, they applied similar study on 15C. They found that a much longer tail will give much closer to the experimental reaction cross section.

In both 14C and 15C, the HO-type density ended at around 4 to 5 fm at density of 10-3 fm-3. For 14C, the tail extended to 12 fm at density of 10-7 fm-3, while for 15C, the tail extended to 20+ fm at density of 10-7 fm-3.

My comment : I am curious that can they reproduce the single-particle energy using the HO-type + Yukawa tail density. As in this post, I fitted the single-particle energy of 15C using ordinary Woods-Saxon potential, and also give the spatial and momentum distribution of the wave-function. The width (FWHM) of the momentum distribution for the 1s1/2 orbital is only ~ 60 MeV/c, smaller than 70 MeV/c. The p-wave width is about 180 MeV/c. Using 70% 1s-wave and 30% p-wave from the experiment, the simple mean 60 \times 0.7 + 180 \times 0.3 = 96 MeV/c, then the width will be larger than the experiment.

I guess one fundamental issue is, how the HO-type + Yukawa tail (or a Woods-Saxon) correctly describe the nuclear density for halo nuclei?