Nilsson Orbital using diagonalization method

Leave a comment

Long time ago, I tried to tackle the Nilsson orbital by solving the Hamiltonian analytically. However, the Hamiltonian is without LS coupling. This times, I redo the calculation according to the reference B. E. Chi, Nuclear Phyiscs 83 (1966) 97-144.


The Hamiltonian is

\displaystyle H = \frac{P^2}{2m} + \frac{1}{2}m\left( \omega_\rho^2 (x^2+y^2) + \omega_z^2 z^2 \right) + C L\cdot S + D L\cdot L

using

\displaystyle \omega_\rho = \omega_0 \left(1+\frac{2}{3}\delta\right)^{\frac{1}{2}} = \omega \left(\frac{3+2\delta}{3-4\delta}\right)^{1/6}

\displaystyle \omega_z = \omega_0 \left(1-\frac{4}{3}\delta\right)^{\frac{1}{2}} = \omega \left(\frac{3-4\delta}{3+2\delta}\right)^{1/3}

\displaystyle \beta = \frac{4}{3}\sqrt{\frac{\pi}{5}}\delta

\displaystyle r^2 Y_{20}(\theta, \phi) = \frac{1}{4}\sqrt{\frac{5}{\pi}} (3z^2-r^2)

The Hamiltonian becomes

\displaystyle H = -\frac{\hbar^2}{2m}\nabla^2 +\frac{1}{2} m \omega_0^2 r^2 - \frac{1}{2} m\omega_0^2 r^2 \beta Y_{20} + C L\cdot S + D L\cdot L

Set x_i^2 \rightarrow  x_i^2 \frac{\hbar}{m \omega_0} , and r^2 \rightarrow \rho^2 \frac{\hbar}{ m \omega_0}

\displaystyle \frac{H}{\hbar\omega_0} = \frac{1}{2}(-\nabla^2 + \rho^2)  - \rho^2 \beta Y_{20} - 2 \kappa L\cdot S - \mu \kappa L\cdot L

Set

\displaystyle \frac{H_0}{\hbar\omega_0} = \frac{1}{2}(-\nabla^2 + \rho^2) - 2 \kappa L\cdot S - \mu \kappa L\cdot L

and the perturbation is

\displaystyle \frac{H_p}{\hbar\omega_0} =  - \rho^2 Y_{20}


The wavefunction for the spherical harmonic is

\displaystyle |Nljk\rangle = A r^l e^{-\frac{r^2}{2}} L_{\frac{N-l}{2}}^{l + \frac{1}{2}}(r^2) \sum_{m m_s} Y_{lm}(\theta, \phi) \chi_{\frac{1}{2} m_s} C_{lm\frac{1}{2} m_s}^{jk}

\displaystyle A = \sqrt{\frac{(\frac{N-l}{2})!(\frac{N+l}{2})! 2^{N+l+2}}{\sqrt{\pi} (N+l+1)!}}


The diagonal elements are

\displaystyle \frac{1}{\hbar \omega_0 }\langle Nljk|H_0|Nljk\rangle = N + \frac{3}{2} - \kappa \langle L\cdot S \rangle - \mu \kappa l(l+1)

where \langle L \cdot S \rangle = \frac{1}{2} ( j(j+1) - l(l+1) - \frac{3}{4} )

The off diagonal elements are

\displaystyle \frac{1}{\hbar \omega_0 }\langle Nljk|H_p|Nljk\rangle = - \langle Nljk| r^2 Y_{20}|Nljk\rangle

( I will evaluate this integral in future )


The rest is diagonalization the Hamiltonian

\displaystyle H = H_0 + \beta H_p

Here is the calculation for the 2nd harmonic for \kappa = 0.05, \mu = 0

Screen Shot 2019-07-25 at 18.25.45.png

The component of each orbital can be directly taken from the eigenvalue. Here is the [521]1/2 state. \kappa = 0.05, \mu(N=3) = 0.35, \mu(N=4) = 0.625, \mu(N=5) = 0.63

Screen Shot 2019-07-25 at 18.28.27.png

 

Rotational Band

Leave a comment

For deformed nuclei, it can be rotated in various angular momentum in Laboratory frame. Assume rigid body rotation, the energy is

\displaystyle E_J = \frac{1}{2}I\omega^2 = \frac{1}{2I}J^2

In QM, that becomes

\displaystyle H = \sum_{i=1}^{3} \frac{\hbar^2}{2I_i} J_i^2

For axial symmetry, I_1 = I_2 = I

\displaystyle H = \frac{\hbar^2}{2I} (J^2 - J_3^2) + \frac{\hbar^2}{2I_3}J_3^2

Remember, in deformed nuclei, the projection of J along the symmetry axis in the body-frame is K . The expected value of the Hamiltonian with state |Nn_z m_l K \rangle in the body-frame is proportional to J(J+1) for J^2 and K for J_3. i.e.

\displaystyle E_J = \frac{\hbar^2}{2I} J(J+1) + E_K

From body-frame to Lab-frame, we should apply the Wigner D-Matrix to the intrinsic wave function. ( I am not sure the following equation is correct, but the idea is rotating the body-frame wavefunction with Wigner D-Matrix to get the Lab-frame wave function. In Lab frame the total angular momentum must be a good Quantum number as rotational symmetry restored, so as J_z = M. The problem of the following equation is that the J is not a good Q-number in Nilsson wavefunction )

\displaystyle |JMK\rangle = \sum_{M} D_{MK}^{J} |Nn_zm_lK\rangle

However, the Wigner D-Matrix does not conserve parity transform:

\displaystyle D_{MK}^J \rightarrow (-1)^{J+K} D_{M-K}^{J}

In order to restored the parity, we need to include \pm K in the Lab-frame wave function.

\displaystyle |JMK\rangle = \sum_{M} \left( D_{MK}^J \pm (-1)^{J+K} D_{M-K}^J \right) |Nn_zm_lK\rangle

where + for positive parity, – for negative parity.

From the above equation, for K^\pi = 0^+ (0^-), J must be even (odd). For K > 0 , J = K, K+1, K+2, ... .

rotaional Band of 205Fm.png

rotational band of 253No.png

rotational band of 19F.png

We can see for K = 1/2 , the J = 5/2, 9/2, 11/2 are lower to the main sequence. This was explained by adding an extra term in the rotation Hamiltonian that connect \Delta K = 1 .

\displaystyle \langle JMK | H'(\Delta K = 1) |JMK \rangle

\displaystyle \rightarrow \langle D_{MK}^J | H' | D_{MK}^J \rangle+ \langle D_{M-K}^J |H'| D_{MK}^J \rangle + \langle D_{MK}^J | H' | D_{M-K}^J \rangle+ \langle D_{M-K}^J |H'| D_{M-K}^J \rangle

The term with \Delta K = 0 vanished. And since $\latex K = 1$, the only non-zero case is K = 1/2 .

A possible form of the H' (\Delta K = 1) = \frac{1}{2} \omega (J_+ + J_-) . These are the ladder operator to rise or lower the m-component by 1. In 19F case, we can think it is a single proton on top of 18O core.  A rotation core affect the proton with an additional force, similar to Coriolis force on earth.

 

Axial Harmonic Oscillator – Nilsson Orbit (III)

Leave a comment

One of the problem (or difficulty) is the conversion between cylindrical quantum number |Nn_zm\rangle and spherical quantum number |Nlm\rangle. In the limit of \delta = 0, all state with same N has same energy and the conversion is not necessary.

With out LS coupling, the energy level can be found by the approximation

\displaystyle E(N, n_z) = \hbar\omega \left( \frac{3}{2} + \left(1 + \frac{\delta}{3} \right) N - \delta n_z \right)

