Квантовые траектории и с чем их едят

22.11.2019

Эта небольшая заметка про то, как рисовать красивые картинки, ну и немного про физику, о которой редко говорят, про бомовскую квантовую механику.
Квантовые траектории и с чем их едят

Небольшое введение

Как нам любит рассказывать всякая научная фантастика и псевдонаучная ерунда, типа фильма «Секрет», законы микромира, очень сильно отличаются от привычных нам, классических.
В мире квантовой механики правит всем вероятность, задаваемая волновой функцией $inline$ (интересующиеся деталями могут заглянуть, например, в пост «Мюонный катализ с точки зрения квантовой химии. Часть I: обычный водород vs. мюонный водород»).
Из вероятностных свойств квантмеха растут ноги всяких забавных вещей, типа котов Шрёдингера, принципов неопределённости Гейзенберга и неравенств Белла.

Но все эти картинки со всякими орбиталями электрона не давали ответа на вопрос «как же электрон летит в пространстве». Чтобы прояснить эту ситуацию, физики потратили много времени, но так и не справились с этим. Зато Дэвид Бом (известный многим по эффекту Ааронова — Бома) окончательно создал один из формализмов квантовой механики (имени себя), в котором всё-таки есть траектории, по которым движется квантовая частица. И, в отличие от фейнмановских интегралов по траекториям, эта траектория у каждой частицы ровно одна. Это свойство принципиально позволяет отследить движение частиц, и сравнить движение классических и квантовых частиц, чем мы и займёмся в этой заметке.

не только формализм

Собственно, сам формализм никому особо не интересен, но из этого формализма можно построить одну из интерпретаций квантовой механики, которую, из-за кажущийся простоты классической механики любят некоторые фрики (не многие, ибо въехать в это дело не очень просто, да и известна/интересна эта часть квантмеха только особым извращенцам).
Эту (как и другие) интерпретацию мы обсуждать не будем.

Классические и квантовые траектории

Мы будем рассматривать достаточно скучную систему: один электрон в поле нескольких протонов. Об этой системе, а также о классической и квантовой механике можно почитать в первой и второй частях постов «Мюонный катализ с точки зрения квантовой химии».

Классическая задача о движении частицы в некотором потенциале даётся вторым законом Ньютона:

$display$

где m — масса частицы, x — координата, F — сила, действующая на частицу, а $ddot{x} = frac{d^2 x}{dt^2}$ — вторая производная координаты частицы по времени, или ускорение. Если в системе действуют только потенциальные силы, то силу можно выразить через новую сущность, потенциальную энергию V как

$F = - frac{dV}{dx}$

В нашем случае, электрона в поле нескольких протонов,

где электрон взаимодействует с каждым из протонов по закону Кулона

$display$

, где k — коэффициент, равный 1 в атомных единицах, e — заряд электрона, а R — расстояние от электрона до протона.
В этом случае суммарный потенциал, действующий на электрон будет равен

$V = sum_{n=1}^N V_n(R_n) = - sum_{n=1}^N frac{k e^2}{R_n} ,$

где индекс n нумерует протоны (всего протонов N штук), а R_n — расстояние от электрона до n-го протона.

Численно решить диффур, который представляет из себя второй закон Ньютона, задача халтурная, главное задать начальные положение и скорость. Если электрон будет лететь слишком быстро, то он вырвется за пределы притяжения протона(ов) и улетит на бесконечность, а если энергии будет чуть-чуть, то он будет вечно бултыхаться туды-сюды в поле одного из ядер, никогда не посещая другие.

Лучистое трение

Если бы мы учли лучистое трение, которое возникает из-за того, что при движении с ускорением, электрон будет часть своей энергии движения отдавать электромагнитному полю, излучая его куда-то, то электрон в итоге скатился бы на ядро за какое-то время.

Так что происходит в классике, мы знаем.

