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

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

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

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

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

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

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

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

где m — масса частицы, x — координата, F — сила, действующая на частицу, а  — вторая производная координаты частицы по времени, или ускорение. Если в системе действуют только потенциальные силы, то силу можно выразить через новую сущность, потенциальную энергию V как 

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

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

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

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

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

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

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

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

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

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

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

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

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

Реализация

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

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.10f\n" % 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.10f\n" % tuple(r))
    
        t+=dt
    outf.close()

exit()

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

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

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

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

и для молекулы H2+

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

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

Заключение

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

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

© Habrahabr.ru