Nilsson_noLS.PNG

We can see from the formula, for \delta>0 , the largest n_z has lowest energy. And the conversion to spherical harmonics can be done using projection method.

With LS coupling, the conversion becomes complicated. One issue is the coupling with the spin. Since the total angular momentum J = L + S is not a conservative quantity due to deformation. In the body-fixed frame, the additional good quantum number is the z-projection of J , which denote as K = |m_l \pm 1/2 |, the quantum state becomes

|Nn_z m_l K\rangle

Note that each state will have only 2 degeneracy with negative K . On the spherical basis, the coupling with J is a standard textbook content. the eigen state is

|N l j mj \rangle = |N l j K\rangle

The conversion get complicated, because there is no clear rule on how | Nn_z m_l K\rangle \rightarrow |Nl j m_j \rangle . When we want to calculate the energy level using cylindrical basis with the LS coupling, since the J is not clear, even through we can write down the formula

\displaystyle E( N, n_z, l, j) = \hbar \omega_z \left(\frac{1}{2} + n_z \right) + \hbar \omega_\rho \left(1 + N - n_z \right) \\ + a \frac{1}{2} ( j(j+1) - l(l+1) - s(s+1)) + b l^2

It is difficult to know the j , l .


The go around this conversion, we can diagonal the Hamiltonian using spherical basis. The Hamiltonian is can be written as

\displaystyle H = - \frac{\hbar^2}{2m} \nabla^2 + \frac{m\omega}{2}r^2 \left( 1 - \frac{2}{3}\sqrt{\frac{16\pi}{5}} Y_{20}(\theta, \phi) \delta \right) + a L\cdot S + b L^2

The first 2 terms is the spherical harmonic oscillator. The energy is \hbar\omega( 3/2 + N) . The Y_20 restricts the coupling between l, l \pm 2 (see here.)

For example, when using the basis on N = 1, which are

\displaystyle |N l j m_j \rangle = \left( |1 1 \frac{3}{2} \frac{3}{2}\rangle , |1 1 \frac{1}{2} \frac{1}{2}\rangle, |1 1 \frac{3}{2} \frac{1}{2}\rangle \right)

at \delta = 0.5, the matrix is

\begin{pmatrix} 2.633 & 0 & 0 \\ 0 & 2.233 & 0 \\ 0 & 0 & 2.233\end{pmatrix}

The eigen values and eigen state are

\displaystyle  2.633 , |1 1 \frac{3}{2} \frac{3}{2}\rangle

\displaystyle  2.233 , |1 1 \frac{1}{2} \frac{1}{2}\rangle ,  |1 1 \frac{3}{2} \frac{1}{2}\rangle

The lower level has mixed j !

Of course, a more complete diagonalization should be admixture from other major shell as well. When we disable the LS coupling, using 7 major shells, the Nilsson diagram look like this

Nilsson_noLS2.PNG

This is consistence with omega_0 = const. , i.e. no volume conservation. With LS coupling switched on.

Nilsson_LS.PNG

Although it is not so clear, we can see the state avoid “crossing”.  The s-state are unaffected by the LS coupling.

One deflect in the plot is that, the slope seem to be incorrect. And also, I tried to add volume conservation, but the result does not consistence….. (sad).

 

 

 

Axial Harmonic Oscillator – Nilsson Orbit (II)

Leave a comment

This time we will project the deformed orbit into spherical orbit. The wavefunctions are stated in here again.

\displaystyle \Psi^D_{N n_z m}(z, \rho, \phi) = |Nn_z m\rangle_D \\ = \sqrt{\frac{1}{\alpha_z \alpha^2}}\sqrt{\frac{ n_\rho !}{2^n_z  n_z! (m + n_\rho)!\sqrt{\pi^3}}} H_{n_z}\left(\frac{z}{\alpha_z}\right) \\ \exp\left(- \frac{1}{2}\left(\frac{z^2}{\alpha_z^2}+\frac{\rho^2}{\alpha^2}\right)\right) \left(\frac{\rho}{\alpha}\right)^{m} L_{n_\rho}^{m} \left(\frac{\rho^2}{\alpha^2}\right) \exp(i m \phi) 

