УДК 519.633, 538.9
Численное моделирование эволюции состояний
полярона
И. В. Амирханов*, Е. В. Земляная*, В. Д. Лахно^, Д. З. Музафаров*, И. В. Пузынин*, Т. П. Пузынина*, З. А. Шарипов*
* Лаборатория информационных технологий Объединённый институт ядерных исследований ул. Жолио-Кюри, д.6, Дубна, Московская область, 141980, Россия ^ Институт математических проблем биологии РАН ул. Институтская, д. 4, г. Пущино, Московская область, 142290, Россия
В работе исследуется эволюция полярона в однородной среде в зависимости от параметров модели и начальных условий, которые выбираются в виде различных комбинаций стационарных поляронных состояний. Представлены вычислительная схема и результаты численного моделирования.
Ключевые слова: эволюция полярона, конечно-разностная схема, численные методы.
1. Введение
Поляроны определяют многие процессы в ионных кристаллах [1], полупроводниках, полярных жидкостях и биологических системах [2-4]. Представление о биполяронных состояниях играет важную роль при объяснении такого явления, как высокотемпературная сверхпроводимость. Поляронные состояния используются в современной наноэлектронике при описании переходов в квантовых точках. Поляронными эффектами объясняются полосы поглощения центров окраски в ионных кристаллах. В полярных средах сольватированные состояния электронов представляют собой поляронные состояния и определяют химические реакции, выполняя роль сильнейшего восстановителя. В полимерах поляроны являются основными носителями тока. Их проводящие свойства используются при создании сверхлёгких проводников и аккумуляторов. В биологии поляроны или солитоны объясняют возможность переноса энергии на большое расстояние. Их изучение даёт основу для создания таких качественно новых устройств нанобио-электроники, как нанобиочипы и электронные нанобиосенсоры.
Модель эволюции произвольного начального состояния полярона описывается системой связанных квантово-классических динамических уравнений [5,6]. Это система нелинейных интегро-дифференциальных уравнений, общими характеристиками которой являются многопараметричность (физическими параметрами задачи являются: V — скорость полярона, П — частота оптических колебаний ионов, 7 — параметр трения, т* — масса полярона и т.д.) и многомерность конфигурационного пространства. Стационарные решения этой системы исследованы многими авторами (см., в частности, [7] и цитируемую литературу).
Основной задачей нашего исследования является изучение временной эволюции различных начальных состояний полярона в отсутствии и при наличии трения. Здесь мы ограничиваемся случаем неподвижного полярона (V = 0). В работе исследуется численная схема для решения системы нелинейных дифференциальных уравнений, описывающей эволюцию полярона в однородной среде, и представлены результаты численного моделирования для конкретного набора значений физических параметров модели.
Работа выполнена при финансовой поддержке РФФИ (гранты 07-07-00313, 08-01-00800, 07-01-00738, 09-01-00770).
2. Постановка задачи
В работе [5] предложена система нелинейных уравнений для описания эволюции полярона. В частном случае, для сферически симметричного неподвижного полярона с учётом трения [6] эта система записывается в следующем виде:
•г. - д д2
г 2т— + + 2т— at дх2 х
Ф = 0,
а2 у
дх2 д2
= в,
(1)
д
Ж2 +7a¡+ и
в=
и
£ X
где ф — волновая функция, ^ — потенциал, в, т, 7, ш, £ — безразмерные параметры модели. Система (1) дополняется следующими начальными и граничными условиями:
Ф(х, í)Uo = Ф*(cos(Afc^/4) + г sin(Afc^/4)),
lí=0 1Ф| а
в(х, í)|t=o = ~ -f, ^ в(х, t)
е х
= 0, у(0) = 0, У (то) = 0.
(2)
í=0
Здесь А& и Ф& — собственные значения и собственные функции соответствующей стационарной задачи:
d2 0 Ф(х)
—- 2тА + 2т——
d х2 х
d2 ^ 1ф2(х) -т-2 Ф(х) = — -——,
dх2 £ х
Ф(х) = 0,
0 < х < то,
(3)
с граничными условиями и с условием нормировки:
сю
Ф(0) = 0, Ф(0) = 0, Ф(то) = 0, Ф'(то) = 0, J Ф2(х)dх = 1. (4)
Решая систему (3),(4) c помощью непрерывного аналога метода Ньютона [7], находим решения {Ф&, А&}, где к = 0,1, 2,..., {Ф0 ,А0} — собственная функция и собственное значение основного состояния (безузлового решения), {Фх,Ах} — собственная функция и собственное значение первого возбуждённого состояния и т.д. На рис. 1 показаны первые три собственные функции системы (3), (4). Соответствующие собственные значения равны А0 = —0,16277, Ai = —0, 0308, А2 = —0, 0125 (т = 1, е = 1).
2
3. Схема численного исследования
Введём равномерную сетку о шагами Нш, Н соответственно по переменным х
и ¿:
{хт = шНж(ш = 0,1,..., ¿), ¿п = = 0,1,...)} .
Для решения системы (1) с начальными и граничными условиями (2) будем использовать следующую неявную конечно-разностную схему порядка аппроксимации 0(Н4 + Н2) [8]:
> /
Рис. 1. Собственные функции системы (3)
фп+1 _ ф гт г:
а
Г„++\ - 2^ + С+-1
2т к1
тп+1 + 1
ткх т
+
+ (1 - а)
1Фт+1 - + ^т-1 + ^
2т №
тк
п
гт
- 1 + ^ = ©п+1 м ©т '
2© п + © п_ 1 ^ т. ^т. * ^т.
(5)
к2
+ 7
9эт+1 ©то + ш2ф"+1 — — 1
к.
е тк х
■Фт — Ф к (сон(\к п/4) + г вт(Ак ^/4)); 0т1 — --
1 К
I 2
9и — ©
£ ткх
^ — 0; р? — рР-1; Ш — 1, 2,...,/; п — 0,1, 2,...
где а — 0, 5, Ф к, Ак соответственно собственные функции и собственные значения стационарной задачи (3), (4).
Для решения задачи (1), (2) по схеме (5) на каждом слое с номером п использовался следующий алгоритм:
1) решается третье уравнение при известном фп относительно ©"+1;
2) решается второе уравнение для найденного ©"+1, определяется п+1;
3) решается первое уравнение и вычисляется гфп+1 на следующем временном слое;
4) переход кп.1 для следующего значения п.
Тестирование вычислительной схемы (5) проводилось с помощью модельных расчётов для уравнения Шрёдингера с кулоновским потенциалом, которое совпадает с первым уравнением системы (1) при (р — 1. В этом случае для уравнения Шрёдингера известны аналитические решения. По результатам проведённого сравнительного анализа численных и аналитических решений уравнения Шрё-дингера были выбраны параметры дискретной сетки к х — 0, 01 и к г — 0, 001, для которых отклонение численного решения от аналитического не превышает 4-10-4 на промежутке времени 0 ^ I < 105, соответствующем характерному периоду колебаний в системе «электрон-деформированное поле».
Для визуализации численных результатов вычислялась энергия Ш (¿) по формуле:
дф(х, Ь) дх
ёж —
ёж.
к
х
Отметим, что поскольку расчёты велись в безразмерных единицах, энергия Ш(¿) также является безразмерной величиной.
4. Численные результаты и выводы
Численные эксперименты [9] показали, что если в качестве начального условия (2) взято стационарное состояние полярона, полученное путём численного решения задачи (3), (4), форма полярона со временем не меняется.
В данной работе исследовалась эволюция начальных поляронных состояний, заданных в форме комбинаций стационарных состояний полярона:
Ф(®, i)|i=o = W Ф(®, i)|t=0 Ф(®, i)|t=0
Ar
4
Ar
Ф0 exp ъъ— + Ф1 exp ¿ъ—
Ai
4
Фо exp гъ-т- + Ф2 exp яъ—-
4
A-
4
Ф1 exp ( j + Ф2 exp i
(6)
(7)
(8)
Здесь N — нормировочная константа, Фо — волновая функция основного состояния, Фх,2 — волновая функция первого и второго возбуждённых состояний.
На рис. 2, 3, 4 представлены результаты численного решения задачи (1), (2) о начальными условиями (6), (7), (8) соответственно, при значениях параметров т = 1, ш = 1, е = 1, N = 0, 5, 7 = 0 и 4.
0,20
0,15
0,10
0,05
0,00
......t = 0
----1 = 200
-1 = 1000
10
20 30
a)
40
50
1000 2000 3000 4000
б)
0,18 0,16 0,14 0,12 0,10 0,08 0,06 0,04 0,02 0,00
Ivl
t = 1000 t = 5000 t = 10000
-0,06
-0,08
-0,10
-0,12-
-0,14
W(t)
2000 4000
6000 8000
Рис. 2. Эволюция полярона из состояния (6) и соответствующая энергия электрона
Ш(£) при 7 = 4 (а, б) и 7 = 0 (в, г)
0,20
Ы
......1 = 0 -0,04
----1 = 1000
-1 = 6000 -0,06
-0,08
-0,10
-0,12
-0,14
10 20 30 40 50 60 70 80
а)
0,07 0,06 0,05 0,04 0,03 0,02 0,01 0,00
М
А
4 = 1000 4 = 6000 I = 10000
6000 8000
б)
6000 8000
Рис. 3. Эволюция полярона из состояния (7) и соответствующая энергия электрона
Ш(г) при 7 = 4 (а, б) и 7 = 0 (в, г)
Ы
......1 = 0 ----1 = 500 -0,02-
-1 = 9000 -0,04-0,06-
-0,08-
-0,10-
-0,12-
\ -0,14-
60 ) 80
-0,02 -0,03 -0,04 -0,05 -0,06 -0,07
Щ*)
6000 8000
б)
6000 8000
Рис. 4. Эволюция полярона из состояния (8) и соответствующая энергия электрона
Ш(г) при 7 = 4 (а, б) и 7 = 0 (в, г)
На основании проведённого численного моделирования можно заключить, что начальные распределения заряда, заданные суперпозициями (6), (7), (8), при наличии в системе затухания (7 = 0), с течением времени эволюционируют в основное состояние. При отсутствии в системе затухания в промежутке времени 0 ^ £ < 105 эволюция в основное состояние не наблюдается.
Время эволюции полярона в основное состояние при 7 = 0 зависит от типа комбинации начального состояния. Так, для начального условия в форме (6) время эволюции в основное состояние составляет £ ~ 3000, для (7) и (8) — соответственно Ь ~ 6000 и Ь ~ 9000.
Литература
1. Пекар С. И. Исследования по электронной теории кристаллов. — М.: Гостех-издат, 1951.
2. Давыдов А. С. Солитоны в молекулярных системах. — Киев: Наукова Думка, 1988.
3. Polarons and Applications / Ed. by V. D. Lakhno. — Wiley, Chichester, 1994.
4. Компьютеры и суперкомпьютеры в биологии / под ред. В. Д. Лахно, М. Н. Устинин. — Москва-Ижевск: Институт компьютерных исследований, 2002.
5. Давыдов А. С., Энольский В. З. Трёхмерный солитон в ионном кристалле // ЖЭТФ. — 1981. — Т. 81, № 3(9). — С. 1088-1098.
6. Lakhno V. D. Dynamical Polaron Theory of the Hydrated Electron // Chemical Physics Letters. — 2007. — Vol. 437. — Pp. 198-202.
7. Пузынин И. В. и др. Обобщённый непрерывный аналог метода Ньютона для численного исследования некоторых нелинейных квантово-полевых моделей // ЭЧАЯ. — 1999. — Т. 30, № 1, вып. 21. — С. 210-262.
8. Самарский А. А. Теория разностных схем. — М.: Наука, 1989. — С. 296-299.
9. Амирханов И. В. и др. Численное исследование динамики поляронных состояний // Вестник ТвГУ, серия «Прикладная математика». — 2009. — № 2[13]. — С. 5-14.
UDC 519.633, 538.9
Modeling of the Evolution of the Polaron States
I. V. Amirkhanov*, E. V. Zemlyanaya*, V. D. Lakhno1", D. Z. Muzafarov*, I. V. Puzynin*, T. P. Puzynina*, Z. A. Sharipov*
* Laboratory of Information Technologies Joint Institute for Nuclear Research Joliot-Curie 6, 141980 Dubna, Moscow region, Russia t Institute of mathematival problems of biology Russian Academy of Sciences Institutskaja str., 4, 142290, Pushchino, Moscow Region, Russia
The evolution of polaron in a homogeneous environment is analyzed depending on parameters of the model and initial conditions which are selected in the form of various combinations of stationary polaron states. A computational scheme and results of the numerical modelling are presented. Work supported by RFBR grants 07-07-00313, 08-01-00800, 07-01-00738, 0901-00770.
Key words and phrases: evolution of polaron, finite-difference scheme, numerical methods.