А что же будет в бомовской динамике?
В этом случае частица тоже будет двигаться по второму закону Ньютона с потенциалом $V = V_mathrm{C} + V_mathrm{Q}$ , где $V_mathrm{C}$ — классический потенциал из обычного закона Ньютона, который в нашем случае имеет вид, данный выше.
Т.е. помимо классического потенциала, на неё будет действовать ещё и другая сущность: квантовый потенциал $V_mathrm{Q}$ , имеющий (в 1D случае) вид

$V_mathrm{Q} = - frac{hbar^2}{2mA} frac{ d^2A}{ dx^2}$

где A — это амплитуда (модуль) волновой функции $inline$ ( $inline$ , где $inline$ — фаза волновой функции).
Значит, чтобы получить уравнение движения квантовой частицы, нам всё равно нужно знать что-то о волновой функции.

Про скрытые параметры

Формализм Бома является теорией со скрытыми параметрами. Но поскольку скрытый параметр (волновая функция) — нелокален, то результаты вычислений в этом формализме всё равно удовлетворяют вышеупомянутым неравенствам Белла.

В случае одного протона нам известно (см., например, тут) точное выражение волновой функции электрона в основном (1s) состоянии:

$display$

О нормировке

В формуле для квантового потенциала, нормировка числителя сократится со знаменателем, поэтому мы о ней париться и не будем.

В случае нескольких протонов, в рамках метода молекулярных орбиталей как комбинаций атомных орбиталей (МО ЛКАО, см. тут) основное состояние с достаточной степенью точности будет даваться суммой 1s-орбиталей каждого из атомов:

$psi approx sum_{n=1}^N psi_n(R_n) = sum_{n=1}^N exp(-R_n)$

Теперь, чтобы узнать квантовый потенциал, нужно всего-то воспользоваться этим выражением.

Нуджные вычисления

Функция

$inline$ как сумма 1s орбиталей — действительная, поэтому

$inline$ .
Поскольку электрон может двигаться в трёх измерениях, нужно одномерную производную

$A_{xx}'' = frac{d^2 A}{dx^2}$ заменить на её трёхмерное обобщение:

$Delta A = A_{xx}'' + A_{yy}'' + A_{zz}''$ . Оператор

$inline$ можно представить как квадрат оператора набла:

$inline$ . Также можно представить расстояние

$inline$ как

$R_n = sqrt{mathbf{R}_n^2}$ , где

$mathbf{R}_n$ — радиус-вектор электрона относительно n-го протона.

Тогда

$Delta A = nabla^2 psi = sum_{n=1}^N nabla^2 psi_n(R_n)$

Первая производная считается легко:

$nabla psi_n(R_n) = nabla exp(-R_n) = exp(-R_n) cdot (-1) cdot frac{1}{2 underbrace{sqrt{mathbf{R}_n^2}}_{R_n}} cdot 2 mathbf{R}_n = -exp(-R_n) cdot frac{mathbf{R}_n}{R_n}$

Вторая же производная уже несколько сложнее:

$nabla (nabla exp(-R_n)) = -frac{mathbf{R}_n}{R_n}nabla exp(-R_n) -exp(-R_n) nabla frac{mathbf{R}_n}{R_n} = exp(-R_n) - frac{2exp(-R_n) }{R_n}$

где $-frac{mathbf{R}_n}{R_n}nabla exp(-R_n) = exp(-R_n) cdot underbrace{ left(-frac{mathbf{R}_n}{R_n} right)^2}_{1} = exp(-R_n)$ и
$-exp(-R_n) nabla frac{mathbf{R}_n}{R_n} = -exp(-R_n) cdot left( frac{overbrace{nabla mathbf{R}_n}^{3}}{R_n} - frac{2mathbf{R}_n^2}{2 R_n^3} right) = - frac{2exp(-R_n)}{R_n}$ .
В итоге остаётся:
$Delta psi = overbrace{sum_{n=1}^{N} exp(-R_n)}^{psi} - sum_{n=1}^{N} frac{2exp(-R_n)}{R_n}$
Поделив всё на $inline$ и домножив на $-frac{hbar^2}{2m}$
получим
$V_mathrm{Q} = -frac{hbar^2}{2m} left( 1 - sum_{n=1}^N frac{2 exp(-R_n)}{R_n}right)$
Единица при дифференцировании для получения силы исчезнет, поэтому можно смело оставить только второе слагаемое.