\displaystyle \Psi^S_{Nlm}(r, \theta, \phi) = |N l m_l\rangle_S\\ =\sqrt{ \frac{1}{\sqrt{\pi}\alpha^{2l+3}} \frac{(\frac{N-l}{2})! (\frac{N+l}{2})! 2^{N+l+2}}{(N+l+1)!}} r^l \exp\left(-\frac{r^2}{2\alpha^2}\right) L_{k}^{l+\frac{1}{2}}\left( \frac{r^2}{\alpha^2} \right) Y_{lm}(\theta, \phi)

Since both functions span the entire space and are basis, thus, we can related them as

\displaystyle|Nn_z m \rangle_D = \sum_{N' l' m'} C_{N n_z m}^{N' l' m'} |N' l' m' \rangle_S

where

C_{N n_z m}^{N' l' m'} = \langle (N'l'm')_S | (N n_z m)_D \rangle

First thing we notice is that m' = m . Because the phi components are the same in both wave function. i.e. \int \exp(- i m' \phi) \exp(i m \phi) d\phi = 0 if m' \neq m. We can omit the m, m' , so that C_{N n_z}^{N' l'}

Second thing is the parity must be the same, thus when N is even (or odd), N' must be even (or odd).


For non-deformed \delta = 0 , N' = N ,, here are some results

\displaystyle|0 0 0 \rangle_D = |0 0  0\rangle_S

\displaystyle|1 1 0 \rangle_D = |1 1 0 \rangle_S

\displaystyle|1 0 1 \rangle_D = |1 1 1 \rangle_S

\displaystyle|2 2 0 \rangle_D = \sqrt{\frac{2}{3}} |2 2 0 \rangle_S - \sqrt{\frac{1}{3}} |2 0 0 \rangle_S

\displaystyle|2 1 1 \rangle_D = |2 2 1 \rangle_S

\displaystyle|2 0 2 \rangle_D = |2 2 2 \rangle_S

\displaystyle|2 0 0 \rangle_D = \sqrt{\frac{1}{3}} |2 2 0 \rangle_S + \sqrt{\frac{2}{3}} |2 0 0 \rangle_S

\displaystyle|3 3 0 \rangle_D = \sqrt{\frac{2}{5}} |3 3 0 \rangle_S - \sqrt{\frac{3}{5}} |3 1 0 \rangle_S

\displaystyle|3 2 1 \rangle_D = -\sqrt{\frac{4}{5}} |3 3 1 \rangle_S + \sqrt{\frac{1}{5}} |3 1 1 \rangle_S

\displaystyle|3 1 2 \rangle_D = |3 3 2 \rangle_S

\displaystyle|3 0 3 \rangle_D = |3 3 3 \rangle_S

\displaystyle|3 1 0 \rangle_D = \sqrt{\frac{3}{5}} |3 3 0 \rangle_S + \sqrt{\frac{2}{5}} |3 1 0 \rangle_S

\displaystyle|3 0 1 \rangle_D = -\sqrt{\frac{2}{5}} |3 3 1 \rangle_S - \sqrt{\frac{4}{5}} |3 1 1 \rangle_S


For \delta = 0.3

\displaystyle|0 0 0 \rangle_D = 0.995|0 0 0 \rangle_S + 0.099|2 2 0 \rangle_S - 0.015|2 0 0 \rangle_S + 0.009|4 4 0 \rangle_S + ...

\displaystyle|1 1 0 \rangle_D = 0.987|1 1 0 \rangle_S + 0.132|3 1 0 \rangle_S - 0.014|5 3 0 \rangle_S + 0.014|5 5 0 \rangle_S + ...

