Метод функций Грина для расчета нестационарных спектров люминесценции неравновесных молекулярных систем

DOI:

Интерес к сверхбыстрым фотохимическим процессам обусловлен их ролью в живой природе и возможностями использования в устройствах преобразования солнечной энергии. В настоящее время значительная часть экспериментальных исследований в этой области проводится с помощью спектроскопических методов и основана на регистрации динамического отклика системы на возбуждение коротким лазерным импульсом. Данный метод является базовым для фотохимии и применяется для исследования как небольших неорганических молекул в жидкостях или смесях, так и сложных биологических объектов, например, фотосинтетических реакционных центров бактерий и растений [1; 5; 9–12; 14; 15; 17; 18].

Время-разрешенные спектры поглощения и люминесценции содержат информацию о заселении электронных и колебательных состояний системы в ходе фотореакции, что дает возможность исследования ее механизмов. Вместе с тем вопрос о корректной интерпретации экспериментальных данных в случае сверхбыстрых реакций на пикосекундных временах остается актуальным. Традиционный подход, основанный на разложении нестационарных спектров на динамические компоненты, каждая из которых «привязывается» к некоторому электронному состоянию (state-associated analysis), оказывается неточным. В условиях неравновесности спектральная динамика системы зависит не только от химических превращений (то есть изменений населенностей состояний), но и от эволюции спектрального отклика каждого из регистрируемых состояний. В частности, релаксация высокочастотных внутримолекулярных колебаний в сверхбыстрых реакциях вызывает не только сдвиг спектра люминесценции (стоксов сдвиг), но и значительные изменения самого спектрального профиля [2]. Предположение о неизменности формы отдельных компонент в спектре в этом случае нельзя считать оправданным, а корректный анализ спектральной динамики требует учета неравновесности внутримолекулярных степеней свободы донорно-акцепторной системы и ее окружения. Разработка теоретических моделей неравновесных процессов (совместно с инструментами численного моделирования спектральной динамики и методов фитирования результатов эксперимента) может быть полезна с точки зрения повышения точности анализа экспериментальных данных.

Одним из ключевых методов экспериментального исследования релаксационных процессов в жидкостях является оптическое зондирование фотовозбужденных молекулярных соединений с переносом заряда. Эта группа соединений характеризуется широкой полосой поглощения, связанной с переносом электрона между функциональными центрами. Возбуждение полосы с переносом заряда формирует неравновесное франк-кондоновское состояние, которое можно рассматривать как электрический диполь в неравновесном окружении. Прямая связь между релаксацией среды и спектральной динамикой позволяет исследовать релаксационные процессы путем регистрации флуоресцентной динамики или с помощью пробных оптических импульсов. В настоящее время данный метод считается одним из самых надежных для изучения релаксации недебаевских растворителей, однако требует использования сложной техники для генерации сверхкоротких лазерных импульсов и специализированных методов обработки экспериментальных данных [9].

Данное исследование посвящено разработке математической модели спектральной динамики макромолекулярных систем, в которых фотовозбуждение запускает серию сверхбыстрых реакций электронного переноса в участием нескольких окислитель- но-восстановительных центров. Основное внимание в работе уделено учету неравновесных состояний среды и внутримолекулярных степеней свободы, формируемых как на этапе фотовозбуждения, так и в ходе многостадийной реакции. Учет релаксационных процессов внутри макромолекулы и в ее окружении позволяет существенно расширить круг исследуемых явлений, в частности, использовать модель для интерпретации экспериментальных данных о многостадийном переносе заряда в биологических объектах. Математической основой разрабатываемой модели является полуклассическая теория многостадийного переноса электрона в многокомпонентном растворителе [4; 8] и метод функций Грина для классических (поляризационных) и квантовых (внутримолекулярных) координат системы.

• 1.    Метод