В итоге наш квантовый потенциал мы можем записать как

$V_mathrm{Q} approx frac{hbar^2}{m} sum_{n=1}^N frac{ exp(-R_n)}{R_n}$

и с этим выражением мы можем уже гонять бомовскую динамику электрона в поле многих протонов.

Реализация

Для этого всего безобразия был написан код на питоне, он доступен тут:

Код на питоне

from math import *
import numpy as np

cutoff=5.0e-4
Quantum=True

def dist(r1,r2):
    return np.dot((r1-r2), (r1-r2))

def Vc(r, r0):
    if dist(r, r0)>=cutoff:
        return -1.0/dist(r, r0)
    else:
        return -1.0/cutoff


    
rH=[]
#h1
#rH.append(np.array([ 0.0, 0.0, 0.0]))

#h2
rH.append(np.array([-1.0, 0.0, 0.0]))
rH.append(np.array([+1.0, 0.0, 0.0]))




def Vat(r):
    V=0.0
    for rh in rH:
        V+=Vc(r,rh)
    return V

def PsiA(r):
    psi=0.0
    for rh in rH:
        if dist(r, rh)<1.0e3:
            psi+=np.exp(-dist(r, rh))
    return psi

def Vq(r):
    vq=0.0
    for rh in rH:
        if dist(r, rh)>=cutoff:
            vq-=2.0*np.exp(-dist(r, rh))/dist(r, rh)
        else:
            vq-=2.0*np.exp(-cutoff)/cutoff    

    vq*=(-0.5) # -0.5*hbar**2/me
    return vq

def GradF(F, r):
    grad=np.zeros(3)
    dx=0.1
    for i in range(0,3):
        dr=np.zeros(3)
        dr[i]=dx
        #print(dr)
        #print(F(r+dr)-F(r-dr))
        grad[i]+=(F(r+dr)-F(r-dr))/(2.*dx)
    return grad


    

dt=0.001
tmax=2.0e1
DR=1.0
dx=0.001

MaxR=10.0

t=0.0

cent=np.zeros(3)

Ntrj=30
m=1.0

def GenRvBox(DX):
    return np.random.uniform(-DX,+DX,3)

def GenRvSph(DX):
    r=np.random.uniform(0.0,DX)
    phi=np.random.uniform(0.0,2.0*np.pi)
    theta=np.random.uniform(0.0,np.pi)
    x=r*np.sin(theta)*np.cos(phi)
    y=r*np.sin(theta)*np.sin(phi)
    z=r*np.cos(theta)
    return np.array([x,y,z])

for ntrj in range(0,Ntrj):
    if Quantum:
        outf=open("bmd_%05i" % (ntrj) + ".trj", "w")
    else:
        outf=open("cmd_%05i" % (ntrj) + ".trj", "w")
    
    nat=np.random.randint(0,len(rH))
    r=rH[nat]+GenRvSph(DR)
    rprev=r+GenRvBox(dx)
    outf.write("%15.10f %15.10f %15.10fn" % tuple(r))
    t=0.0
    while t<=tmax and dist(r,cent)<=MaxR:
        F= -GradF(Vat, r)  
        if Quantum:
            F-= GradF(Vq, r)
        rnew= 2.*r - rprev + (F/m)*dt**2
        rprev=r
        r=rnew
        outf.write("%15.10f %15.10f %15.10fn" % tuple(r))
    
        t+=dt
    outf.close()

exit()

Мы же только обсудим несколько моментов.
Интегрируется второй закон Ньютона при помощи алгоритма Верле:

$x(t + Delta t) = 2 x(t) - x(t - Delta t) + frac{F(t)}{m} Delta t^2$

Начальное положение генерируется путём случайного выбора одного из протонов, вокруг него случайно выбирается направление (используя сферические координаты). Чтобы задать начальную скорость, нужно задать ещё одно, предыдущее положение. Оно выбирается ещё при помощи одного маленького случайного вектора.

Включая/выключая квантовый потенциал, мы переходим в квантовый/классический режимы движения.

Ну а дальше, можно построить красивые картинки, используя Gnuplot, для атома водорода

и для молекулы H₂⁺

Как видно, классические траектории (верхние, синие) или очень локализованы, или, если слишком быстро заставить двигаться электрон, убегают от ядер. В квантовом же случае (нижние, розовые), квантовый потенциал позволяет электронам гулять существенно дальше от ядра, а в случае молекулы H₂⁺, позволяет бегать от одного протона к другому, что является прямой визуализацией химических связей.

Пару слов о построении картинок: чтобы создать эффект неона, каждая траектория рисуется много раз, от тонкой белой, до толстой чёрной, через тени интересующего цвета. Для удобства выбора такой палитры можно, например, воспользоваться сайтом https://www.color-hex.com/
Пример скрипта приведён ниже.

Скрипт для Gnuplot

unset key

set xyplane relative 0
unset box 
set view map
set size ratio -1
unset border 

unset xtics

unset ytics
set terminal pngcairo size 2160,4096 backgr rgb "black"

set output "tmp.png"

yshift=-5.0 maxiC=29 maxiQ=29 splot for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 30.0 lc rgb "#030d19" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 18.0 lc rgb "#071b33" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 17.0 lc rgb "#0a294c" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 16.0 lc rgb "#0e3766" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 15.0 lc rgb "#11457f" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 14.0 lc rgb "#155399" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 13.0 lc rgb "#1861b2" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 12.0 lc rgb "#1c6fcc" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 11.0 lc rgb "#1f7de5" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 10.0 lc rgb "#238bff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 9.0 lc rgb "#3896ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 8. lc rgb "#4ea2ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 7. lc rgb "#65adff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 6. lc rgb "#7bb9ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 5. lc rgb "#91c5ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 4. lc rgb "#a7d0ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 3. lc rgb "#bddcff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 2. lc rgb "#d3e7ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 1. lc rgb "#e9f3ff" not, for [i=0:maxiC] sprintf("cmd_%05i.trj", i) w l lw 0.5 lc rgb "#ffffff" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 30.0 lc rgb "#190613" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 18.0 lc rgb "#330c27" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 17.0 lc rgb "#4c123b" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 16.0 lc rgb "#66184f" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 15.0 lc rgb "#7f1e63" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 14.0 lc rgb "#992476" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 13.0 lc rgb "#b22a8a" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 12.0 lc rgb "#cc309e" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 11.0 lc rgb "#e536b2" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 10.0 lc rgb "#ff3dc6" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 9.0 lc rgb "#ff50cb" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 8. lc rgb "#ff63d1" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 7. lc rgb "#ff77d7" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 6. lc rgb "#ff8adc" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 5. lc rgb "#ff9ee2" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 4. lc rgb "#ffb1e8" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 3. lc rgb "#ffc4ed" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 2. lc rgb "#ffd8f3" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 1. lc rgb "#ffebf9" not, for [i=0:maxiQ] sprintf("bmd_%05i.trj", i) u 1:($2+yshift):3 w l lw 0.5 lc rgb "#ffffff" not

Заключение

Бомовские траектории, хоть и физически неинтересны, позволяют рисовать красивые картинки, которые показывают, насколько квантовая механика веселее и богаче классической.

Если есть комментарии, вопросы, предложения: пишите. 🙂

Источник

Gnuplot, python, sciart, водород, квантовая механика, классическая механика