\displaystyle|1 0 1 \rangle_D = 0.994|1 1 1 \rangle_S + 0.108|3 3 1 \rangle_S + 0.017|3 1 1 \rangle_S + 0.011|5 5 1 \rangle_S + ...

\displaystyle|2 2 0 \rangle_D = 0.790|2 2 0 \rangle_S - 0.564|2 0 0 \rangle_S - 0.153|4 2 0 \rangle_S + 0.139|4 4 0 \rangle_S + ...

\displaystyle|2 1 1 \rangle_D = 0.986|2 2 1 \rangle_S + 0.158|4 4 1 \rangle_S - 0.052|4 2 1 \rangle_S - 0.012|6 4 1 \rangle_S + ...


For \delta = 0.6

\displaystyle|0 0 0 \rangle_D = 0.958|0 0 0 \rangle_S + 0.258|2 2 0 \rangle_S - 0.084|2 0 0 \rangle_S + 0.061|4 4 0 \rangle_S + ...

\displaystyle|1 1 0 \rangle_D = 0.885|1 1 0 \rangle_S + 0.319|3 3 0 \rangle_S - 0.274|3 1 0 \rangle_S - 0.119|5 3 0 \rangle_S + ...

\displaystyle|1 0 1 \rangle_D = 0.955|1 1 1 \rangle_S + 0.281|3 3 1 \rangle_S + 0.078|5 5 1 \rangle_S - 0.043|5 3 1 \rangle_S + ...

\displaystyle|2 2 0 \rangle_D = 0.598|2 2 0 \rangle_S - 0.45|2 0 0 \rangle_S - 0.379|4 2 0 \rangle_S + 0.299|4 4 0 \rangle_S \\ - 0.259|0 0 0 \rangle_S + 0.198|4 0 0 \rangle_S + 0.173|6 2 0 \rangle_S - 0.168|6 4 0 \rangle_S + ...

\displaystyle|2 1 1 \rangle_D = 0.882|2 2 1 \rangle_S + 0.381|4 4 1 \rangle_S - 0.191|4 2 1 \rangle_S - 0.129|6 6 1 \rangle_S \\ - 0.114|6 4 1 \rangle_S + ...

We can see, more deform, more higher angular momentum states are involved.

Also, for a pure state when non-deform, the mixing is still small.

decompo.PNG

Axial Harmonic Oscillator – Nilsson Orbit

Leave a comment

The Hamiltonian is

\displaystyle H = -\frac{\hbar}{2m}\nabla^2 +\frac{m}{2}(\omega_\rho^2(x^2+y^2)+\omega_z z^2 )

Use cylindrical coordinate, the Schrodinger equation is

\displaystyle \left(-\frac{\hbar}{2m}\left(\frac{d^2}{dz^2} + \frac{1}{\rho}\frac{d}{d\rho}\left(\rho\frac{d}{d\rho}\right) + \frac{1}{\rho^2}\frac{d^2}{d\phi^2}\right) +\frac{m}{2}(\omega_\rho^2(\rho^2)+\omega_z z^2 ) \right) \Psi = E \Psi

Where the energy is

\displaystyle E = \hbar\omega_\rho(1+ n_x + n_y) + \hbar\omega_z\left(\frac{1}{2}+ n_z\right)

Note that n_x+n_y \neq n_\rho. One of the reason is there are 2 degree of freedom, it cannot be solely expressed into 1 parameter.

Set \Psi = Z P \Phi ,

\displaystyle -\frac{\hbar}{2m}\left(\frac{d^2Z}{dz^2} P \Phi  + \frac{1}{\rho}\frac{d}{d\rho}\left(\rho\frac{dP}{d\rho} \right) Z \Phi + \frac{1}{\rho^2}\frac{d^2\Phi}{d\phi^2} Z P\right)  \\ +  \frac{m}{2}\omega_\rho^2\rho^2 ZP\phi+\frac{m}{2}\omega_z z^2 ZP\Phi= E ZP\Phi