Рассмотрим механизмы динамического влияния окружения макромолекулы на спектральный профиль люминесценции на временах, сравнимых с временем релаксации. Инициированное электронным туннелированием перераспределение плотности заряда внутри макромолекулы сопровождается реорганизацией среды и возбуждением внутримолекулярных колебаний. В силу более слабого взаимодействия окружения с переносимым электроном, моды растворителя, как правило, низкочастотны (Q ^ к в Т ), поэтому при моделировании традиционно рассматриваются как классические степени свободы. Внутримолекулярная реорганизация, наоборот, чаще связана с высокочастными колебаниями (Q >  к _в Т ), поэтому внутримолекулярные координаты описываются квантово-механически. Количество задействованных в фотореакции колебательных состояний может достигать 10 ⁴ –10 ⁶ при суммарной энергии реорганизации Л _у = 0,2 - 0, 6 эВ [7].

С учетом большого числа активных внутримолекулярных и межмолекулярных мод, полуклассические модели электронного переноса традиционно используют «коллективные» координаты реакции, в частности, энергоразностную координату AB (energy gap coordinate), вычисляемую путем проецирования множества ядерных степеней свободы на энергетический зазор между состояниями реагентов и продуктов [13; 19]. Для многостадийных реакций более характерно использование целого набора независимых поляризационных растворителей [4; 8]. В дальнейшем в этой работе мы используем математический аппарат общей теории многостадийных реакций в недебаевской среде, предполагающий построение многомерного конфигурационного пространства.

Рассмотрим макромолекулярное соединение со структурой, включающей I различных окислительно-восстановительных центров. Состояние с локализацией электрона на г-м центре обозначим как |ф « ) (г = 1, ••• , I), при этом будем предполагать, что в общем случае электронные переходы возможны между любыми состояниями |ф « ) и ф). Введем энергию реорганизации Х^ , характеризующую отклик среды на смещение электронной плотности |ф г ) ^ |ф ^ ). Поскольку Хц = 0 и X _tj = Х^ , значения X _tj образуют симметричную квадратную матрицу с нулевой главной диагональю. Введем набор независимых координат среды Q _k , где к пробегает значения от 1 до К = I — 1. В приближении линейного взаимодействия электрона со средой диабатическая поверхность свободной энергии (ПСЭ) электронного состояния |ф » ) в координатах Q _k может быть записана в виде [4; 8]:

К

'■' ' ( Q ) = E (q₆

k=1

-

₂          a

Q k,i^ + ^G i +     K^ a Tl Q. .

a =1

Здесь Q _k,_i — координаты минимума ПСЭ в К -мерном пространстве; G _i — равновесное значение свободной энергии в г-м состоянии; ^ _а — частота a-й высокочастотной квантовой моды (а = 1, • • • , Л). Колебательное состояние системы описывается набором квантовых чисел n = (71 1 ,71 2 , ... ,71 a ) , задающих уровни возбуждения внутримолекулярных мод. Отметим, что при построении системы (1) значения G _i , Л и Q _a считаются известными параметрами модели, а величины Q _k,_i требуют определения на основе матрицы λ. Общие алгоритмы расчета координат минимума ПСЭ описаны в работах [4; 8].

Энергия реорганизации электронного перехода |фi) ^ ф) может быть определена соотношением Aij = Gj0)(Qi) — G(0)(CQi), тогда для энергетического зазора между состояниями |фi, n) и ф, m) в любой точке пространства справедлива формула кa

^G^j ’ ) = \j + AGij — 2 ^E QQk — Qk,i^ (Qk,j — Qk,i^ +    Tl9,a (ma — На) .

k=1a

Здесь AG _ij = G j — G _i — свободная энергия электронного перехода. Перепишем выражение (2) в более компактном виде, используя вектора

Q i = (Q1

^{— Q} 1 ,i , ^Q 2 ^{— Q} 2 ,i , . . . , ^Q K

^— Q K,i?) ,

D ij = ( Q 1,j ^— Q 1,i , ^Q 2,j ^{— Q} 2,i , .

, ^Q K,j

^{— Q}K,/) ,

определяющие положение точек Q и Q j в локальной системе координат, связанной с точкой Q _i . С учетом этих обозначений

А

^AG^ j , ) = A ij + ^{A G} ij ^{— 2 Q} i • D ij + ^^ h^ a ^{( m} a ^— ^-a ) ,               ⁽³⁾

