In [1]:
# SCHRODINGER EQUATION FOR 1-ELECTRON ATOM
# approximations: fixed nucleus, 1 electron, spherical symmetry...
# 
# Schrodinger equation takes the usual form: Hamiltonian[Psi] = E*Psi
# in this case we KNOW the Hamiltonian and Psi = orbitals from textbooks
# spherical symmetry => solutions must be spherically symmetric = 1s,2s,3s...
# 
# All we do below is:
# (a) take/define known Hamiltonian and Psi
# (b) build Schrodinger equation: Hamiltonian[Psi] = E*Psi
# (c) solve Schrodinger equation for E => get energies (E) of orbitals (1s,2s,3s..)
#
# Derivation of Hamiltonian for spherically symmetric 1-electron atom
# => textbooks and also my handwritten notes in JPG-files
#
# Using known solutions = Psi = orbitals 1s,2s,3s... is not a nonsense
# => theoretically, we could guess to solutions
# => guessing the solution of differential equation is entirely Ok
#    ...although in this case the guess would be quite difficult
# => simple analogy with common equation:
#    ...our standard equation: x+2 = 3
#    ...solutions of standard equations are numbers
#    ...we guess intuitively that the solution is: x=1
# => the same can be done with differential equations
#    ...our differential equation: d/dx[f(x)] = f(x)
#    ...solutions of differential equations are functions
#    ...we guess intuitively that the solution is: f(x)=exp(x)
In [2]:
import sympy as sp
In [3]:
sp.init_printing()
In [4]:
# Properties of symbols/functions are very important
# (r + h,m,e,Z: real and positive
# (E: real and NEGATIVE
# (F: real
(r) = sp.symbols('r', real=True, positive=True)
(h,m,e,Z) = sp.symbols('h m_e e Z', real=True, positive=True)
(E) = sp.symbols('E', real=True)
(F) = sp.Function('\psi', real=True)
In [5]:
# (1) DEFINE SCHRODINGER EQUATION = SchrEq
# (SchrEq: Hamiltonian[Psi] = E*Psi
# (we use known Hamiltonian for: fixed nucleus, 1 electron, spherical symmetry
In [6]:
SchrEq = sp.Eq(-h**2/(2*m) * (1/r**2) * sp.diff( r**2 * sp.diff(F(r),r), r) - Z*e**2/r * F(r),E*F(r))
In [7]:
display(SchrEq)
$\displaystyle - \frac{Z e^{2} \psi{\left(r \right)}}{r} - \frac{h^{2} \left(r^{2} \frac{d^{2}}{d r^{2}} \psi{\left(r \right)} + 2 r \frac{d}{d r} \psi{\left(r \right)}\right)}{2 m_{e} r^{2}} = E \psi{\left(r \right)}$
In [8]:
# (2) DEFINE SOLUTIONS of SchrEq
# (solutions of our Schr.eq. = functions Psi = orbitals 1s,2s,3s...
# (the formulas for orbitals 1s,2s,3s... can be found in QM textbooks and www
In [9]:
orbital_1s = sp.exp(-Z*e**2*m*r/h**2)
In [10]:
orbital_2s = (2 - Z*e**2*m*r/h**2) * sp.exp(-Z*e**2*m*r/(2*h**2))
In [11]:
a = sp.symbols('a')
a = h**2/(m*e**2)
orbital_3s = (27 - 18*(Z*r/a) + 2*(Z**2*r**2/a**2)) * sp.exp(-(Z*r)/(3*a))
In [12]:
display(orbital_1s)
$\displaystyle e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}$
In [13]:
display(orbital_2s)
$\displaystyle \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}$
In [14]:
display(orbital_3s)
$\displaystyle \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}$
In [15]:
# (3) INSERT KNOWN SOLUTIONS into SchrEq
In [16]:
SchrEq_1s = SchrEq.subs(F(r),orbital_1s)
In [17]:
SchrEq_2s = SchrEq.subs(F(r),orbital_2s)
In [18]:
SchrEq_3s = SchrEq.subs(F(r),orbital_3s)
In [19]:
display(SchrEq_1s)
$\displaystyle - \frac{Z e^{2} e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}}{r} - \frac{h^{2} \left(r^{2} \frac{\partial^{2}}{\partial r^{2}} e^{- \frac{Z e^{2} m_{e} r}{h^{2}}} + 2 r \frac{\partial}{\partial r} e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}\right)}{2 m_{e} r^{2}} = E e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}$
In [20]:
display(SchrEq_2s)
$\displaystyle - \frac{Z e^{2} \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}}{r} - \frac{h^{2} \left(r^{2} \frac{\partial^{2}}{\partial r^{2}} \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}} + 2 r \frac{\partial}{\partial r} \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}\right)}{2 m_{e} r^{2}} = E \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}$
In [21]:
display(SchrEq_3s)
$\displaystyle - \frac{Z e^{2} \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}}{r} - \frac{h^{2} \left(r^{2} \frac{\partial^{2}}{\partial r^{2}} \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}} + 2 r \frac{\partial}{\partial r} \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}\right)}{2 m_{e} r^{2}} = E \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}$
In [22]:
# (4) PERFORM DERIVATIVES in SchrEq
# (do all derivatives in SchrEq's for 1s,2s,3s...
In [23]:
SchrEq_1s = SchrEq_1s.doit(); display(SchrEq_1s)
$\displaystyle - \frac{Z e^{2} e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}}{r} - \frac{h^{2} \left(\frac{Z^{2} e^{4} m_{e}^{2} r^{2} e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}}{h^{4}} - \frac{2 Z e^{2} m_{e} r e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}}{h^{2}}\right)}{2 m_{e} r^{2}} = E e^{- \frac{Z e^{2} m_{e} r}{h^{2}}}$
In [24]:
SchrEq_2s = SchrEq_2s.doit(); display(SchrEq_2s)
$\displaystyle - \frac{Z e^{2} \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}}{r} - \frac{h^{2} \left(\frac{Z^{2} e^{4} m_{e}^{2} r^{2} \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 6\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}}{4 h^{4}} + 2 r \left(- \frac{Z e^{2} m_{e} \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}}{2 h^{2}} - \frac{Z e^{2} m_{e} e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}}{h^{2}}\right)\right)}{2 m_{e} r^{2}} = E \left(- \frac{Z e^{2} m_{e} r}{h^{2}} + 2\right) e^{- \frac{Z e^{2} m_{e} r}{2 h^{2}}}$
In [25]:
SchrEq_3s = SchrEq_3s.doit(); display(SchrEq_3s)
$\displaystyle - \frac{Z e^{2} \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}}{r} - \frac{h^{2} \left(\frac{Z^{2} e^{4} m_{e}^{2} r^{2} \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{9 h^{4}} - \frac{14 Z e^{2} m_{e} r}{3 h^{2}} + 19\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}}{h^{4}} + 2 r \left(- \frac{Z e^{2} m_{e} \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}}{3 h^{2}} + \left(\frac{4 Z^{2} e^{4} m_{e}^{2} r}{h^{4}} - \frac{18 Z e^{2} m_{e}}{h^{2}}\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}\right)\right)}{2 m_{e} r^{2}} = E \left(\frac{2 Z^{2} e^{4} m_{e}^{2} r^{2}}{h^{4}} - \frac{18 Z e^{2} m_{e} r}{h^{2}} + 27\right) e^{- \frac{Z e^{2} m_{e} r}{3 h^{2}}}$
In [26]:
# (5) SOLVE SchrEq
# (i.e. solve SchrEq's for E for 1s,2s,3s...
# (i.e. find energy of electrons in orbitals 1s,2s,3s...
# (note that the solution of rather complex equations above is quite simple...
In [27]:
SchrEq_1s_sol = sp.solve(SchrEq_1s,E)
In [28]:
display(SchrEq_1s_sol)
$\displaystyle \left[ - \frac{Z^{2} e^{4} m_{e}}{2 h^{2}}\right]$
In [29]:
SchrEq_2s_sol = sp.solve(SchrEq_2s,E)
In [30]:
display(SchrEq_2s_sol)
$\displaystyle \left[ - \frac{Z^{2} e^{4} m_{e}}{8 h^{2}}\right]$
In [31]:
SchrEq_3s_sol = sp.solve(SchrEq_3s,E)
In [32]:
display(SchrEq_3s_sol)
$\displaystyle \left[ - \frac{Z^{2} e^{4} m_{e}}{18 h^{2}}\right]$
In [33]:
# (6) PRINT ALL SOLUTIONS AND GENERALIZE...
# (the solutions in (5) suggest that E = konst * Z**2/n**2
In [34]:
# 1st step: Prepare calculated solutions for energy of orbitals 1s,2s,3s... 
# (take just the first (and the only) solutions of all equations
if type(SchrEq_1s_sol) == list: SchrEq_1s_sol = SchrEq_1s_sol[0]
if type(SchrEq_2s_sol) == list: SchrEq_2s_sol = SchrEq_2s_sol[0]
if type(SchrEq_3s_sol) == list: SchrEq_3s_sol = SchrEq_3s_sol[0]
In [35]:
# 2nd step: Define [generalized solution] according to results in (5)...
# ...we note that the solutions are very similar and can be generalized
# ...we define symbol n = principal quantum number (n=1 for 1s, 2 for 2s...)
n = sp.symbols('n', positive=True, integer=True)
# ...and we define formula for general equation for energy E(n)
def E(n): return -(Z**2 * e**4 * m) / (2 * n**2 *h**2)
# ...display the general solution
display([E(n)])
$\displaystyle \left[ - \frac{Z^{2} e^{4} m_{e}}{2 h^{2} n^{2}}\right]$
In [36]:
# 3rd step: Compare the calculated solutions for orbitals 1s,2s,3s...
# with our general relation for E for n=1,2,3...
print('Calculated solutions for 1s,2s,3s ...')
display([SchrEq_1s_sol, SchrEq_2s_sol, SchrEq_3s_sol])
print('Generalized solution E = -const * Z**2/n**2 ...')
display([E(1), E(2), E(3)])
if ([SchrEq_1s_sol, SchrEq_2s_sol, SchrEq_3s_sol] == [E(1), E(2), E(3)]):
    print('Calculated and generalized solutions are equivalent.')

