Study of proton-proton Scattering using Phase Function Method

Shikha Awasthi¹, Anil Khachi² and O.S.K.S. Sastri¹ ¹Department of Physics

\&

Astronomical Sciences
Central University of Himachal Pradesh, Dharamshala 176215, India.
Email: [email protected]
² Department of Physics, St. Bedes College, Himachal Pradesh, India

Abstract

Background: The study of $np$ and $pp$ scattering, central to understanding nuclear force, still remains an optional topic in many undergraduate nuclear physics curriculum.
Purpose: The main thrust of this paper is to study $pp$ scattering using phase function method to obtain the observed S-wave phase shifts and cross section at various energies.
Methods: The $pp$ interaction has been modeled by choosing Malfliet-Tjon potential for nuclear part along with screened Coulomb potential. The phase equation has been solved to obtain scattering phase shifts using fourth order RK-method (RK-4).
Results: The interaction potential obtained from optimised parameters matches well with the realistic Argonne V18 potential for ${}^{1}S_{0}$ state of $pp$ scattering and the scattering phase shifts as well as the cross section for energies ranging from $1-350$ MeV have been found to be in good agreement with expected data.
Conclusion: Introducing phase function method for S-wave ( $\ell$ =0) could bring this interesting study of nucleon-nucleon scattering to the undergraduate classroom.

March 2024

Keywords: Proton-Proton scattering, Scattering Phase Shifts (SPS), Phase function method (PFM), RK-4 method, Gnumeric/Excel, Scilab.

1 Introduction

Study of scattering in very light nuclei such as nucleon-nucleon (N-N), nucleon-nucleus (N-n) and nucleus-nucleus (n-n), is an important way to obtain understanding of nuclear force and gain insight into nuclear structure. Nuclear force is hypothesized to be charge independent. That is, the interaction strength involving neutron-proton and neutron-neutron is similar or equivalent to that of proton-proton, excluding its Coulomb repulsion. Hence, it is important to explore the topic of neutron-proton ( $np$ ) and proton-proton ( $pp$ ) scattering through the observed expected phase shift data phase shifts and cross sections at various energies to confirm or refute the hypothesis. The scattering phase shift (SPS) is a vital quantity to understand the nucleon-nucleon scattering. It signifies the changes in phase of the scattering wave occuring as a result of the interaction betwen the nucleons. Proper estimation of the scattering phase shift is crucial to compare experimental results with theoretical predictions and to understand the true nature of nuclear forces.
Typically, $np$ scattering is discussed even today, by modeling the interaction to be constant over the inter-nuclear distance as in square well potential [1], even though Yukawa [2] got Noble prize for his interaction model over 80 years ago. The reason being that it is easy to solve time independent Schr $\ddot{\text{o}}$ dinger equation (TISE) analytically for square well as opposed to lack of one, even though very involved, for the Yukawa potential till recently. Although Yukawa potential plays a fundamental role in explaining the strong nuclear force, it is unable to describe the tensor and spin-orbit components of the nuclear force. Because at low energies, the attractive component of nuclear force dominates [1] and the Yukawa potential provides an adequate description of nuclear interactions.
But at higher laboratory energies, the repulsive core; which results from the strong interaction at extremely short distances; becomes more significant due to the Pauli exclusion principle. Also, at higher laboratory energies, the interacting particles have more kinetic energy which allow them to overcome the attractive forces and come into closer proximity with each other[3]. At such small distance, the repulsive part of the potential becomes very important and plays important role in reduction of scattering cross section at large scattering angles. Therefore it becomes necessary to include both the repulsive and attractive parts of the nuclear force at higher energies, in order to describe nuclear interactions accurately.
This can be achieved by considering complex potential models that contain both attractive and repulsive components of the interaction potential. While potential models like the Woods-Saxon potential and realistic nucleon-nucleon potentials derived from Effective Field Theories (EFT) provide valuable insights into the behavior of light nuclear systems, it’s important to note that their applicability and accuracy may be limited when extended to heavier nuclei [3]. This can be due to the difference in structural properties, shell effects, effective interactions, and the limitations imposed by available experimental data as compared to heavier nuclei. Therefore, specialized models and approaches may be necessary to accurately describe the behavior of light nuclei.
There are many realistic nucleon-nucleon (NN) potentials available which are capable of properly describing the interaction between the nucleons [4, 5]. But in order acquire such an accuracy, they need to adjust large number of parameters (up-to nearly 40 parameters), to calculate the scattering phase shift accurately. It is quite challenging for an undergraduate or graduate student to calculate such a large number of parameters in order to calculate the scattering phase shift. Malfliet-Tjon potential [6] on the other hand is a three parameter potential, which is based on a non-relativistic approach to study nucleon-nucleon scattering, and provides a comprehensive description of nucleon-nucleon interactions, including both central and tensor components.
Recently, we have considered both Yukawa and Malfliet-Tjon potentials to explain the $np$ interaction [7] for laboratory energies up-to 350 MeV. We observed in our work [7] on $np$ scattering, that although Yukawa potential is a very good theoretical model used to describe the nuclear force between nucleons for low energy scattering. But at large value of incident projectile energy i.e above 50 MeV the Yukawa potential is not able to predict phase shifts for neutron-proton scattering. Malfliet-Tjon potential on the other hand is a modified Yukawa potential with repulsive term and is able to explain the interaction properly above 50 MeV as well.
Most often the scattering phase shifts are obtained from the wave function [3, 8, 9, 10] using R-matrix method [11], S-matrix method [12], Jost function method [13] etc. The cross section and phase shift calculations present in textbooks (if any) are only by considering wave function approach only. Even the formulae to obtain scattering phase shifts given in standard text book by Krane [1] is far beyond the scope of under graduate and post graduate courses. So, it is important to provide a simple procedure to calculate scattering phase shifts which will be helpful for under graduate and post graduate students.
We have alternatively preferred to use much simpler method known as phase function method (PFM) [14, 15], wherein second order time independent Schr $\ddot{\text{o}}$ dinger equation for wave function is cast into a first order Ricatti equation for phase shifts, called phase equation. The interaction potential and scattering energy in center of mass frame are given as input to this phase equation. It is solved numerically to obtain scattering phase shifts which can in turn be utilised to determine scattering cross sections.
Proton-proton (pp) scattering which is the scattering of identical nucleons is also an essential topic in nuclear and particle physics course. Although the proton-proton scattering is mentioned in some standard text books [1, 8, 9, 10] at under graduate and post graduate level. However, barely any of the books describe the methods and procedure for calculating scattering phase shifts. The formulae given in standard text book [1] to calculate SPS is far beyond the scope of under graduate and post graduate course.
In current work, we are studying proton-proton $pp$ interaction by considering Malfliet-Tjon (MT) potential [6] as interaction model. Proton-proton scattering is very important since it is the interaction which involves the interaction of two fundamental forces of nature i.e the electromagnetic and the nuclear forces, at subatomic scale. Studying $pp$ interaction using phase function method is extremely challenging due to the presence of long range Coulomb interaction alongside the short range nuclear interaction [16, 17, 18, 9]. Hence, one replaces the bare Coulomb potential with a screened one [17, 18], considering the observation that an isolated charge is typically surrounded by residual particles due to polarisation in actual experiments. Proton-proton scattering has also been studied using different potentials e.g., square-well potential [1]; to calculate cross section; and Yukawa potential [19] to calculate scattering parameters and cross section.
In this paper, our goal is to study $pp$ scattering by utilising Malfliet-Tjon (MT) potential as model for nuclear interaction and atomic Hulthén potential [17, 18] as an ansatz for screened Coulomb interaction. Using this combination, we have solved the phase equation numerically to obtain scattering phase shifts (SPS) and cross sections of $pp$ scattering for energies rannging from $1-350$ MeV. The detail of the implementation have been discussed in step by step simulation methodology for pedagogic purpose.