а =1

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

Для скорости излучательных переходов вида |ф _i , n ) ^ |ф ^- , m ) используем выражение, полученное в рамках теории возмущений по величине оптической связи между начальным и конечным состояниями

К П ^т>( _Ш ) = 2П | v ;j" ^m) | ²,_: A-: '-m' — Ц.                   (4)

В отсутствие механизмов уширения излучательные переходы между этими состояниями характеризуются спектральной плотностью S _ij (ш), которая может быть рассчитана из уравнения [2]

S ij (ш) = Л о E | ^ j n ^{, m)} | ² / ^d Q Q i n^} ( Q ,t) 5 ( AG ; n ^{, m)} — Ь ш ) ,            (5)

n , m

B (n , m)

этом выражении vij    — энергия связи; Ло — нормировочная константа; y(n) — плотность частиц в состоянии |фi, n).

Вычислим далее y ⁽ⁿ ) ( Q ,t) с учетом релаксационных процессов в молекулярной системе, используя технику функций Грина. Эволюция функции плотности в состоянии с локализацией электрона на г-м молекулярном центре задана уравнением

где C — лиувиллиан системы. В отсутствие безызлучательных переходов C может быть представлен в виде суммы

ка где £^к и С^ соответствуют отдельным классическим (низкочастотным) и квантовым (высокочастотным) координатам. Предполагая независимость координат Qk и па, представим функцию Грина для оператора Лиувилля Д в виде произведения

R(Q, n, t \ Qo, no) = П F(Qk, t\Qko) П P(na, tW,(8)

кa где F(Qk, t\Qko) и P(na, t\nao) — функции Грина для классического и квантового операторов Лиувилля соответственно dF(Q,t\Qo) -№(Q,t\Qo) = 5(Q -Qo)5(t),(9)

9tP(n,t\no) -C^P(n,t\no) = 5и,по5(t).(10)

Для упрощения записи в последних уравнениях опущены индексы классических и квантовых координат, а в правой части использованы стандартные обозначения для символа Кронекера и δ-функции Дирака. Отметим, что квантово-классическое расщепление (8) возможно только при определенных условиях, накладываемых на лиувиллиан молекулярной системы. Далее мы рассмотрим несколько известных моделей динамики квантовых и классических степеней свободы, которые допускают такое расщепление. Модели являются стандартными и широко применяются для описания релаксационных процессов в электрон-донорно-акцепторных соединениях [7].

• 1.1.    Модели внутримолекулярной релаксации и релаксации среды

Найдем выражения для F(Q,t\Q _o ) и P(n,t\n _o ), отождествляя координаты Q и n с поляризационными модами среды и внутримолекулярными колебаниями. В модели Зусмана [19] эволюция системы вдоль координаты Q носит характер диффузии в гармоническом потенциале (1) с коэффициентом диффузии D = к _в Т/ т _к , где т _к - время релаксации классической моды. Вид функции Грина для соответствующего оператора известен [7]

F (Q,t\Q o ) = я(д - Q, (Q o - Q)e ^{- t/ T} " , 2к в Т (1 - e ^{-2 t/ T} ^ )) ,

где JV(х, а, о^-2 ) — функция Гаусса

Я (х,а, о ² ) = ^ 2по ² ) ^-1 / ² _е ^- ( ^{ж -} " ) ^{2 / 2о} 2 .

Уравнение (11) описывает релаксацию к равновесному состоянию системы, и при t ^ то воспроизводит распределение Больцмана вдоль координаты Q на параболической поверхности. Функция F (Q,t|Q ₀ ) также используется для программной генерации стохастических траекторий частиц на диабатических поверхностях свободной энергии в рамках численного метода броуновского моделирования [3; 16].

Найдем далее функцию Грина Р (n,t|n ₀ ) для квантовой подсистемы, принимая во внимание возбуждение и релаксацию высокочастотных колебаний в ходе электронного переноса. Рассмотрим стандартный одноквантовый механизм внутримолекулярной релаксации, когда переход из колебательно-возбужденного состояния возможен только на ближайшее нижележащее состояние, переходы на все другие уровни запрещены. В этой модели Р (n,t|n ₀ ) удовлетворяет дифференциальному уравнению вида

Р ^Ы = —Р ^Ы + ,„ Р (п + 1,tlno),

—У где -Уп) — скорость переходов из п-го колебательного состояния. Решение системы ОДУ

• (13)    может быть найдено методом вариации постоянных. Представим это решение в виде

  Р (п, t|n o ) = C „ (t)exp

  
  

  (^- й ’

  
  
  
  

где C _n (t) — неизвестная функция, подлежащая нахождению. Подставим (14) в (13) и получим общую рекуррентную формулу для определения C _n (t)

C _n (t) =

- (n +l) ^Т У

f

I C_n +! (t')exp <  0

1         1 J ,‘

- ■    -        ' /

В модели Лина скорость колебательной релаксации пропорциональна номеру колебательного подуровня (фактически, амплитуде внутримолекулярных колебаний). Время распада населенности п-го уровня - У^{п )} в этой модели обратно пропорционально п и может быть представлено в виде

^{( / )}

т ( ^{п )} = . (16) п

Скорость внутримолекулярного перераспределения энергии колебаний здесь полностью определяется параметром т^ ) , имеющим смысл времени жизни 1-го колебательно-возбужденного состояния относительно безызлучательного перехода в основное состояние. С учетом зависимости (16) из уравнения (15) находим

C _n (t) =

п + 1

^Т У

[ C n+^'K¹'/ ^- !'" dt.

Это уравнение может быть решено путем последовательного интегрирования для всех п от п0 до 1. В результате получим j no n               к

Рц (n,t|n o ) =

^п 0 ! ^        ^{( - 1)}        _е - ( п + к ) 1/т <^г)

(п = 1, 2,...,п 0 ).     (18)

п! ^ к\(п ₀ — п -к).

Pitting)