Calculated solutions for 1s,2s,3s ...
$\displaystyle \left[ - \frac{Z^{2} e^{4} m_{e}}{2 h^{2}}, \ - \frac{Z^{2} e^{4} m_{e}}{8 h^{2}}, \ - \frac{Z^{2} e^{4} m_{e}}{18 h^{2}}\right]$ 
Generalized solution E = -const * Z**2/n**2 ...
$\displaystyle \left[ - \frac{Z^{2} e^{4} m_{e}}{2 h^{2}}, \ - \frac{Z^{2} e^{4} m_{e}}{8 h^{2}}, \ - \frac{Z^{2} e^{4} m_{e}}{18 h^{2}}\right]$ 
Calculated and generalized solutions are equivalent.
In [37]:
# (7) FINALIZE THE FORMULA FOR GENERAL SOLUTION AND CALCULATE CONSTANT
In [38]:
# Write final general formula for energy of electron in 1-eln atom 
# (introduce undefined function for E: for equation, we need...
# (...LHS: undefined function: declared below as general SciPy by sp.Function
# (...RHS: defined function E: declared above as standard Python function
Energy     = sp.Function('E', real=True)
Energy_fin = sp.Eq(Energy(n),E(n))
display(Energy_fin)
$\displaystyle E{\left(n \right)} = - \frac{Z^{2} e^{4} m_{e}}{2 h^{2} n^{2}}$
In [39]:
# Rewrite the formula so that all constant values were combined in one const 
const = sp.symbols('const')
Energy_fin.subs(((e**4 * m) / (2 * h**2)),const)
Out[39]:
$\displaystyle E{\left(n \right)} = - \frac{Z^{2} const}{n^{2}}$
In [40]:
# Calculate the value of constant
# (the constant should be Rydberg unit of energy = 13.6 eV
# (more info: https://en.wikipedia.org/wiki/Rydberg_constant
In [41]:
# Prepare constants
# (eps = permitivity of free space
# (e0 = charge of electron including permitivity
# (me = mass of electron
# (hr = Planck constant / (2*pi)
from math import pi,sqrt
eps = 8.854e-12
e0  = 1.602e-19/sqrt(4*pi*eps)
me  = 9.109e-31
hr  = 6.626e-34/(2*pi)
In [42]:
# Calculate numerial value of constat = Rydberg unit of energy = Ry
# (Ry = symbolic formula
# (Ry_value = numerical value in J
# (Ry_value_eV = Ry_value / electron_charge
# (E[J] = QU[CV] => Q/e => charge in electron units => E[CV] -> E[eV]
Ry = (e**4 * m) / (2 * h**2)
Ry_value = Ry.subs([(e,e0),(m,me),(h,hr)])
Ry_value_eV = Ry_value/1.602e-19
sp.Eq(Ry,Ry_value_eV)
Out[42]:
$\displaystyle \frac{e^{4} m_{e}}{2 h^{2}} = 13.601485784325$
In [43]:
# (8) CALCULATE ENERGY OF EDX TRANSITION USING ABOVE-DERIVED RELATIONS
# (we take one typical/frequent EDX transition: CuKa = 8.04 keV
# (objective: verify that our solution gives reasonable results
In [44]:
def E(Z,n): return -Ry_value_eV * Z**2 / n**2
In [45]:
E_CuKa_eV  = E(29,2)-E(29,1)
E_CuKa_keV = E_CuKa_eV / 1000
In [46]:
print('E(CuKa), tabulated : %.2f [keV]' % 8.04)
print('E(CuKa), calculated: %.2f [keV]' % E_CuKa_keV)

E(CuKa), tabulated : 8.04 [keV]
E(CuKa), calculated: 8.58 [keV]
In [47]:
# SUMMARY
# 1) We have verified solution of Schrodinger equation
#  ..simplified version of Schr.eq: for 1-eln atom, spherically symmetric
#  ..verificatin by insertion of known spherically symmetric solutins: orbitals 1s,2s,3s...
# 2) The verification yielded the following...
#  ..confirmation that our simplified form of Schr.eq is correct (otherwise it would not work)
#  ..expression for energy of orbitals 1s,2s,3s... (for 1-eln atom it holds: E(2s)=E(2p))
#  ..generalized expression for energy of orbitals (for 1-eln atom it holds: E = -const*Z**2/n**2)
#  ..numerical value of the constant: E = -const*Z**2/n**2 = -13.6[eV]*Z**2/n**2
# 3) Finally, we have used the relation E = -13.6[eV]*Z**2/n**2...
#  ..to calculate energy of real EDX transition = E(CuKa)
#  ..despite all aproximations, the tabulated and calculated values were in good agreement