CC BY
175
АЛГОРИТМЫ И ВЫЧИСЛИТЕЛЬНЫЕ МЕТОДЫ
УДК 519.6+539.194
Г. В. Квитко, Э. Л. Кузин, Д. В. Шоть
ЧИСЛЕННОЕ РЕШЕНИЕ УРАВНЕНИЯ ШРЁДИНГЕРА С ПОЛИНОМИАЛЬНЫМИ ПОТЕНЦИАЛАМИ (Часть II)
Описана технология численной реализации фундаментального решения задачи Коши для одномерного уравнения Шрёдингера с полиномиальными потенциалами. Вариационным методом Ритца решается стационарное уравнение Шрёдингера для протона в потенциальной яме с двумя минимумами, разделенными барьером; на дискретном спектре этих решений формируется функция Грина, после чего численно рассчитывается фундаментальное решение нестационарного уравнения Шрёдингера. Конкретные расчеты применены к задачам с модельными адиабатическими потенциалами, параметры которых характерны для соединений с внутримолекулярными водородными связями.
The algorithm of the numerical realization of a fundamental solution of a Cauchy problem for the one-dimensional Schrodinger equation with polynomial potentials is described. The stationary Schrodinger equation for a proton in the double potential well is solved by using Ritz variation method. The Green's function is formed on a discrete spectrum of these solutions and after that the fundamental solution of the non-stationary Schrodinger equation is numerically calculated and applied to compounds with intramolecular hydrogen bonds.
Ключевые слова: уравнение Шрёдингера, протон, задача Коши, функция Грина, спектр, метод Ритца, адиабатический потенциал, численное решение.
Key words: Schrodinger equations, proton, Cauchy problem, Green's function, specter, Ritz method, adiabatic potential, fundamental solution.
1. Введение
В работе [1] была определена постановка задачи численного решения одномерного уравнения Шрёдингера с полиномиальными потенциалами, отмечены основные технологии для получения такого решения и приведен пример по расчету конкретно поставленной задачи. Речь, в частности, шла о предложенной авторами методике по численной реализации известного, но весьма формализованного фундаментального решения задачи Коши и о широко известном классическом разностном методе, основанном на неявной разностной схеме Кранка — Николсон.
В настоящей работе подробно описывается и исследуется предложная авторами методика численной реализации фундаментального решения у(x, t) задачи Коши для уравнения Шрёдингера:
Вестник Балтийского федерального университета им. И. Канта. 2011. Вып. 10. С. 79 — 93.
80
.н дц^) = _ 2_ 5 *) + V(х) у(х, *), у(х, 0) = /(х), (1)
где т = 6, т — масса квантовой частицы (протон); Й — постоянная Планка; і — мнимая единица.
Фундаментальное решение у(х, *) задачи (1) имеет вид:
+да
у(х, *) =| 0(х, х', *)/(х') йх'. (2)
_да
Оно справедливо для задач с произвольным потенциалом V(х) и для достаточно широкого класса функций /(х) [2]. Стоящая под интегралом в (2) функция Грина Є(х, х', *) для задач с дискретным спектром (именно они нас интересуют), имеет следующий вид:
х, ^ *) = |] ф„(х)ф„(х')ехр ^_Й-Еп *| (3)
Здесь фп (х), Еп — решения спектральной задачи
НФп (х) = Еп Фп (х) (4)
для оператора Гамильтона
h2 д2
H = —— —f + V (x) (5)
2m dx
с полиномиальным потенциалом
n k
V(x) = ^ vkx (n = 4; 6), (6)
k=0
моделирующим силовое поле в виде ямы с двумя минимумами, разделенными потенциальным барьером. В качестве начальной функции f(x) (для удобства) взят так называемый гауссов волновой пакет
/(х|ст, q, p) = (h стл) 1/4 exp (-(x - q)2 / 2 Йст + ip(x - q)/й) с полушириной
ст > 0 и центром (q, p), определяемым средними значениями операторов координаты x = x и импульса p = -й -3 / 3x.
2. Решение спектральной задачи
Из соотношений (3, 4) ясно, что для решения нестационарного уравнения Шрёдингера (1) необходимо предварительно решить задачу на собственные значения (4), то есть стационарное уравнение Шрёдингера. Такая задача для потенциалов вида (6) представляет собой самостоятельную и достаточно трудоемкую вычислительную проблему. Для ее решения в настоящей работе был использован известный метод приближенного решения стационарного уравнения Шрёдингера — вариационный метод Ритца [3], являющийся одним из вариантов так называемых проекционных методов [4]. Основные требования к его
применимости следующие: функции фп (х) должны быть определены в
гильбертовом пространстве, а оператор H должен быть симметричным и положительно определенным.
В соответствии с технологией вариационного метода Ритца, приближенное решение фп (х) задачи на собственные значения (4) ищется в виде разложения по некоторому конечному набору полной системы функций, образующих базис в рассматриваемом пространстве. В настоящей работе было использовано следующее представление:
N
Фп ю) = Е CnkUk (Х, юХ (7)
k=0
где uk (х, ю) — собственные функции некоторого гармонического осциллятора с собственной частотой (параметром) ю:
ик(х, ю) = (2k к!Х(ю)л/л)—1/2 exp[-(х / Х(ю))2 / 2] Hk (х /Х(ю)). (8)
Здесь Hк (^) — полиномы Эрмита, ^ = х / Х(ю), а Х(ю) = (ю m / h) 1/2.
Свойства базисных функций ик (х, ю) хорошо изучены: они орто-нормированы и удовлетворяют условию полноты. Коэффициенты разложения cnk в (7) находятся из решения системы линейных уравнений:
N
Е (нк— 8 jkE) • Cnk = 0
k=0 J J
+да
где 8 jk — символ Кронекера, Hjk = (HUj, uk) = J UjH uk йх, а сама эта
—да
система получается варьированием некоторого функционала E = (Hф, ф) / (ф, ф) (средняя энергия) по всем параметрам cnj. Условие
ее разрешимости следующее: det( Hjk — E 8 jk) = 0.
Совершенно очевидным становится факт, что для решения спектральной задачи (4) прежде всего требуется в базисе собственных функций гармонического осциллятора определить матричные элементы Hkj оператора Гамильтона H.
3. Расчет матричных элементов гамильтониана в базисе собственных функций гармонического осциллятора
Для эффективного отыскания матричных элементов оператора Гамильтона H в данной работе был применен известный операторный метод, разработанный Дираком [5] для гармонического осциллятора. Для удобства задача решалась в так называемом {QJ-представлении [6]. В этом случае оператор H Гамильтона (5) с полиномиальным потенциалом У(х) принимает следующий вид:
H = P2 / 2 + V(Q), (9)
81
где обобщенные операторы энергии H, координаты Q и импульса P связаны с прежними операторами H, X и p соотношениями: H = Й/Йю, Q = ^jmю /h • X и P = p/yfmoh. Выражение (9) легко преобразуется к следующему виду:
н = H +z (Xk-5k2/2)Qk, (l0)
k=0
где H0 = (P2 + Q2) / 2 — гамильтониан гармонического осциллятора. Здесь коэффициенты X k (k = 0, l,..., б) связаны с коэффициентами v ^
полинома V(x) соотношением: Xk = (h / mra)k/2(Йю)-) vk. Добавка к оператору H0 в (10) подчеркивает тот факт, что исходная задача о движении протона в полиномиальном потенциале может рассматриваться как ангармонический осциллятор.
Метод, предложенный Дираком, позволяет решать задачу на собственные значения в общем виде без ссылок на какое-либо представление, опираясь исключительно на постулаты пространства Гильберта и коммутационное соотношение Гейзенберга [Q, P] = i. Для этого вводят так называемые операторы рождения а + = (Q - i P) / V2 и уничтожения а = (Q + i P)/ V2 состояний, которые образуют алгебру Ли динамической группы Лоренца 0(2, 1), связывающей все уровни гармонического осциллятора в одно неприводимое представление [7]. Эти операторы эрмитово сопряжены друг с другом, удовлетворяют коммутационному соотношению [ й, й+ ] = 1 и действуют в пространстве, базис которого — функции un(x, ю) или векторы |un) = (nі)-1/2(а+ )n|u0) (в обозначениях Дирака), построенные с помощью операторов рождения а + из основного состояния |u„) = Х(ю)л-1/2 exp [-(x / Х(ю))2 /2].
Для решения поставленной задачи требовалось, во-первых, исследовать результат действия операторов Qm, входящих в потенциал V(Q),
и оператора P2 на собственные векторы гармонического осциллятора
и, во-вторых, вычислить матричные элементы гамильтониана (5) в базисе этих векторов. Для операторов Q и P имеют место следующие правила действия на собственные векторы | un) оператора H:
Q|un) = (MU«-,) + ^п + 11 un+l))/л/2,
P I Un) = -i • [Jn I Un-l)-Vn + 1 I Un+1)) / V2, (11)
\Un+l) =УІ2/(n + 1) Q|Un)-4n/(n + 1) Iun-1.
Последовательно и многократно действуя слева операторами Q и P на собственные векторы гармонического осциллятора в соответствии с правилами (11), получим все требуемые результаты типа Qm\u^) и
Pk | u^j. Затем, после несложных, но достаточно громоздких вычисле-
ний, удалось найти в явном вице матричные элементы оператора Гамильтона Н уравнения Шрёдингера (4). Для полинома степени т = 6 они имеют следующий вид:
Нп,п =Х0 +(Х2 +1/2)(п+1/2)+3/4-Х4(2и2 +2п+1) + 5/8-Х6(4и3 +6п2 +8п+3),
нп,п+г = Н+і,п = і + 3/2 X з(я +1) + 5/4X5(2 п2 + 4 п + 3)][(п +1)/2]1/2,
Н„,„+2 = Нп+2,п = [(X1 -1/2)+Х4(2п+3)+5/4X(п2 + 3п+3)] [(п+1)(п+1)/4]1/2,
4 5
Н п+з = Нп+з,„ = [X з + 5/2 X з(п + 2)][(п + 1)(п + 2)(п + 3)/8]1/2,
Нп+4 = Н„+4,„ = [X4 + 3/2-Xб(2п + 5)] [(п + 1)(п + 2)(п + 3)(п + 4)/16]1/2,
Н п+5 = Нп+5п = X 5 [(п + 1)(п + 2)(п + 3)(п + 4)(п + 5)/32 ]1/2,
Нп+6 = Н„+б,„ = X 6 [(п + 1)(п + 2)(п + 3)(п + 4)(п + 5)(п + 5)/64]1/2.
Все остальные матричные элементы Нк■ равны нулю. Видно, что матрица ||Нк)|| имеет ленточный вид и является симметрической: Нкі = Н^к. Если положить, что X5 =X6 = 0, то получаем матричные элементы для полинома степени т = 4.
Отметим важный момент в рассмотренной выше технологии. Поскольку мы имеем возможность определить любые результаты типа Qm |ип) и Рк\ип), а следовательно, и разные комбинации QmPk\un) или Рк0>т |ип), то описанная выше технология дает возможность находить
матричные элементы, а значит, и решать спектральные задачи для дифференциальных уравнений более широкого класса, чем уравнение Шрёдингера (4) с гамильтонианом (5).
После того как матрица ||Нк)||, представляющая оператор Гамильтона в базисе собственных функций гармонического осциллятора, была определена, стало возможным и решение собственно самой спектральной задачи.
4. Оптимизация решений по частоте базисных функций гармонического осциллятора
Для проведения численных расчетов в среде современных математических пакетов была разработана и отлажена вычислительная программа. В ходе вычислений варьировались параметры, задающие геометрию потенциала У(х), моделирующего потенциальную яму с двумя минимумами, разделенными барьером. Для задачи Коши (1) варьировались еще и параметры начальной функции / (х). Поскольку нас интересует решение прикладной задачи о движении протона в эффективном поле атомов в молекулах с водородными связями (прототропных таутомерах), то конкретные числовые значения основных параметров потенциалов, положения экстремумов {х1, х2, хз} и их значения (У1, У2, Уз} в расчетах привязывались к имеющимся экспериментальным данным, взятым из [8]. Координаты выражались в ангстремах, а характери-
стики потенциалов — в специфических единицах типа см-1, а также в электрон-вольтах (эВ) (1 см-1 -1,23984 -10-4 эВ ). Были рассмотрены как симметричные, так и несимметричные потенциалы.
Вариационный метод Ритца обеспечивает получение приемлемых решений ф(х) стационарного уравнения Шрёдингера (4) за счет соответствующего оптимального выбора коэффициентов спк в разложении (7) искомых решений по базису собственных функций ик (х, ю) гармонического осциллятора. Однако сами базисные функции зависят от параметра ю — частоты колебаний осциллятора. Понятно, что при ограниченном (обрезанном) базисе N величина этого параметра может существенно повлиять на решения ф(х) стационарного уравнения Шрё-дингера и, как следствие, в силу формулы (2) — на решения у(х, Ґ) нестационарного уравнения Шрёдингера (1). Поэтому было проведено исследование степени влияния параметра ю на получаемые приближенные решения. Для этого был предложен и реализован способ дополнительного уточнения численных решений задачи (4). Он достаточно прост. На варьируемом интервале частот [юо, ю^ ] с шагом Дю
для всех значений ю ^ (і = 0, 1, ..., й) методом Ритца решалась спектральная задача (4). По найденным решениям вычислялся функционал (критериальная функция) Е(ю), минимум которого соответствовал «наилучшему» решению ф(х) задачи.
Были рассмотрены два варианта задания Р(ю). В первом, наиболее простом в практической реализации, в качестве такого функционала использовалось минимальное собственное значение (нулевая энергия) Ео(ю). Во вто-
2 2 1/2
ром — среднеквадратичное отклонение ДН^ = [< Н >■ - < Н >2] га-
мильтониана Н, вычисленное для его низших собственных значений Є і так, чтобы выполнялось условие: (< Н >. -АН.) $ є. (< Н >. + АН-).
+да
Г *
Входящие в него интегралы вида < Н > = І фі (х, ю)Нфі (х, ю)йх опреде-
—да
ляют среднее значение энергии протона в і-м состоянии. Наибольшую трудность в практической реализации второго варианта представляла процедура отыскания матричных элементов для оператора Н2 (квадрат гамильтониана). Так, для потенциала У6(х) приходилось вычислять матрицы степеней операторов координаты Q вплоть до Q12.
Расчеты показали, что обе технологии обеспечивают практически совпадающие результаты. Но поскольку с вычислительной точки зрения и по временным затратам первая технология была существенно проще, то именно она и применялась в дальнейшем для проведения численных расчетов задачи.
Методика поиска Е0(ю) = тіп варьировалась от применения простого пошагового расчета до использования разных вариантов метода наискорейшего спуска. Расчеты показали, что в большинстве рассмат-
риваемых задач с потенциалами У(х) различной геометрии этот минимум присутствовал. Для потенциалов У4(х) такая «энергетическая яма» в основном имела вид плато (рис. 1). На нем значения энергии Е0 (ю і ) почти не менялись или слабо осциллировали в некотором интервале [ю', ю''] с [ю0, юа ] .
Для потенциалов У6(х), моделируемых полиномами 6-й степени, характер поведения в общих чертах был таким же, что и для У4(х). Однако плато уже не было постоянным, а имело вид некоторой монотонно меняющейся функции. Построенные на рисунке 1 кривые соответствуют потенциалам У4(х) и
о
У6(х) с геометрией {-0,3; 0,0; 0,3} А и {0; 2000; 0} см 1. Для задачи с потенциалом У6(х) были проведены сравнительные расчеты уравнения Шрёдингера для найденной по описанной выше процедуре оптимальной частоты ю0 = 3 • 1014 с 1 и частот, удаленных от ю0 (ниже в частотах ю і множитель 1014 п-1 опускаем).
Графики на рисунках 2 — 5 показывают, насколько важно при проведении расчетов стационарного уравнения Шрёдингера с ограниченным базисом функций осциллятора учитывать правильный выбор параметра м. Были проведены расчеты с частотами ю 1 = 0,5<ю0 и
ю2 = 20>ю0. Для трех функций ф0(х, юх), ф0(х, ю0) и ф0(х, ю2) получена картина, которая представлена на рисунке 2. Хорошо видно, что решение ф0(х, юх) имеет совсем иной вид — у него один максимум в центре, тогда как решения с другими частотами имеют два максимума, разделенных минимумом. Изменение значения м2 приводит к тому, что при сохранении формы с ростом м2 функция ф0(х, ю2) все сильнее расплывается по пространству координат, а ее амплитуда понижается. Что касается решения ф0(х, юх), то, как видно из увеличенного фрагмента на рисунке 5, а, оно имеет четыре корня, что неприемлемо для функции нулевого состояния, которая не должна иметь узлов.
Решения (2) задачи Коши (1) определяются через собственные функции стационарного уравнения Шрёдингера (4). Поэтому, как видно из рисунка 4, а, б, в стационарные решения в плохом базисе (с неудачно выбранной частотой осциллятора) сильно исказили основанные на них решения ^6(х, £ I м) нестационарного уравнения Шрёдин-
о
гера с потенциалом Уб(х): {-0,3; 0,0; 0,3} А и {0; 2000; 0} см-1. На рисунке 4, в, г визуализированы решения ^6(х, £ I м0 = 3) и ^4(х, £ I ш0 = 3) с одинаковой частотой ю0 = 3 • 1014 с-1. В этом случае различие в графиках
Рис. 1. Зависимость низших собственных значений (нулевых энергий) от частоты ш базисных осцилляторов
85
связано главным образом с различием в форме (крутизне) потенциалов. В силу этого фактора меняется частота колебаний самой функции у(х, Ґ). Видно, что за один и тот же промежуток времени в потенциале У4(х) функция совершила почти в 2 раза больше колебаний, чем в У6(х).
86
Рис. 2. Низшие волновые функции для У6(х) при разных частотах ш
. а и 1 х, А
- 0.4 - оУ о.о4 1 ^0.2 0.4
Рис. 3. Фрагмент графика функции нулевого состояния для потенциала У6(х), полученного для частоты ша < ш0
в
г
Рис. 4. Модуль волновой: а — функции Ф6(х, Ь1 Ш1 = 0,5); б — функции Ф6(х, Ь I ш2 = 20); в — функции Ф6(х, ^ I ш0 = 3); г —функции Ф6(х, t I ш0 = 3)
На рисунке 5 показано, как меняется ширина энергетического плато-интервала наименьших значений Е0(ю) для различных размеров базиса N. Расчеты проводились для потенциалов У4(х) и У6(х) с одинаковыми параметрами, зафиксированными выше. Полученные результаты подтверждают, вообще говоря, очевидный факт: при неограниченном увеличении размера базиса (N ) безразлично, какое значение ю ис-
пользовать в расчетах — все одинаково хороши. Однако при фиксиро-
ванном и существенно ограниченном базисе этот параметр начинает играть весьма существенную роль.
а б
Рис. 5. Зависимость нулевых энергий Е0(ш) (в эВ) от частоты ш (в пс-1) при различных размерах базиса N для потенциалов: а — У4(х), б — У6(х)
87
S. Результаты численных расчетов
При решении стационарной задачи в рамках вычислительной программы удалось получить (без введения дополнительных программ) собственные функции возбужденных состояний Фи(х, ю) вплоть до
N = 1SC. Основная проблема при получении решений с более высокими возбужденными уровнями обусловлена тем, что в выражение (S) для базисной функции ик (х, ю) входит нормировочный множитель
P( к) = (2 kk І Цю)7^ Г1/2 , который при росте к очень быстро убывает. Вычислительные возможности математических пакетов типа Maple, в которых была реализована расчетная программа, позволяют напрямую производить расчеты примерно для указанного выше ограничения по размеру базиса. С другой стороны, вести борьбу за расширение базиса в рамках вариационного метода Ритца вряд ли целесообразно. Хотя этот метод и является весьма эффективным для приближенного решения стационарного уравнения Шрёдингера, однако хорошо известны проблемы, связанные с его применением [9]. Расчет возбужденных состояний вариационным методом сложен и неточен. Нижний корень представляет собой значительно более хорошее приближение для оценки основного состояния, чем большие по величине корни для возбужденных состояний. Трудно однозначно ответить на вопрос о том, позволит ли последовательное расширение базиса достаточно точно аппроксимировать собственные значения и собственные функции. Кроме того, вариационный принцип практически не позволяет оценить погрешности результата. Вместе с тем применение метода Ритца для решения спектральной задачи уравнения Шрёдингера (4) вполне оправданно и целесообразно при решении огромного количества практических задач квантовой физики, химии и биологии, особенно тех, где высоковозбужденные состояния не играют существенной роли. Расчеты стационарного уравнения подтвердили известный в физике
ss
факт. В потенциалах рассматриваемого типа — в вице ямы с двумя минимумами, разделенными барьером, при определенной геометрии, когда потенциальный барьер достаточно высок, — возможно образование (под потенциальным барьером) так называемых расщепленных энергетических уровней. Иначе говоря, среди решений возникают пары очень близких, иногда почти вырожденных состояний, когда разница между соседними нижними уровнями энергии ничтожно мала.
На рисунке 6 прецставлен график, наглядно иллюстрирующий этот факт, с результатами численного расчета стационарного уравнения Шрёдингера для потенциала V4(x) с параметрами: ю0 = 5 • 1014 с 1 и Утах =5000 см 1 * 0,62 эВ, {-0,5; 0; 0,5} А
и {0; ^тах;°}.
С учетом соотношений (3) и (7) формальное выражение (2) для решения у(х, £) задачи Коши (1) можно записать в следующем, более удобном для проведения численных расчетов, виде:
(12)
Рис. 6. Энергетический спектр (в эВ) для потенциала P4(x): є0 = 0,13898, є2 = 0,13899, є2 = 0,39459, є3 = 0,39671, є4 = 0,58703, є5 = 0,62645
ЦІ(x, t) = Хфп (x) ЄХР (—i En t / ЙУп (N),
n=0
N +да N
где In(N) = SCnk J uk(|, ra)/(I)dI = SCnkJ(k).
k=0 _да k=0
+да
Вычисление несобственных интегралов типа J(k) = J Uk(I, ra)/(I)dI
—да
сопряжено с определенными вычислительными трудностями, связанными со специфическим видом базисных функций. Созданная вычислительная программа позволяла напрямую производить расчеты интегралов ](k) примерно дня значений k 70. При больших номерах k > 70 вычислительные системы не справляются с задачей и приводят к аво-стам. Было проведено исследование по вычислению определенных ин-
+L
тегралов вида J(k, L) = J Uk(I, ra)/(I)dI при L , которое показало, что
—L
встроенные процедуры интегрированных математических пакетов достаточно хорошо справляются с расчетом несобственных интегралов. Сравнение велось с результатами, полученными на основе специально разработанных программ интегрирования в системе VisUal Fortran. Численные исследования показали, что решения задачи Коши, получаемые при работе с базисом N $ 50, практически не отличались от решений для максимально реализуемого в расчетной программе базиса N = 70. Поэтому все дальнейшие расчеты проводились для N = 50.
При формировании решений y(x, t) нестационарного уравнения Шрёдингера по формуле (12) был исследован на сходимость ряд
^ Фи(х)ехр{-ЇЕпї / Й)4(И). Было показано, что он довольно быстро
п=0
сходится. Приемлемые решения получаются уже при М = 20. Последнее значение и было использовано авторами при проведении расчетов.
Численные расчеты показывают, что внесение слабой асимметрии в форму потенциала У4(х) оказывает влияние в основном лишь на нижние состояния протона. Из рисунка 7 (параметры те же, что и на рис. 6) хорошо видно, что для полностью симметричного потенциала {0; 5000; 0} см-1 низшее энергетическое состояние (рис. 7, а) симметрично и имеет два ярко выраженных максимума, локализованных в районе левой и правой ям потенциала каждый, в то время как первое возбужденное состояние (рис. 7, в) антисимметрично.
89
в
Рис. 7. Низшие (а), (б) и первые возбужденные состояния (в), (г) для симметричного (слева) и асимметричного (справа) потенциалов Р4(х)
Такая картина отвечает свойствам решений стационарного уравнения Шрёдингера. Но для асимметричного потенциала вица {0; 5000; 50} ситуация уже иная. Для низшего состояния (рис. 7, б) волновая функция ф[0“](х) практически полностью локализована в области левой (более глубокой) ямы потенциала с максимумом в центре ямы. Первое возбужденное состояние (рис. 7, г), наоборот, локализовано в правой яме, и максимум функции приходится на правый минимум потенциала.
Анализ показывает (рис. 8), что это состояние, как и требуется, имеет один узел (при х = -0,4884 А). Спектр энергий практически одинаков для обоих случаев. Решения у(х, Ь) уравнения Шрёдингера для обоих вариантов потенциала также практически одинаковы.
б
а
г
90
Один из вариантов для У4(х) представлен на рисунке 9. Причем
рисунок 9, б — это контур графика
модуля функции у(х, Ї), то есть его
вид сверху. Хорошо видна пространственно-временная симметрия функции | у(х, £)|. Графики представлены в узлах следующей сетки:
^Ьт = ^к Х ^1, где =
= {хэ = эк, б = 0, 1, ..., Их, ЬИх = Ь, Ь = хИх - х0 х0 = -1,0А, хИх = 1,0А, Их = 100},
йх = {ік = кі, к = 0, 1, ., Ыи = Т, Т = N - ^0, ї0 = 0, N = 5-10—13 с, N = 200}.
Рис. 8. Фрагмент рис. 6, г для ф^х), показывающий наличие узла
\ % у М г
а б
Рис. 9. Модуль решения у(х, {) уравнения Шрёдингера (1) для потенциала У4(х) с параметрами N = 100, У1 = 0, У2 = 5000, У3 = 0, х1 = -0,5, х2 = 0, х3 = 0,5, ш0 = 5, ц = 0, р = 0,1, ст = 0,005: а — общий вид; б — вид сверху
На рисунках 10 — 14 представлены результаты численных расчетов для потенциалов У4(х) и У6(х) с небольшой асимметрией при следующих расчетных параметрах: базис N = 100; значения потенциалов У1 = 0, У2 = 5000, У3 = 50 (в см-1) заданы в точках экстремумов хг = -0,5, х2 = 0 и х3 = 0,5 А; частота базисного осциллятора Ш0 = 5 • 1014 с-1; параметры начальной функции (волнового пакета) ц = 0, р = 0,1, ст = 0,005.
На рисунке 10 на фоне разных одномерных потенциалов показаны нижние уровни энергетического спектра. є0=0,21179,
Рис. 10. Потенциалы У4(х) и У6(х) и энергии в эВ
є2=0,51248,
є4=0,70052,
для Уб(х):
Є!=0,02408,
є3=0,5524,
є5=0,83030.
Для У4(х): Є!=0,21834, є3=0,53421, є5=0,83686; є0=0,01055, є2=0,53642, є4=0,69229,
Ф04]( х) 2 А
—V- -1-І 0 °х(А)
а б
Рис. 11. Низшее (а) и первое возбужденное (б) состояния в потенциале У4(х)
Фі(х)
На рисунках 11—12 даны нулевое и первое возбужденное стационарные состояния. В целом общий характер поведения решений в этих потенциалах схожий, но
численно их характеристики заметно различаются. Первые два низших энергетических уровня в потенциале У6(х) меньше, чем в У4(х), на порядок, тогда как два следующих уровня почти не отличаются. Амплитуды состояний, как и следовало ожидать, отвечают минимумам потенциалов, но в У6(х) они достаточно высокие в обеих ямах, тогда как в У4(х) их различие более заметно.
На рисунках 13 — 14 в разных ракурсах представлены графики модуля |ф(х, Ґ) I решения задачи Коши. Как уже было отмечено ранее (см. рис. 4, в, г), заметное различие между приведенными графиками связано с различием в форме (крутизне) сравниваемых потенциалов. Частота колебаний функции ф(х, Ґ) в потенциале У6(х) оказалась существенно выше, чем в У4(х). Поэтому за один и тот же расчетный интервал времени функция ф(х, Ґ) в У6(х) совершает гораздо большее число колебаний, проходит полный цикл и фактически возвращается к своему первоначальному виду.
х(А)
б
Рис. 12. Низшее (а) и первое возбужденное (б) состояние в потенциале У6(х)
91
Рис. 13. Модули решения ф(х, І) уравнения Шрёдингера (1) для потенциала У4(х)
Рис. 14. Модули решения y(x, t) уравнения Шрёдингера (1) для потенциала V6(x)
Эффективность работы вычислительной программы с потенциалами в виде полиномов разной степени предварительно была протестирована на решении более простых задач с потенциалами, для которых имелись известные аналитические решения задачи Коши (1): свободное движение, движение в линейном и квадратичном поле и в потенциале Морзе. Проверка показала хорошее согласие аналитических решений с решениями, получаемыми путем численной реализации фундаментального решения (2).
Как было отмечено в работе [1], помимо описанной выше технологии для отыскания фундаментального решения y(x, t) задачи Коши (1) авторами использовалась широко применяемая на практике классическая двухслойная неявная схема Кранка — Николсон [10]. Проведенные по этой схеме сравнительные численные расчеты ряда задач с потенциалами разной геометрии показали хорошее согласие с результатами, численно реализующими фундаментальное решение (2).
В перспективе планируется в рамках квантово-статистической модели [11; 12] использовать предложенные выше методики численного решения уравнения Шрёдингера для исследования процессов, обусловливающих внутримолекулярные таутомерные превращения.
Работа выполнена при финансовой поддержке РФФИ по проекту № GS-G1-GG-431.
Список литературы
1. Квитко Г. В., Кузин Э. Л., Шоть Д. В. Численное решение уравнения Шре-дингера с полиномиальными потенциалами (Ч. 1) // Вестник Балтийского федерального университета им. И. Канта. Калининград, 2011. Вып. 5. С. 115 — 119.
2. Polyanin A. D. Phandbook of linear partial differential equations for engineers and scientists, Chapman & Hall/CRC, 2002. URL: http://eqworld.ipmnet.ru/en/ so-lutionslpde/ lpde108.pdf.
3. Михлин С. Г. Вариационные методы в математической физике. М., 1970.
4. Марчук Г. И., Агошков В. И. Введение в проекционно-сеточные методы. М., 1981.
5. Дирак П. А. М. Принципы квантовой механики. М., 1960.
6. Мессиа А. Квантовая механика. Т. 1. М., 1978.
7. Малкин И. И., Манько В. И. Динамические симметрии и когерентные состояния квантовых систем. М., 1979.
8. Brickmann J., Zimmermann H. Lingerig time of proton in well of double-minimum potential of hidrogen bonds // The Journal of Chemical Physics. 1966. Vol. 50, N 4. P. 1608-1618.
9. Цюлике Л. Квантовая химия. Т. 1. М., 1976.
10. Самарский А. А. Теория разностных схем. М., 1977.
11. Квитко Г. В., Кузин Э. Л., Новиков В. И. Квантовая статистическая модель внутримолекулярного таутомерного превращения // Теоретическая и экспериментальная химия. 1975. Т. 11, № 6. С. 754 — 761.
12. Квитко Г. В., Кузин Э. Л., Шоть Д. В. Математическая модель внутримолекулярного таутомерного превращения и процессы релаксации протона // Вестник Российского государственного университета им. И. Канта. 2009. Вып. 10. С. 104 — 111.
Об авторах
Геннадий Васильевич Квитко — канд. физ.-мат. наук, доц., Балтийский федеральный университет им. И. Канта, e-mail: gkvitko.univ@ gmail.com.
Эдуард Леонидович Кузин — канд. хим. наук, доц., Балтийский федеральный университет им. И. Канта, e-mail: eduard_kuzin@mail.ru.
Дмитрий Владимирович Шоть — асп., Балтийский федеральный университет им. И. Канта, e-mail: d. schott@triumph-adler.ru.
Authors
Dr Gennady Kvitko — assistant prof., I. Kant Baltic Federal University, e-mail: gkvitko.univ@gmail.com
Dr Eduard Kuzin — assistant prof., I. Kant Baltic Federal University, e-mail: eduard_kuzin@mail.ru
Dmitry SAott — PhD student, I. Kant Baltic Federal University, e-mail: d. schott@triumph-adler.ru