n o - 1

= E k=g

(-1) ^fc n g !

(k + 1)!(ng — к — 1)

^ 1 - _g - ( ^ +1) ^/^ » )

Индекс l здесь используется для обозначения релаксационной модели Лина. Выражения (18), (19) описывают внутримолекулярную релаксацию из n g -го начального колебательного состояния путем последовательных переходов частиц на нижележащие колебательные уровни с накоплением в состоянии n = 0. Ниже мы рассмотрим методы расчета нестационарных спектров люминесценции с использованием этих выражений.

Проанализируем далее другую известную модель одноквантовой релаксации, в которой скорость колебательных переходов «вниз» не зависит от номера n т^) = Т,) = const.

В литературе такая модель часто упоминается как «лестничная» модель колебательной релаксации (stepladder model). Применяя описанный выше метод интегрирования системы ОДУ (13) с учетом (20), получим следующий результат

P _s (n,t|n g ) =

g-t/т^

(n g — n)

( \ no-n k)

(n = 1, 2, . . . ,n g ).

n o - 1

P _s (0,t|n g ) = 1 — E k=g

e^t/t^

^{( n} g ^{— k —} 1)

( \ no - n

• 1.2.    Спектральная функция Грина

Рассчитаем профиль спектра люминесценции S(^,t|Qg,ng), формируемый излучательными переходами между выделенными электронными состояниями |ф^) и ф), в случае, когда динамика распределения частиц в начальном электронном состоянии задана функцией вида q(Q, n, t|Qg, ng) = F(Q,t|Qg)P(n,t|ng). Здесь и далее для упрощения записи мы опускаем индексы электронных состояний, а также квантовых и классических мод. Профиль S(w,tlQg,ng) в этом случае может интерпретироваться как спектральная функция Грина, так как представляет собой отклик выделенного оптического перехода на резонансное возбуждение в состояние (Qg,ng). Согласно (5), функция S(ш, t|Qg, ng) рассчитывается как интеграл no

