In [1]:

    

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from starkhelium import *
from tqdm import trange, tqdm
import os
from scipy.constants import h, hbar, c, alpha, m_e, e, epsilon_0, atomic_mass, pi, physical_constants
a_0 = physical_constants['Bohr radius'][0]
E_h = physical_constants['Hartree energy'][0]


    

In [3]:

    

# Transition dipole moments
n_1, l_1, m_1 = 70, 69, 69
n_2, l_2, m_2 = 71, 70, 70
transition_dipole_moment = np.abs(
    stark_int(n_1, n_2, l_1, l_2, m_1, m_2, field_orientation='parallel', dm_allow=[-1,0,+1])) * e*a_0 # in atomic units e a_0
print('Transition dipole moment: ', transition_dipole_moment/(e*a_0), ' (a.u.)')

# Spontaneous decay rate
nmin = np.min([n_1, n_2])
nmax = np.max([n_1, n_2])
S = 1
m = 0
n_vals, L_vals = get_nl_vals(nmin, nmax, m)
J_vals = get_J_vals(S, L_vals, 1)
# energy levels
En = W_n(S, n_vals, L_vals, J_vals)

E_upper = En[np.intersect1d(np.where(n_vals == n_1), np.where(L_vals == l_1))[0]] * E_h
E_lower = En[np.intersect1d(np.where(n_vals == n_2), np.where(L_vals == l_2))[0]] * E_h
omega_transition = 2*pi*np.abs(E_upper - E_lower)/h # E=hf=h(w/2pi), w=2pi*E/h
einstein_A = ((2*omega_transition**3)/(3*epsilon_0*h*c**3)) * np.abs(transition_dipole_moment)**2
print('Spontaneous emission rate: ', einstein_A, ' /s')
print('Radiative lifetime, T1: ', (1/einstein_A) * 10**3, ' ms')


    

    




Transition dipole moment:  3476.95686325  (a.u.)
Spontaneous emission rate:  6.02070366054  /s
Radiative lifetime, T1:  166.093542613  ms

In [ ]: