Oscylacje neutrin to postawowe zjawisko, w którym możemy obserwowć szerg interesujących włąsności tych tajemniczych cząstek. Zakładamy, że możliwy jest fenomenolobicxzny opis zaczepnięty z prac F.Benattiego i R, Floreanniniego:
$ H=H_0+ H_A$
gdzie część $H_0$ opisuje oscylacje sowobodne, zaś $H_A$ oddziaływanie z materią. Założymy, że możemy się ograniczyć do rozważania neutrin mionowych i elektronowych, wówczas:
$H_0=\left(\begin{array}{cc}E-\omega_0-\omega_3 & \omega_1-i\omega_2 \\ \omega_1 +i\omega_2 & E+\omega_0+\omega_3\end{array}\right)$
gdzie $E$ to średnia energia neutrin, $\omega_0\sim \Delta m^2$ opisuje różnicę mas różnych typów neutrin. Pozostaę parametry opisują oddziaływanie z otoczeniem, ale dla uproszczenia założymy ich znikanie.
Część hamiltonianu opisująca oddziaływanie z materią:
$H_A= \frac{A}{2}\left(\begin{array}{cc} 1+\cos(2\theta) & e^{-i\phi}\sin(2\theta)\\ e^{i\phi}\sin(2\theta) & 1-\cos(2\theta)\end{array} \right) $
zależy od $A$ reprezentujacego efektywną "gestość materii" ośrodka, w którym neutrina propagują, zaś $\theta$ to tzw. kąt mieszania, parametr znany z coraz liczniejszych eksperymentów.
Początkowo zakładamy neutrina elektronowe w stanie:
$\rho_{\nu_e}=\left(\begin{array}{cc} \cos^2(\theta) & e^{-i\phi}\cos(\theta)\sin(\theta) \\ e^{i\phi}\cos(\theta)\sin(\theta) & \sin^2(\theta)\end{array}\right) $
odpowiadający im stan neutrin mionowych otrzymujemy z zależności:
$\rho_{\nu_\mu}=1-\rho_{\nu_e}$
Paramert $\phi$ jest niezerowy w przpadku, gdy neutrina miałyby być cząstkami typu Majorany.
from qutip import *
from pylab import *
%pylab inline
Populating the interactive namespace from numpy and matplotlib
WARNING: pylab import has clobbered these variables: ['power', 'linalg', 'draw_if_interactive', 'random', 'save', 'load', 'info', 'fft'] `%pylab --no-import-all` prevents importing * from pylab and numpy
Wprowadźmy oznaczenia:
up=basis(2,0)
dn=basis(2,1)
rpp=up*up.dag()
rpm=up*dn.dag()
rmp=rpm.trans()
rmm=dn*dn.dag()
Określmy parametry oscylacji:
E=1.0
w0=1.0
w3=1.0
w2=0.0
w1=0.0
Hamiltonian:
H0=(E-w0-w3)*rpp+(E+w0+w3)*rmm+(w1-1j*w2)*rpm+(w1+1j*w2)*rpm
Dalsze parametry:
A=1.0
theta=pi/3.0
phi=0.0
HA=A/2.0*((1.0+cos(2.0*theta))*rpp+(1.0-cos(2.0*theta))*rmm+exp(-1j*phi)*sin(2.0*theta)*rpm+exp(1j*phi)*sin(2.0*theta)*rmp)
i cały hamiltonian:
H=H0+HA
a także stan początkowy:
psi0=cos(theta)**2.0*rpp+sin(theta)**2.0*rmm+cos(theta)*sin(theta)*(exp(-1j*phi)*rpm+exp(1j*phi)*rmp)
psi_mu=qeye(2)-psi0
Nie wykluczajmy również obecności dekoherencji wynikającej z wpływu otoczenia:
D=sqrt(0.1)*sigmaz()
tlist=linspace(0,5.0,10.0)
outlist=zeros(len(tlist))
obliczmy dynamikę neutrin:
qq=mesolve(H, psi0, tlist, [D], [])
Charakterystyką badaną będzie prawdopodobieństwo oscylacji:
$P_{\nu_{e}\rightarrow \nu_{\mu}}=\mbox{Tr}(\rho_{\nu_{e}}(t)\rho_{\nu_{\mu}})$
for i1 in range(len(tlist)):
rho=qq.states[i1]
outlist[i1]=abs((rho*psi_mu).tr())
Przedstawmy wizualizację otrzymanego wyniku:
plot(tlist,outlist,'-',label='', linewidth=4)
xlabel('t',fontsize=20)
ylabel('P(nu->mu)',fontsize=20)
show()