\displaystyle -\frac{\hbar}{2m}\left(\frac{d^2Z}{dz^2} \frac{1}{Z}  + \frac{1}{\rho}\frac{d}{d\rho}\left(\rho\frac{dP}{d\rho} \right) \frac{1}{P} + \frac{1}{\rho^2}\frac{d^2\Phi}{d\phi^2} \frac{1}{\Phi}\right)  \\ +  \frac{m}{2}\omega_\rho^2\rho^2 +\frac{m}{2}\omega_z z^2 = E 

The angular part, we can set

\displaystyle \frac{d^2\Phi}{d\phi^2} \frac{1}{\Phi} = -m_\phi^2

The solution is \Phi(\phi) = \exp(i m_\phi \phi)

The z-part is usual 1D Harmonic oscillator

Thus the rest is

\displaystyle -\frac{\hbar}{2m}\left( \frac{1}{\rho}\frac{d}{d\rho}\left(\rho\frac{dP}{d\rho} \right) \frac{1}{P} - \frac{m_\phi^2}{\rho^2}\right) +  \frac{m}{2}\omega_\rho^2\rho^2 =  \hbar\omega_\rho(1+ n_x + n_y)

rearrange

\displaystyle -\frac{\hbar}{2m}\left( \frac{1}{\rho}\frac{d}{d\rho}\left(\rho\frac{dP}{d\rho} \right) - \frac{m_\phi^2}{\rho^2}P\right) +  \frac{m}{2}\omega_\rho^2\rho^2 P - \hbar\omega_\rho(1+ n_x + n_y) P = 0

Set a dimensionless constant \alpha^2 = \frac{\hbar}{m \omega_\rho} , and x = \rho/\alpha

\displaystyle \frac{1}{x}\frac{d}{dx}\left(x\frac{dP}{dx}\right) + \left( 2(1+n_x+n_y) - \frac{m_\phi^2}{x^2} - x^2 \right) P = 0

Using the normalisation formula, \int P^2 x dx = 1, set u = P \sqrt{x}

\displaystyle x^2\frac{d^2u}{dx^2}+ \left( 2(1+n_x+n_y) - \frac{m_\phi^2-\frac{1}{4}}{x^2} - x^2 \right) u = 0

Because of the long-range and short-range behavior, set

\displaystyle  u(x) = f(x) x^{m_\phi+\frac{1}{2}} \exp\left(-\frac{x^2}{2}\right)

\displaystyle x\frac{d^2f}{dx^2} + (1+2m_\phi-2x^2)\frac{df}{dx} + 2(n_x+n_y-m_\phi)x f= 0

Set y = x^2

\displaystyle y \frac{d^2f}{dy^2} + (m_\phi+1-y)\frac{df}{dy} + \frac{n_x+n_y-m_\phi}{2} f= 0

Define n_\rho = \frac{n_x+n_y-m_\phi}{2} \rightarrow n_x+n_y = 2 n_\rho + m_\phi = N - n_z.

\displaystyle y \frac{d^2f}{dy^2} + (m_\phi+1-y)\frac{df}{dy} + n_\rho f= 0

This is our friend again! The complete solution is

\displaystyle \Psi_{n_z n_\rho m_\phi}(z, \rho, \phi) \\ = \sqrt{\frac{1}{\alpha_z \alpha^2}}\sqrt{\frac{ n_\rho !}{2^n_z  n_z! (m_\phi + n_\rho)!\sqrt{\pi^3}}} H_{n_z}\left(\frac{z}{\alpha_z}\right) \\ \exp\left(- \frac{1}{2}\left(\frac{z^2}{\alpha_z^2}+\frac{\rho^2}{\alpha^2}\right)\right) \left(\frac{\rho}{\alpha}\right)^{m_\phi} L_{n_\rho}^{m_\phi} \left(\frac{\rho^2}{\alpha^2}\right) \exp(i m_\phi \phi)