2 Methodology

Methodology is necessary to model the physical system under consideration and to simulate the problem [20], and it provides a step-by-step procedure for obtaining the solution. The simulation methodology [20] can be divided into four steps described below:

1.

Modelling the physics problem: This involves description of objects, interactions and process followed by formulation of a mathematical model, which typically consists of a dynamical equation with interaction embedded into it and specification of various initial and boundary conditions [21, 22].
2.

Preparation of the system for numerical solution: Here, one rephrases the mathematical model into an appropriate choice of units emicable for numerical implementation and sets up the revised iterative equations based on the chosen numerical technique.
3.

Implementation of the numerical method in a computer: By starting with an algorithm, one has to write a code in a chosen software and optimise the algorithm parameters to obtain a working code for simulation.
4.

Simulation and Discussion of Results: The system is studied by varying its physical parameters and one discusses their effect in understanding the underlying physics of the problem.

Now, we detail the above four steps for studying proton-proton scattering.

2.1 Modeling the physical system

2.1.1 Description stage:

The type of objects involved in a scattering system in the laboratory consists of an incoming projectile, that is accelerated to a certain energy, made to collide with a target that is at rest. In this case, our system is composed of a proton beam, accelerated to lab energies in the range of 1-350 MeV, to be incident on to a target consisting of liquid Hydrogen which can be thought of as collection of protons (mass $m_{p}$ , spin $1/2$ ). One can neglect the interaction due to electrons which are typically too light to contribute.
The interactions between the incoming protons and protons at rest are of two types:

1.

Short range attraction due to nuclear force and
2.

Long range repulsion due to Coulomb force

Even though, originally Yukawa potential was successful in explaining the nuclear force at low energies, later it was realised that an extra repulsive term similar to that of Yukawa needs to be added to account for observations at higher energies.
The process of scattering is more easily understood theoretically, by transofrming the two body system in lab frame to one that consists of a single particle in center of mass system. The transformation involves position vectors of two particles $m_{1}$ and $m_{2}$ with $(\vec{r_{1}},\vec{r_{2}})$ to be related to $(\vec{r},\vec{R})$ , where the first co-ordinate is the relative distance between the two particles and $\vec{R}$ is the center-of-mass (CoM) coordinate. This reduces the two body problem to a single particle system with a mass given by

\mu=\frac{m_{1}m_{2}}{m_{1}+m_{2}}=\frac{m_{p}}{2}

(1)

called the reduced mass and its position vector is $\vec{r}$ w.r.t to CoM co-ordinate, chosen as the origin. Since, both nuclear and Coulomb forces are only dependent on the relative distance co-ordinate, they both are of central character and spherical polar co-ordinates are a good choice for reference system.
Typically the dynamics, at the microscopic level, are governed by time dependent Schrodinger equation and the state of the system is given by it’s wave function: $\Psi(\vec{r,t})$ . But, in this work, we are utilising phase function method where in Schrodinger equation is replaced by Ricatti equation in terms of phase shifts. Hence, scattering phase shifts $\delta(r)$ can be considered to be describing state of the system.

2.1.2 Formulation stage:

Mathematical Model for Interactions:
The Yukawa model [2], which introduced the idea of meson exchange, gave an early theoretical foundation for understanding the strong nuclear force. The short-range nature of the nuclear force is effectively explained by the Yukawa model.
The inclusion of repulsive interaction of the core along with the attractive part is necessary in order to take into account observed scattering phase shifts at higher energies. Such a potential is proposed by Malfliet and Tjon (MT [6], given by:

V_{N}(r)=V_{MT}(r)=-V_{A}\Big{(}\frac{e^{-\mu_{A}r}}{r}\Big{)}+V_{R}\Big{(}\frac{e^{-\mu_{R}r}}{r}\Big{)}

(2)

where, $\mu_{R}$ is taken equal to $2\mu_{A}$ in units of $fm^{-1}$ . Thus, the MT potential is composed of three parameters containing attractive and repulsive parts.
For the charged systems like $pp$ scattering, the Coulomb potential is added to the nuclear interaction potential. But when the particles interact through the Coulomb potential, they do not behave as free particles [17]. Therefore, for the charged particle scattering, the conventional scattering theory needs to be revised. Also due to the presence of Coulomb forces, the asymptotic condition which is the basis of traditional scattering theory, is not valid and hence the conventional definitions of scattering phase shift are invalid [23]. Therfore, within the formalismof non-relativistic scattering theory the standard Coulomb potential is replaced by a short range potential such that at large values of distance $r$ the potential should decrease rapidly [24].
The ambiguity of long range Coulomb potential can be avoided by considering a screened potential which behave as Coulomb potential at short distances and dies down exponentially at large distances [16, 17], say for $r$ $>$ $R$ . The screened Coulomb potential is commonly utilised in numerous fields of physics like nuclear physics, atomic physics, plasma physics and in many dynamic calculations [17].
In this work, we have replaced Coulomb interaction due to protons by Hulthén potential, which is a screened potential [17] and is a modified form of Yukawa potential containing exponential term. The Hulthén potential behaves like a Coulomb potential for small values of distance $r$ , and compared to the Coulomb potential, it rapidly declines exponentially at large values of $r$ ; given by

V_{C}(r)=V_{H}(r)=V_{0}\frac{e^{-r/a}}{1-e^{-r/a}}

(3)

Here $V_{0}$ is the strength of the potential and parameter $a$ is chosen in such a way that $aV_{0}=2k\eta$ . Here $\eta$ is a constant quantity known as Sommerfeld parameter given by

\eta=\frac{\alpha}{\hbar v}

(4)

Where, $\alpha=Z_{1}Z_{2}e^{2}$ and $v=\sqrt{\frac{2E}{\mu}}$ is the relative velocity of reacting particles at large separation.
Hence, $\eta$ will finally be given by

\eta=\frac{Z_{1}Z_{2}e^{2}\mu}{\hbar^{2}k}

(5)

Therefore, $aV_{0}$ after multiplying and dividing by $c^{2}$ will be equal to

aV_{0}=2k\eta=\frac{2Z_{1}Z_{2}e^{2}\mu c^{2}}{\hbar^{2}c^{2}}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ fm^{-1}

(6)

Which is a fixed value for a particular interaction, since $\hbar c=197.327$ MeV-fm, $e^{2}=1.44$ MeV-fm and $\mu c^{2}$ is the reduced mass of the system in the units MeV/ $c^{2}$ .
Hence, the total interaction potential for proton-proton (pp) system is:

V(r)=-V_{A}\Big{(}\frac{e^{-\mu_{A}r}}{r}\Big{)}+V_{R}\Big{(}\frac{e^{-\mu_{R}r}}{r}\Big{)}+\frac{V_{0}}{a}\frac{e^{-r/a}}{1-e^{-r/a}}

(7)

Phase equation for scattering phase shifts:
Phase equation was originally given by Morse [25] and has been briefly discussed here for the sake of completeness. The radial time independent Schr $\ddot{\text{o}}$ dinger wave equation (TISE), for partial wave $\ell=0$ (S-state) is given by

\frac{d^{2}u_{0}(r)}{dr^{2}}+\frac{2\mu}{\hbar^{2}}\big{[}E-V(r)\big{]}u_{0}(r)=0

(8)

Dividing Eq. 8 by $u_{0}(r)$ , we get

\frac{1}{u_{0}(r)}\frac{d^{2}u_{0}(r)}{dr^{2}}+\frac{2\mu}{\hbar^{2}}(E-V(r))=0

(9)

The wavefunction must satisfy condition $u_{0}(0)=0$ at r = 0. Further, at a distance $r_{f}$ beyond which V(r) is zero, one obtains the asymptotic wavefunction $u_{a}(r)$ , as

u_{a}(r)=\alpha(r)\leavevmode\nobreak\ \sin(kr+\delta(r))

(10)

where $k=\sqrt{\frac{2\mu E}{\hbar^{2}}}$ . $\alpha(r)$ and $\delta(r)$ represents amplitude and phase in the wavefunction [3] respectively; depends on distance $r$ and carries information about the region inside the interaction potential $V(r)$ [26].
Also, the wavefunction and its derivative both need to be continuous at a distance $r_{f}$ . That is, choosing $u_{a}(r)$ to be asymptotic solution of Eq. 8 for $r>r_{f}$ , we must have

u_{0}(r)\big{|}_{r=r_{f}}=u_{a}(r)\big{|}_{r=r_{f}}

(11)

similarly

\frac{du_{0}(r)}{dr}\bigg{|}_{r=r_{f}}=\frac{du_{a}(r)}{dr}\bigg{|}_{r=r_{f}}

(12)

These two conditions are combined together into a single equation by considering their logarithmic derivatives to satisfy boundary condition [3, 27], obtained as

\frac{1}{u_{0}(r)}\frac{du_{0}(r)}{dr}\bigg{|}_{r=r_{f}}=\frac{1}{u_{a}(r)}\frac{du_{a}(r)}{dr}\bigg{|}_{r=r_{f}}

(13)

Let us define the logarithmic derivative of wavefunction, based on equation Eq. 13, as A(r), equals to

A(r)=\frac{1}{u_{0}(r)}\frac{du_{0}(r)}{dr}=\frac{1}{u_{a}(r)}\frac{du_{a}(r)}{dr}

(14)

Now differentiating $A(r)$ by using first part of Eq. 14, i.e within the potential region, one obtains

\frac{dA(r)}{dr}=\frac{d}{dr}\bigg{(}\frac{1}{u_{0}(r)}\frac{du_{0}(r)}{dr}\bigg{)}=\frac{1}{u_{0}(r)}\frac{d^{2}u_{0}(r)}{dr^{2}}-\frac{1}{u_{0}^{2}(r)}\bigg{(}\frac{du_{0}(r)}{dr}\bigg{)}^{2}

(15)

From which, we obtain

\frac{1}{u_{0}(r)}\frac{d^{2}u_{0}(r)}{dr^{2}}=\frac{dA}{dr}+A^{2}(r)

(16)

Putting Eq. 16 in Eq. 9, we get

\frac{dA(r)}{dr}+A^{2}(r)=\frac{2\mu}{\hbar^{2}}\big{(}V(r)-E\big{)}

(17)

Eq. 17 is the Riccati type non linear differential equation which can be linearized to get radial time independent Schr $\ddot{\text{o}}$ dinger equation (TISE) back.
Now substituting Eq. 10 in Eq. 14 by considering asymptotic wavefunction $u_{a}(r)$ , we obtain

A(r)=\frac{1}{\alpha\leavevmode\nobreak\ \sin(kr+\delta)}\frac{d}{dr}[\alpha\leavevmode\nobreak\ \sin(kr+\delta)]

(18)

which on simplification gives

A(r)=k\leavevmode\nobreak\ \cot(kr+\delta)

(19)

It is important to keep in mind that while using asymptotic solution, $\delta$ is a constant quantity that needs to be determined at distance $r_{f}$ .
Now representing A(r) within the potential well; where the phase shift changes with interaction potential; in terms of wavefunction $u_{0}(r)$ and taking derivative of Eq. 19, we get

\frac{dA(r)}{dr}=-\frac{(k^{2}+k\frac{d\delta}{dr})}{\sin^{2}(kr+\delta)}

(20)

In this case, $\delta$ is not treated as a constant but rather as a function of $r$ . Further substituting Eq. 20 in Eq. 17, we obtain

-\frac{(k^{2}+k\frac{d\delta}{dr})}{\sin^{2}(kr+\delta)}+k^{2}\frac{\cos^{2}(kr+\delta)}{\sin^{2}(kr+\delta)}=\frac{2\mu}{\hbar^{2}}(V(r)-E)

(21)

Using relation $\cos^{2}\theta=1-\sin^{2}\theta$ , it gets simplified to result in following phase equation

\boxed{\frac{d\delta(r)}{dr}=-\frac{2\mu}{\hbar^{2}}\frac{V(r)}{k}\leavevmode\nobreak\ \sin^{2}(kr+\delta(r))}

(22)

Equation Eq. 22 is known as phase equation and the function $\delta(r)$ , is known as “phase function”. Phase function is a function of interaction potential $V(r)$ at a particular interaction distance $r$ and also at a particular value of laboratory energy in terms of wave number $k$ . The phase function $\delta(r)$ at any given location $r=R$ yields the value of “phase shift” corresponding to interaction potential $V(r)$ . A similar expression can be found in work by Morse and Viterbo [25, 26]. Being a nonlinear equation, Eq. 22 has no analytical solutions and hence one has to take the approach of numerical methods.

2.2 Preparation of the system for numerical solution

Rephrasing the problem in appropriate units
Choice of units: In nuclear physics, the correct choice of units for energy and distance are MeV and fm respectively. So, $\hbar c$ gets converted from J-m to $197.329$ MeV-fm. Reduced mass is chosen in MeV/c². Hence, phase equation on the right is multiplied and divided by c² to obtain

\boxed{\frac{d\delta(r)}{dr}=-\frac{2\mu c^{2}}{\hbar^{2}c^{2}}\frac{V(r)}{k_{a}}\leavevmode\nobreak\ \sin^{2}(k_{a}r+\delta(r))}

(23)

with $k_{a}=\sqrt{\frac{2\mu c^{2}E}{\hbar^{2}c^{2}}}$ .
The physical system which is to be solved analytically or numerically should be expressed in appropriate units. If all the parameters are not in appropriate units, there will be round off errors in the solution obtained. Therefore the physical problem should be rephrased in suitable set of units in order to obtain acceptable solutions.
Observing phase equation Eq. 23, it is clear that $k$ must have inverse units of $r$ , so that argument inside sine function remains dimensionless. We know $\hbar$ has units of J-s and there is no time appearing in the equation. So, idea is to multiply $\hbar$ by speed of light $c$ with units of $m/s$ so that we have energy-meter units. The suitable choice of units in nuclear physics are MeV for energy and $fermi$ $(fm)$ for distance. Hence by considering value of $\hbar c$ equal to $197.329$ MeV-fm, J-m will be converted to MeV-fm. The trick is to multiply and divide factor $\frac{\hbar^{2}}{2\mu}$ by $c^{2}$ and choose value of reduced mass in MeV/ $c^{2}$ . So, factor will have units of MeV- $fm^{2}$ . This results in $k$ having

k=\sqrt{\frac{E}{\hbar^{2}/2\mu}};\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \bigg{(}\text{units\leavevmode\nobreak\ inside\leavevmode\nobreak\ square\leavevmode\nobreak\ root}:\leavevmode\nobreak\ \frac{(MeV)}{(MeV-fm^{2})}\bigg{)}

(24)

units of $fm^{-1}$ as needed. Finally, observe that LHS has units of $fm^{-1}$ and factor multiplying $\sin^{2}$ term also has same units by choosing potential to have units of MeV. This requires, attractive and repulsive strengths, $V_{A}$ and $V_{R}$ in Malfliet-Tjon potential to be in units of MeV-fm. Of course, $\mu_{A}$ has to be in units of $fm^{-1}$ .
Region of Interest (RoI): One has to understand the nature of interaction to chose appropriate of region of interest in which the phase equation can be solved. Since, nuclear force is of short range, it is important to study the long range Coulomb for various screening radii to decide on the final integration distance, say $r_{f}$ , where Coulomb potential approaches to zero (i.e $V_{C}$ $\sim$ $10^{-3}$ MeV). We have plotted the screened Coulomb potentials for various values of screening parameter $a$ in Fig.1.

Refer to caption — Figure 1: Variation of screened Coulomb potential w.r.t distance ‘ $r$ ’ for different values of screening parameter ‘a’ (fm). The actual Coulomb potential is also plotted for comparison.

One can observe that with increasing screening parameter value $a$ , the potential falls off slowly with distance but has smaller repulsive strength at $r$ close to 0. The cutoff values for $a=2,5,10$ and $20$ are $r_{f}=15,30,50$ and $90$ respectively. So, RoI would have to be accordingly chosen and sampled with a certain step size, say $h$ .
It is clear from the Fig.1 that the screened Coulomb potential plotted with value of screening parameter $a=5$ following the same trend as that of the variation of original Coulomb potential for proton-proton interaction. Also, the distance up-to which effect of potential is felt for a particular interaction is limited to a certain region. Hence final distance ( $r_{f}$ ) considered for calculations should be little greater than the total interaction radius of system under consideration. In current work for proton-proton system, we have taken $r_{f}=30$ fm since at this distance, the screened Coulomb potential reduces to zero approximately.

2.

Choice of Numerical Technique:
The first order phase equation can be solved by using Runge Kutta methods. While both RK-2 and RK-4 methods are stable, they have accuracies given by $O(h^{2})$ and $O(h^{4})$ respectively. Considering that $pp$ expected scattering phase shift data is available up to three decimal places, it is preferable to utilise RK-4 technique.
In addition, when compared to the RK-4 method, the RK-2 method requires significantly more iterations to obtain the required accuracy. This could be due to the fact that the RK-2 method requires only two slopes to accomplish the solution, whereas the RK-4 method requires more slopes, resulting in early convergence of the results [28].
Considering an initial value of function $\delta(r_{0})$ = $\delta_{0}$ , Runge kutta methods involve calculating value of function $\delta(r_{1})$ , by using the previous values of the function at $r_{i}$ (i=0,1,2,…n-1). Runge Kutta fourth order method (RK-4) [29] consists of four steps with solution given as:

$\delta(r_{i+1})=\delta(r_{i})+\frac{h}{6}(F_{1}+2F_{2}+2F_{3}+F_{4})$ (25)

where $h$ is the step size and $F_{i}^{\prime}s$ represent the slopes of the function $f$ at various points within the interval $[0,h]$ . Since, the distance $r_{f}=r_{n-1}$ is a distance at which potential becomes negligible, the corresponding $\delta(r_{f})$ will reach a constant value.

2.3 Implementation of the numerical method in a computer

The algorithm for obtaining scattering phase shifts using RK-4 method is given in supplemental material provided separately. Once algorithm is carefully tested in a worksheet environment such as Gnumeric/Excel, one can translate it into a code or program in any software such as Scilab, Matlab or Python.
Mean square error (MSE) for obtained SPS, is calculated as

MSE=\frac{1}{N}\sum_{i=1}^{N}(|\delta^{e}_{i}-\delta^{o}_{i}|)^{2}

(26)

Where $\delta^{e}_{i}$ and $\delta^{o}_{i}$ are the expected phase shifts and obtained phase-shifts for different values of energies $E_{i}$ .
An optimisation routine based on Variational Monte-Carlo (VMC) tehnique[30, 31, 32] is utilised, where in the model parameters are randomly varied in each iteration, so that MSE value converges on to a minimum value.
Of course one could use other cost functions such as mean absolute percentage error (MAPE) or chi-squared error ( $\chi^{2}$ ). The optimisation code is explained briefly in supplemental material provided separately.

2.4 Experimental observable quantities

Cross section
After obtaining the SPS $\delta_{\ell}(E)$ , one can calculate the partial cross section $\sigma_{\ell}(E)$ , which for S-wave ( $\ell=0$ ) is given by:

\sigma_{0}(E)=\frac{4\pi}{k_{a}^{2}}\sin^{2}({\delta_{0}(E)})

(27)

Scattering Parameters
Scattering parameters i.e scattering length ‘ $a_{0}$ ’ and effective range ‘ $r_{0}$ ’, obtained from scattering phase shifts, can be calculated for low energy $pp$ scattering by using effective range theory [1] given by:

k\cot(\delta_{0})=-\frac{1}{a_{0}}+\frac{1}{2}r_{0}k^{2}

(28)

The slope and intercept from straight line plot of $\frac{1}{2}k^{2}$ vs $k\cot(\delta_{0})$ gives scattering parameters ‘ $a_{0}$ ’ and ‘ $r_{0}$ ’.

3 Results and Discussion

In order to calculate the scattering phase shifts, we have integrated phase equation Eq. 23 from origin to asymptotic region by using RK-4 method. Mean square error (MSE) is then optimised by using Variational Monte Carlo technique w.r.t the expected phase shift values of [33]. The advantage of considering PFM is its simplicity that only interaction potential is required for calculating scattering phase shifts.
Model parameters for net interaction potential were determined by minimizing the mean square error (MSE) for different values of screening parameter $a$ . The value of $a$ for which mean square error is minimum is finalised. After finalising the value of $a$ , the corresponding final integration distance $r_{f}$ is taken as the distance greater than the distance $r_{C}$ where Coulomb potential ( $V_{C}$ ) $\approx$ $10^{-3}$ MeV i.e $r_{f}>r_{C}$ .
For $pp$ interaction, we have obtained screened Coulomb potential for different values of screening parameter $`a^{\prime}$ , e.g $a=2,5,10,20$ fm etc. And for each value of $a$ , the depth of potential and corresponding mean square error for simulated scattering phase shifts is calculated. Out of all the value of $a$ , the mean square error comes out to be least for $a=5$ fm i.e mean square error= $0.29$ . The depth of the interaction potential which matches with that of realistic potential for ${}^{1}S_{0}$ state of $pp$ interaction ( $\approx 100$ MeV) also corresponds to value of $a=5$ fm. For $a=5$ fm, the Coulomb potential approaches to zero near $r_{C}=29$ fm. Hence, we conclude that ‘a’ $=5$ fm is the appropriate value for calculating screened Coulomb potential for $pp$ scattering with minimum mean square error = $0.29$ and with final interaction distance $r_{f}=30$ fm.
The optimised model parameters obtained by VMC technique, for MT potential with $a=5$ and $r_{f}=30$ fm are given in Table 1. We have plotted the final interaction potential for ${}^{1}S_{0}$ state of $pp$ scattering by using the model parameters for MT potential given in Table 1 and is shown in Fig. 2(a).

Table 1: Model parameters for MT potential and mean square error (MSE) w.r.t expected data for singlet (

{}^{1}S_{0}

) state of pp interaction.

$V_{r}$	$V_{a}$	$\mu_{A}$	$a$	MSE
(MeV-fm)	(MeV-fm)	( $fm^{-1}$ )	(fm)
7171.2823	1570.9068	2.45	5	0.29

The MT potential obtained by using our optimised parameters is in good agreement with that of realistic Argonne AV- $18$ potential [5, 35] given by S. Ohkubo [34] for ${}^{1}S_{0}$ state of proton-proton interaction thereby validating our approach. It is found that the interaction potential so obtained with parameters has depth $V_{d}=-99.44$ MeV with equilibrium distance $r_{d}=0.85$ fm, closely match with that of realistic Argonne AV- $18$ potential [5, 34, 35] having depth $\approx$ -103.5 MeV at equilibrium distance $r_{d}=0.89$ fm. Hence we can validate that the MT potential parameters which we have obtained for $pp$ scattering, describes the interaction well.
The Coulomb barrier height $V_{CB}$ is also calculated for interaction potential obtained. $V_{CB}$ is found to be $=0.11$ MeV at distance $r_{C}=4.3$ fm. The variation beyond matches with screened Coulomb potential for $a=5$ fm. This is because nuclear force being short range, becomes $\sim$ 0 for distance of about $4-5$ fm itself.
Simulated phase shifts obtained by considering MT potential parameters from Table 1 along with expected phase shifts for ${}^{1}S_{0}$ state of pp scattering are given in Table 2. It is evident that that all the phase shifts have discrepency less than $6$ ( $\%$ ). The variation of scattering phase shifts for lab energies $1-350$ MeV are shown in Fig. 2(b) along with expected data points from Granada [33]. From the figure, one can observe an excellent match at all energies.

Table 2: Scattering phase shifts (in degrees) for

{}^{1}S_{0}

states of

pp

interaction, w.r.t expected data [33], for lab energies up-to

350

MeV. The absolute percentage error is also given for each data point in fourth column and the mean absolute percentage error(MAPE) is given at the bottom.

Lab energy	${}^{1}S_{0}$ (Expected)[33]	${}^{1}S_{0}$ (Malfliet-Tjon)	Absolute percentage
$E_{lab_{i}}$ (MeV)	$\delta_{e_{i}}$	$\delta_{o_{i}}$	error ( $\%$ )
1	32.677	34.285	4.92
5	54.895	53.348	2.82
10	55.32	54.383	1.69
25	48.848	48.662	0.38
50	39.182	39.137	0.11
100	25.357	25.199	0.62
150	15.229	15.074	1.02
200	7.076	7.057	0.27
300	-5.694	-5.353	5.99
350	-10.828	-10.397	3.98
MAPE ( $\%$ )			2.18

One of the main advantage of phase function method is that, it provides the information of how phase shifts progressively builds-up due to interaction potential which can be observed by plotting phase shift w.r.t interaction distance $r$ . This feature helps students to understand the effect of interaction potential on phase shift as the interaction distance increases. We have plotted phase shifts w.r.t distance for different laboratory energies. The phase shifts $\delta(r)$ w.r.t distance $r$ for lab energies $1,5,10,25,50,100,200,300$ and $350$ MeV are plotted in Fig. 3 by considering the parameters from Table 1. It is clear from the figure that the phase shifts remains negative up-to approximately $1$ fm, till the repulsive core of the potential and then starts increasing towards positive side when the potential becomes attarctive. Phase shift reaches to a constant value at $r=r_{f}=30$ fm where the potential becomes negligible. The phase shift values vary for each value of incident energy as well.

Partial cross section ( $\sigma_{0}$ ) for ${}^{1}S_{0}$ state of pp scattering is calculated by using Eq. 27 and is given in figure 5. Since we have calculated cross section corresponding to $\ell=0$ partial wave (i.e S-state) only, therefore there is a slight deviation of the calculated $\sigma_{0}$ from experimental data taken from tables (Table 4, 5) published by R. A. Arndt $et.al.$ in their work [36].

Table 3: Scattering length ‘

a

’ and effective range ‘

r_{0}

’ parameters calculated for

{}^{1}S_{0}

state of pp scattering along with the expected values, the errors are given within brackets

States	$a(fm)$ [37]	$a(fm)$ (calc.)	$r_{0}(fm)$ [37]	$r_{0}(fm)$ (calc.)
${}^{1}S_{0}$ - $pp$	-7.8063 ( $\pm$ 0.003)	-7.5857 ( $\pm$ 0.0008)	2.794 ( $\pm$ 0.001)	1.8272 ( $\pm$ 0.021)

The Scattering parameters i.e scattering length ‘ $a_{0}$ ’ and effective range ‘ $r_{0}$ ’ parameters for ${}^{1}S_{0}$ state of pp scattering are also calculated. These scattering parameters have been calculated from slope and intercept of $kcot(\delta_{0})$ $vs$ $(0.5\times k^{2})$ plots given in Fig. 5(a). The obtained scattering parameters along with expected data are tabulated in Table 3. It can be seen that the calculated scattering parameters are in good agreement with the results obtained by Miller $et.al.$
Based on the above discussions, we conclude that:

1.

The scattering phase shifts calculated using simple phase function method along with our VMC optimisation procedure shows very good results with mean square error of 0.29 w.r.t expected data of Granda group [33].
2.

Following the pedagogical approach, the algorithm as a pseudo code for our VMC optimisation procedure and the implementation in Gnumeric worksheet to obtain scattering phase shift has been given in the supplemental material provided separately for easy access of students.
3.

The interaction potential which we have obtained by using our optimisation procedure for ${}^{1}S_{0}$ state of pp scattering matches well with the realistic Argonne Av- $18$ potential given by Wiringa $et.al$ [35, 34], thereby validating our approach.
4.

The scattering phase shifts calculated using simple phase function method by using only four parameters of total interaction potential, three for Malfliet-Tjon potential and one for the screened Coulomb potential, unlike most of the relistic potentials which require approximately 40 parameters for the calculation.
5.

The scattering parameters, i.e scattering length ‘ $a$ ’ and effective range ‘ $r_{0}$ ’ calculated by using obtained scattering phase shifts also matches with the expected values given by Miller $et.al$ [37].
6.

Scattering cross section so calculated shows matching trend with the experimental cross section data for $pp$ scattering taken from [36]. Discrepencies in the calculated cross section are due to the fact that in this current work, we have not considered higher partial waves for $pp$ scattering which is our future work.

4 Conclusion

The proton-proton ( $pp$ ) interaction modeled, using phenomenological Malfliet-Tjon potential for nuclear interaction and screened Coulomb potential has been found to explain the observed scattering phase shifts for the ${}^{1}S_{0}$ singlet state. The obtained scattering phase shifts match the expected data with a mean absolute percentage error of 2.18 $\%$ and mean square error of $0.29$ . The scattering cross-sections and scattering parameters are found to be reasonably close to the respective expected data. This procedure based on phase function method for understanding $np$ and $pp$ scattering can be easily integrated into classroom and laboratories at the undergraduate level. Further, students can be encouraged to take-up projects by considering other interactions like nucleon-nucleus ( $nd$ , $pd$ , $n\alpha$ , $p\alpha$ ) and nucleus-nucleus ( $\alpha-\alpha$ ) scattering systems.

Declarations:

Funding: We hereby state that there is no funding received for this work.

Conflict of Interest: There is no conflict of interest whatsoever. References

References

[1] Krane, Kenneth S. “Introductory nuclear physics”. John Wiley $\&$ Sons, (1991).
[2] Yukawa, H., Sakata, S., Kobayasi, M. and Taketani, M., 1955. On the Interaction of Elementary Particles. IV. Progress of Theoretical Physics Supplement, 1, pp.46-71.
[3] Bertulani, C.A. and Danielewicz, P., 2019. Introduction to nuclear reactions. CRC Press.
[4] Naghdi, M. “Comparing some nucleon-nucleon potentials.” Physics of Particles and Nuclei Letters 11 (2014): 410-431.
[5] Wiringa, R. B., V. G. J. Stoks, and R. Schiavilla. “Accurate nucleon-nucleon potential with charge-independence breaking.” Physical Review C 51, no. 1 (1995): 38.
[6] Malfliet, R. A., and J. A. Tjon. “Solution of the Faddeev equations for the triton problem using local two-particle interactions.” Nuclear Physics A 127, no. 1 (1969): 161-168.
[7] Awasthi, S. , Khachi, A., Kumar, L., and Sastri, O. S. K. S. “Numerical Simulation Study of Neutron-Proton Scattering using Phase Function Method.” arXiv preprint arXiv:2310.07336 (2023).
[8] Wong, Samuel SM. “Introductory nuclear physics”. John Wiley $\&$ Sons, 2008.
[9] Arya, A.P., “Fundamentals of nuclear physics”. Allyn and Bacon, Inc. Boston, (1966).
[10] Bohr, A.N. and Mottelson, B.R., 1998. Nuclear Structure (in 2 volumes). World Scientific Publishing Company.
[11] Wigner, Eugene Paul, and Leonard Eisenbud. “Higher angular momenta and long range interaction in resonance reactions.” Physical Review 72, no. 1 (1947): 29.
[12] Mackintosh, Raymond S. “Inverse scattering: applications to nuclear physics.” arXiv preprint arXiv:1205.0468 (2012).
[13] Jost, Res, and Abraham Pais. “On the scattering of a particle by a static potential.” Physical review 82, no. 6 (1951): 840.
[14] Babikov, V. V. “The phase-function method in quantum mechanics.” Soviet Physics Uspekhi 10, no. 3 (1967): 271.
[15] Calogero, Francesco. Variable Phase Approach to Potential Scattering by F Calogero. Elsevier, 1967.
[16] M.A. Preston, Structure of the Nucleus. CRC Press, (2018).
[17] Laha, U., 2020. Nuclear Hulthén potential and the scattering phase shifts for l= 3. Current Science (00113891), 118(4).
[18] Bhoi, J., R. Upadhyay, and U. Laha. “Parameterization of Nuclear Hulthén Potential for Nucleus-Nucleus Elastic Scattering.” Communications in Theoretical Physics 69, no. 2 (2018): 203.
[19] Babenko, V. A., and Petrov, N. M. ”Study of the Pion-Nucleon Coupling Constant Charge Dependence on the Basis of the Low-Energy Data on Nucleon-Nucleon Interaction.” arXiv preprint arXiv:1604.02912 (2016).
[20] Aditi Sharma and O S K S Sastri 2020 Eur. J. Phys. 41 055402 (2020).
[21] A. Sharma, S. Gora, J Bhagavathi, and OSK S. Sastri. “Simulation study of nuclear shell model using sine basis.” Am. J. Phys. 88, 7, 576-585 (2020).
[22] O. S. K. S. Sastri, A. Sharma, J. Bhardwaj, S. Gora, V. Sharda, and B. P. Jithin. “Numerical solution of square well potential with matrix method using worksheets.” Phys. Educ. 36 (2019).
[23] Taylor, John R. “A new rigorous approach to Coulomb scattering.” Il Nuovo Cimento B Series 11 23, no. 2 (1974): 313-334.
[24] Semon, M. D., and J. R. Taylor. “Scattering by Potentials with Coulomh Tails.” Il Nuovo Cimento A (1971-1996) 26 (1975): 48-58.
[25] P.M. Morse, and W.P. Allis. “The effect of exchange on the scattering of slow electrons from atoms.” Physical Review, 44 (4), 269 (1933).
[26] V. D. Viterbo, Nelson HT Lemes, and João P. Braga. “Variable phase equation in quantum scattering.” Revista Brasileira de Ensino de Fisica 36 (2014): 1310.
[27] A. N. Khondker, M. Rezwan Khan, and A. F. M. Anwar. “Transmission line analogy of resonance tunneling phenomena: The generalized impedance concept.” Journal of applied physics 63 (10), 5191-5193 (1988).
[28] V. Zhaba.“Calculation of phases of np- scattering up to $T_{lab}=3$ GeV for Reid68 and Reid93 potentials on the phase-function method.”arXiv preprint arXiv:1604.06006 (2016).
[29] J. Butcher, ”Runge-kutta methods.” Scholarpedia 2, no. 9 (2007): 3147.
[30] A. Sharma and O.S.K.S. Sastri. ”Numerical solution of Schrödinger equation for rotating Morse potential using matrix methods with Fourier sine basis and optimization using variational Monte‐Carlo approach.” International Journal of Quantum Chemistry 121, no. 16 (2021): e26682.
[31] Sastri, O. S. K. S., Anil Khachi, and Lalit Kumar. “An Innovative Approach to Construct Inverse Potentials Using Variational Monte-Carlo and Phase Function Method: Application to np and pp Scattering.” Brazilian Journal of Physics 52, no. 2 (2022): 1-6.
[32] S. Gora, O. S. K. S. Sastri, and S. K. Soni. “Optimization of semi-empirical mass formula co-efficients using least square minimization and variational Monte-Carlo approaches.” European Journal of Physics 43 (3), 035802 (2022).
[33] Pérez, R. Navarro, J. E. Amaro, and E. Ruiz Arriola. “The low-energy structure of the nucleon–nucleon interaction: statistical versus systematic uncertainties.” Journal of Physics G: Nuclear and Particle Physics 43, no. 11 (2016): 114001.
[34] S. Ohkubo, “Luneburg-lens-like structural Pauli attractive core of the nuclear force at short distances.” Physical Review C 95, no. 4 (2017): 044002.
[35] V. I. Zhaba, ”The variable phase approach: phase, amplitude and wave functions of the states for np-system for Argonne v18 potential.” World Scientific News 123 (2019): 161-180.
[36] R. A. Arndt, W. J. Briscoe, A. B. Laptev, I. I. Strakovsky, and R. L. Workman. “Absolute Total np and pp Cross-Section Determinations.” Nuclear science and engineering 162, no. 3 (2009): 312-318.
[37] Miller, Gerald A., Bernard Marie Karel Nefkens, and Ivo Šlaus. ”Charge symmetry, quarks and mesons.” Physics Reports 194, no. 1-2 (1990): 1-116.