S ^.t Q ,n I V г ,. ² / dQ q ( Q, n,t|Q g , n g ) E⁵( ^AGk m —

П ^ .     (23)

n=0

Данный интеграл может быть вычислен аналитически с использованием известных свойств дельта-функции, результат вычисления имеет вид no

S (w,t\Q0,n0) = Лд EE |^",m|2 Q^nm + ^Ш, n, 11 Qg, ng).(24)

n =0 m

Здесь Q n_m — координата точки пересечения двух ПСЭ, вычисляемая из соотношения (3). Используем теперь разложение функции q ( Q, n, t|Q _g , n _g ) на отдельные компоненты. Подстановка F (Q,t|Q _g ) и P (n,t|n _g ) вместо q ( Q, n, t|Q _g , n _g ) в уравнение (24)

дает следующий результат no

S( w ,t\Q ₀ ,n ₀ ) = Л д ^F (n,t\n ₀ ) ^2 | ^ n,m |² X n =0           m

X ЯQ nm - Q + M (Q o - QK^{t/X L} , 2к в T (1 - -    ' )) .

Полученное выражение представляет собой искомую функцию Грина. В зависимости от модели колебательной релаксации, для F (n,t\n ₀ ) в уравнении (25) могут использоваться выражения (18)–(19) или (21)–(22). При наличии нескольких оптических переходов в системе функция S( w ,t\Q ₀ ,n ₀ ) вычисляется как сумма соответствующих компонент.

Выражение (25) может использоваться для расчета полного спектрального профиля люминесценции системы s(w,t). Пусть распределение частиц в возбужденном состоянии в момент t = 0 (непосредственно после импульса накачки) задано функцией P _ex (Q,n). Тогда величина s( m ,t) может быть вычислена через интеграл по начальному распределению

s( w ,t) = ^ / ^d Q 0 P ex (Q 0 ,^ 0 )S( w ,t\Q 0 ,n 0 ).

ⁿ 0 -^

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

P ex (Q,n) = 5 „,„ o я ( q,Q , 2к в Т ) , (27)

где величина Q обозначает координату минимума ПСЭ основного состояния, а параметр п0 определяет номер колебательного подуровня, на который происходит оптический переход. Задача сводится к вычислению интеграла (26) с подынтегральной функцией, представляющей собой произведение двух гауссиан. Интеграл может быть рассчитан аналитически, в итоге получим следующее выражение no

/ ^Qn,m + ^^Ш _— Г^" _e-t/r_L

V2kT   V kT

s( w ,t) = Л0'^ F(n,t\n 0 )^2 ^|^ П,т ^| exp n =0           m

Структура выражения позволяет легко интерпретировать полученный результат. Согласно (28) спектральный профиль неравновесной системы представляет собой набор люминесцентных полос, соответствующих излучательным переходам между состояниями \ф г ,n) и |ф j ,m). Ожидаемым результатом здесь является зависимость интенсивности этих полос от кинетики населенности колебательных состояний F(n,t\n ₀ ), что отражает внутримолекулярную релаксацию системы. Другое важное следствие уравнения (28) — это смещение и изменение профиля самих полос во времени. Легко видеть, что это явление связано с релаксационными процессами в окружении макромолекулы, так как динамический компонент спектра s( m ,t) напрямую зависит от величины т ^ .

Заключение

Люминесценция молекулярных систем со сверхбыстрыми реакциями переноса электрона часто характеризуется сложной динамикой, на которую оказывают влияние не только многостадийный внутримолекулярный перенос заряда, но и релаксационные процессы в реагентах и среде. Результаты спектроскопических экспериментов в этих случаях могут быть корректно интерпретированы только в рамках моделей, учитывающих все основные механизмы спектральной динамики. Данное исследование посвящено построению математической модели релаксационных процессов внутри макромолекулы и в окружающей среде, позволяющей оценить влияние этих процессов на спектры люминесценции. В работе рассмотрен ряд стандартных подходов к описанию релаксации классических (низкочастных) и квантовых (высокочастотных) колебательных координат, получены аналитические выражения для функций Грина соответствующих эволюционных операторов. Разработаны эффективные схемы расчета нестационарных спектров, которые могут использоваться для анализа экспериментальных данных, полученных методами фемтосекундной спектроскопии.