where \alpha_z^2 = \hbar/m/\omega_z

The energy is

\displaystyle E = \hbar\omega_z \left(\frac{1}{2} + n_z \right) + \hbar\omega_\rho \left( 2n_\rho + m_\phi + 1 \right)

The quantum number m_\phi has same meaning as m_l in spherical case.

The notation for the state is

|Nn_z m_l \rangle

with the spin, the only good quantum number is the z-component of the total angular momentum J = L+S  long the body axis, i.e. K = m_\phi \pm 1/2 , thus, the state is

|Nn_z m_l K \rangle

Note that the total angular momentum J^2 is not a good quantum number in deformation, as the rotational symmetry is lost. However, the quantum number K is linked with the angular momentum of the Nilsson single particle orbit. It is because when a particle has m_j = K, the angular momentum must at least j \geq K .


The above is a general solution for the harmonic oscillator in cylindrical coordinate. When \omega_z = \omega_\rho , it reduce to spherical case.

According to P. Ring & P. Schuck (2004) (The Nuclear Many-Body Problem, P.68), the Hamiltonian can be expressed as a quadruple deform field by setting

\displaystyle \omega_\rho^2 = \omega^2 \left(1+\frac{2}{3} \delta\right) \\ \omega_z^2 = \omega^2 \left(1-\frac{4}{3} \delta\right)

\displaystyle H = -\frac{\hbar^2}{2m}\nabla^2 + \frac{1}{2} m \omega^2 \left(r^2- \frac{2}{3}\sqrt{\frac{16\pi}{5}} \delta r^2 Y_{20}(\theta,\phi) \right)

This Hamiltonian has similarity with the deformation

\displaystyle R(\theta, \phi) = R_0 \left(1 + \alpha_00 + \sum_{\lambda=1}^{\infty}\sum_{\mu=-\lambda}^{\lambda} \alpha_{\lambda \mu} Y_{\lambda \mu}(\theta, \phi) \right)

take the first quadruple term, and calculate R^2(\theta, \phi)

R^2(\theta, \phi) = R_0^2 (1 + 2 \beta Y_{20}(\theta,\phi) )

Compare, we have

\displaystyle  \beta = \frac{1}{3}\sqrt{\frac{16\pi}{5}}\delta \approx 1.05689 \delta

delta.png

The ratio

\displaystyle \frac{\omega_z}{\omega_\rho} =  \left(\frac{\alpha_\rho}{\alpha_z} \right)  = \sqrt{\frac{1-\frac{4}{3}\delta}{1+\frac{2}{3}\delta}} = 1- \delta + \frac{1}{6}\delta^2 - \frac{5}{18}\delta^3... 

We need volume conservation

\omega_x \omega_y \omega_z = \omega_0^3

Thus,

\displaystyle \omega = \omega_0 \left(\frac{1}{(1+2/3\delta)(1-4/3\delta)} \right)^{\frac{1}{6}} \approx \omega_0 \left(1+\frac{2}{9}\delta^2 + \frac{8}{81}\delta^3 \right)

The energy is

\displaystyle E = \hbar\omega \left( \sqrt{1-\frac{4}{3}\delta}\left(\frac{1}{2} + n_z \right) + \sqrt{1+\frac{2}{3}\delta} \left( N - n_z + 1 \right) \right) \\ \approx \hbar \omega_0 \left( \frac{3}{2} + \left(1 + \frac{1}{3}\delta\right) N - \delta n_z \right)

Nilsson_noLS.PNG

Here is some plots with various \delta = 0, 0.3, 0.5

Nilsson_beta00.PNG

Nilsson_beta03.PNG

Nilsson_beta05.PNG

Next time, I will add the LS coupling and L-square term to recreate the Nilsson diagram. Also I will expand the solution from cylindrical coordinate into spherical coordinate. This unitary transform is the key to understand the single particle-ness.