シュレーディンガー方程式：たまの週末：SSブログ

	ブログをはじめるログイン

シュレーディンガー方程式　[科学、数学] [編集]

シュレーディンガー方程式(Schrödinger equation) $\Psi\left(\overrightarrow{x},t\right)$ は、以下のように表されるらしいです。
$\\ -\cfrac{\hbar~2}{2m}\nabla^2\Psi\left(\overrightarrow{x},t\right)+V\left(\overrightarrow{x}\right)\Psi\left(\overrightarrow{x},t\right)=\mathrm{i}\hbar\cfrac{\partial\Psi\left(\overrightarrow{x},t\right)}{\partial t}$
$\\ \\$

この方程式は、位置 $\left(\overrightarrow{x}\right)$ に依存する式と、時間 $\left(t\right)$ に依存する式にわけることができ、位置に依存する（時間に依存しない:time-independent Schrödinger equation; TISE）式は次のようになるそうです。
$\\ \\$

$-\cfrac{\hbar~2}{2m}\nabla^2\Psi\left(\overrightarrow{x}\right)+V\left(\overrightarrow{x}\right)\Psi\left(\overrightarrow{x}\right)=E\Psi\left(\overrightarrow{x}\right)$
$\\ \\$

方程式を解きやすくするための「変数分離法」というテクニックだそうです。今まで読んだ教科書では、このあたりはアプリオリにサラッと書かれていてよくわからなかったのですが、先日買ったEMANさんの本「趣味で量子力学」やWEBサイトには丁寧な説明がありました。

本を読むだけだと、今一つわかった気がしないので、このページ
https://en.wikiversity.org/wiki/User:Tclamb/FEniCS
を参考に、Pythonで、方程式を解いてみました。

from dolfin import *

#define mesh and function space
mesh = IntervalMesh(100, -10, 10)
V = FunctionSpace(mesh, 'Lagrange', 3)
 
#define boundary
def boundary(x, on_boundary):
    return on_boundary

#apply essential boundary conditions
bc = DirichletBC(V, 0, boundary)
 
#define functions
u = TrialFunction(V)
v = TestFunction(V)

x2 = Expression('x[0]*x[0]')
 
#define problem
a = (inner(grad(u), grad(v)) + x2*u*v)*dx
m = u*v*dx

#assemble stiffness matrix
A = PETScMatrix()
assemble(a, tensor=A)
M = PETScMatrix()
assemble(m, tensor=M)

#create eigensolver
eigensolver = SLEPcEigenSolver(A,M)
eigensolver.parameters['spectrum'] = 'smallest magnitude'
eigensolver.parameters['solver']   = 'lapack'
eigensolver.parameters['tolerance'] = 1.e-15

#solve for eigenvalues
eigensolver.solve()

#extract first eigenpair
r, c, rx, cx = eigensolver.get_eigenpair(0)
print 'eigenvalue:', r

#assign eigenvector to function
u = Function(V)
u.vector()[:] = rx

#plot eigenfunction and probability density
plt  = plot(u, rescale=False, range_min=-1.0, range_max=1.0)
plt.write_png('img_%05d' % 0)

i = 1
while i < eigensolver.get_number_converged():
    #extract next eigenpair
    r, c, rx, cx = eigensolver.get_eigenpair(i)
    #assign eigenvector to function
    u.vector()[:] = rx
    #plot eigenfunction
    plt.update(u)
    plt.write_png('img_%05d' % i)
    i = i+1

シュレーディンガー方程式からさらに球面調和関数（spherical harmonics）を求めると、電子の軌道を求めることができるらしいです。
「ファインマン物理学 [ 5 ] 量子力学」には、次のような図が掲載されています。

ファインマン物理学が書かれた1965年当時と違って、今はパソコン使って手軽（？）に電子の軌道を図にすることができるようになりました。

QuTiP : Quantum Toolbox in Pythonを使った例

mayaviを使った電子軌道表示の例

QuTipの例

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import numpy as np
from qutip import *
from pylab import *
#
# plot angular wave function for l=3
#
phi=linspace(0,2*pi,90)
theta=linspace(0,pi,45)

c2=basis(7,4) #2l+1
figure(figsize=[6,4])
y=orbital(theta,phi,c2)
sphereplot(theta,phi,y)
show()

mayaviで電子軌道の表示

from mayavi import mlab
import numpy as np
from scipy.special import sph_harm

# Create a sphere
r = 0.3
pi = np.pi
cos = np.cos
sin = np.sin
phi, theta = np.mgrid[0:pi:101j, 0:2 * pi:101j]

x = r * sin(phi) * cos(theta)
y = r * sin(phi) * sin(theta)
z = r * cos(phi)

mlab.figure(1, bgcolor=(1, 1, 1), fgcolor=(0, 0, 0), size=(400, 300))
mlab.clf()
for n in range(1, 6):
    for m in range(n):
        s = sph_harm(m, n, theta, phi).real
        s[s < 0] *= 0.97
        s /= s.max()
        mlab.mesh(s * x - m, s * y - n, s * z + 1.3,scalars=s, colormap='Spectral')

mlab.show()