Role of nonlocal electrostatic effects in the stabilization of monovalent cations in an aqueous cavity surrounded by a weakly polar environment

Full Text

ВВЕДЕНИЕ

Ионный транспорт сквозь каналы в биомембранах

Проблема объяснения, каким образом ионы входят внутрь биомембраны из водного раствора электролита, преодолевая ее выталкивающее действие ("входной барьер”) и стабилизируясь в центре полостей внутри ее каналов, заполненных водой и окруженных слабополярной диэлектрической средой, является фундаментально важной для понимания функционирования ионных каналов [1, 2]. Отметим, что термин “барьер” является общепринятым для специалистов, изучающих ионные каналы клеточных мембран [1]. Вход иона в водную полость канала является первым этапом перед его последующим переходом внутрь селективного фильтра канала, а затем внеклеточного раствора. Поясним природу происхождения этого барьера при входе иона в канал. В растворе ион окружен полярной средой, которая сильно ослабляет электрическое поле иона, за счет чего понижается энергия иона на величину энергии сольватации. При переносе иона внутрь водной полости, находящейся в середине канала, часть пространства вокруг иона за пределами канала занята малополярной средой, которая гораздо слабее ослабляет поле. В результате при расчете в рамках локальной электростатики энергия сольватации существенно уменьшается, что означает большое положительное значение свободной энергии для перехода иона из объема раствора внутрь полости [2], следовательно, пренебрежимо малую вероятность нахождения ионов внутри канала. Такой результат находится в противоречии с экспериментальными наблюдениями о вполне измеримых потоках ионов по каналу между средами в контакте с биологической мембраной – раствором электролита и внутриклеточным раствором.

В работе [3] была впервые расшифрована молекулярная структура бактериального калиевого канала KcsA K⁺ и заложены основы понимания его ионной проводимости и селективности. Было обнаружено [3] существование водной полости в центре KcsA K⁺- канала, которая окружена спиралями белков, отделяющими ее от фосфолипидов биомембраны с меньшей диэлектрической проницаемостью. Кроме этих белковых спиралей имеются четыре α-спирали белков, направленные в центр полости, как это было показано на рис. 1 работы [2]. Электростатическое влияние центральной полости и поровых α-спиралей в KcsA K⁺-канале биомембраны на стабилизацию катиона K⁺ в полости было проанализировано теоретически в работе [2] с использованием конечно-разностного метода решения уравнения Пуассона. Было показано, что полость и α-спирали помогают катиону преодолевать выталкивающий барьер биомембраны, обусловленный уменьшением энергии сольватации иона внутри полости по сравнению с раствором. При этом авторам [2] для полного объяснения стабилизации иона K⁺ в полости пришлось провести расчет энергии системы в рамках мульти-ионной модели K⁺-канала. Было выдвинуто предположение об одновременном движении трех катионов K⁺, причем два из них находятся в селективном фильтре канала, в то время как третий ион – в центре водной полости, находящейся в середине канала. В расчетах работы [2] использовалась локальная электростатика (LE), при этом диэлектрическая проницаемость в водной полости канала принималась равной 80. Структура калиевого KcsA K⁺- канала во многом сходна с другими калиевыми каналами, поэтому появился ряд работ [4–9], в которых анализировалась селективность, проводимость и выпрямление в его водной полости.

Наша работа была инициирована вопросом: обязательно ли использовать для объяснения высокой ионной проницаемости биоканалов столь сложную мульти-ионную модель, т.е. не могут ли другие факторы приводить к резкому уменьшению высоты электростатического барьера при входе иона внутрь канала?

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

Свойства растворов электролитов внутри водных полостей нанометрового масштаба

Интерес к диэлектрическим свойствам межфазной воды [10] вызван тем, что они играют определяющую роль в таких явлениях, как гидратация поверхностей твердых электронных проводников и непроводящих материалов, сольватация ионов и межионные взаимодействия, а также молекулярный и ионный транспорт через нанопоры. Экспериментальные исследования [11–13] показали, что межфазная вода может образовывать упорядоченную структуру с малой поляризуемостью. Однако величины диэлектрической проницаемости и толщины слоя такой воды оставались неизвестными из-за сложности измерений. Значительный прогресс в исследовании свойств межфазной воды был достигнут в работе [14]. В ней экспериментально было доказано существование между плоскими стенками, разделенными различными расстояниями вплоть до 1 нанометра, слоя воды с исчезающе малой дипольной поляризацией, так что эффективная диэлектрическая постоянная этого слоя составляет около 2, что близко к оптической диэлектрической постоянной воды за счет поляризуемости электронной подсистемы.

Изменение энергии иона при его переходе из свободного раствора в сферическую водную полость, окруженную диэлектриком: расчет в рамках локальной электростатики (LE)

Ввиду указанных выше данных о существенно пониженных величинах эффективной диэлектрической постоянной воды внутри полостей малого размера, необходимо прежде всего исследовать, как этот параметр влияет на высоту барьера для входа сферического иона радиуса а из объема раствора внутрь полости малого радиуса. Рассмотрим этот эффект, если его описывать в рамках LE, т.е. рассматривать воду как имеющую различные диэлектрические постоянные: ε_sol = 80 в макроскопическом объеме, но иное (значительно меньшее) значение ε_cav внутри полости радиуса R, которая окружена локальным диэлектриком с диэлектрической постоянной ε_p (рис. 1а). Нижний индекс “cav” означает полость от английского “cavity”.

 

Рис. 1. (а) Схематическое изображение системы: ион радиуса а в центре сферической полости радиуса R (указана ее диэлектрическая проницаемость – в общем случае описываемая нелокальной диэлектрической функцией ε_cav(k)), окруженной диэлектрической средой с проницаемостью ε_р. (б) Изменение свободной энергии Гиббса (ΔG = –ΔW_, в единицах kT) при переходе иона K⁺ из раствора в центр полости в зависимости от диэлектрической проницаемости белкового окружения водной полости ε_p. Расчет (ΔG = ΔG_LE) по формулам (1.1)–(1.2) модели локальной электростатики. Диэлектрическая проницаемость внутри полости ε_cav равна 80 (кривая 1) или 4.9 (кривая 2). Радиус полости R = 5Å. Радиус катиона K⁺ (а = 1.49Å) взят по шкале Гурари и Адриана.

 

Тогда изменение свободной энергии иона ΔG_LE, т.е. работа по его переносу из объема раствора в центр полости может быть вычислена как разность энергий сольватации иона в объеме раствора W_sol и в центре полости W_cav:

$\Delta G_{\mathrm{LE}}=-\Delta W_{\mathrm{LE}}=W_{\mathrm{sol}}-W_{\mathrm{cav}}$, (1.1)

для которой расчет для локальной модели дает формулу (1.2):

$\Delta W_{\mathrm{LE}}\left(R,a\right)=\frac{e^2}{2R}\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{p}}}\right)+\frac{e^2}{2a}\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{sol}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right)$, (1.2)

где e – заряд иона, параметры a и R определены на рис. 1а. Нижний индекс LE в формулах (1.1)-(1.2) означает “локальная электростатика”.

Результаты расчетов по этой формуле показаны на рис. 1б.

Результат работы [2] получается, если диэлектрические постоянные воды в растворе (ε_sol) и внутри полости (ε_cav) равны 80, когда в формуле (1.2) остается только первый член. Как показывает кривая 1 на рис. 1б для этого случая, высота энергетического барьера для переноса иона внутрь полости составляет десятки kT при предполагаемых значениях диэлектрической постоянной среды, окружающей полость (ε_p) – в соответствии с работой [2].

Если в соответствии с данными о диэлектрических свойствах воды внутри малых объемов использовать пониженную величину ее эффективной диэлектрической постоянной внутри полости ε_cav (например, ε_cav= 4.9, что отвечает диэлектрической проницаемости воды в промежуточном интервале длин волн k, см. ниже), то высота барьера увеличивается еще в несколько раз (кривая 2 на рис. 1б).

Таким образом, можно сделать вывод, что расчеты высоты электростатического барьера для иона в рамках модели локальной электростатики исключают стабилизацию одиночного иона внутри полости – не только при ε_cav= 80, но и при пониженных значениях этого параметра, когда высота электростатического барьера ΔG увеличивается еще в несколько раз.

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

Такое описание не учитывает эффектов структуры полярной среды, возникающих как за счет конечных размеров ее молекул, так и короткодействующих взаимодействий между ними. В результате поляризация среды P в точке r определяется не только электрическим полем Е в этой же точке r, но и значениями поля в области, окружающей эту точку. В результате пространственные распределения P(r), Е(r) и D(r) (индукция: D = Е + 4π Р) связаны между собой линейными интегральными соотношениями, которые представляют основу нелокальной электростатики (NE).

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

Таким образом, для описания характеристик системы, в частности величин энергий сольватации ионов как в объеме раствора (W_sol), так и внутри полостей, заполненных средой (W_cav), необходимо применять соотношения NE. Такие нелокально-электростатические расчеты ΔW для перехода иона в сферическую водную полость (рис. 1а) выполнены в данной работе.

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

Отметим, что впервые на необходимость использования NE для рассмотрения процесса стабилизации иона K⁺ в полости KcsA K⁺-канала было указано в работе [15], однако детального расчета представлено не было. При NE-описании этого процесса в [16, 17] оказывается, что изменение свободной энергии Гиббса ΔG при переходе иона из раствора в сферическую полость, т.е. высота электростатического барьера, понижается по сравнению с расчетом на основе локальной электростатики. Аналогичный вывод был сделан при использовании NE в работах [18, 19] в канале с плоскопараллельными стенками.

В работе Бардхана (Bardhan) [16] изменение энергии сольватации иона внутри сферической полости радиуса 13Å по сравнению с объемом раствора (ΔW) рассчитывалось численными методами; не было получено какой-либо аналитической формулы для этой характеристики. В нелокальном уравнении Пуассона использовалось ядро ε(r–r′):

$\mathrm{\epsilon }\left(\mathit{r}\mathbf{-}{\mathit{r}}^{\mathbf{\text{'}}}\right)={\mathrm{\epsilon }}^{*}\delta \left(\mathit{r}\mathbf{-}{\mathit{r}}^{\mathbf{\text{'}}}\right)+\left({\mathrm{\epsilon }}_w-{\mathrm{\epsilon }}^{*}\right)\text{\hspace{0.17em}\hspace{0.17em}}\frac{\exp \left(-\left|\left(\mathit{r}\mathbf{-}{\mathit{r}}^{\mathbf{\text{'}}}\right)\right|/\lambda \right)}{4\pi {\lambda }^2\left|\left(\mathit{r}\mathbf{-}{\mathit{r}}^{\mathbf{\text{'}}}\right)\right|}$, (2.1)

где выбирались следующие значения параметров: ε_w = 80 – диэлектрическая постоянная воды; ε^* = 1.8, т.е. равна величине высокочастотной (“оптической”) диэлектрической проницаемости воды (ε₁ = 1.8, см. ниже про трехмодовую модель); λ = 24Ả [16].

Как было отмечено в работе [17], функциональный вид ядра (2.1) соответствует для фурье-образа диэлектрической функции воды ε(k) ее однополюсной/одномодовой модели (1M) ε_1M(k) [20]:

${\mathrm{\epsilon }}_{1\text{M}}\left(k\right)={\mathrm{\epsilon }}_{\text{o}}+\left(\mathrm{\epsilon }-{\mathrm{\epsilon }}_{\text{o}}\right){\left[1+{\left(k\mathrm{\Lambda }\right)}^2\left(\mathrm{\epsilon }/{\mathrm{\epsilon }}_{\text{o}}\right)\right]}^{-1}$, (2.2)

которое можно также записать в виде [20, 21]:

$\frac{1}{{\mathrm{\epsilon }}_{1\mathrm{M}}\left(k\right)}=\frac{1}{{\mathrm{\epsilon }}_{\text{o}}}-\left(\frac{1}{{\mathrm{\epsilon }}_{\text{o}}}-\frac{1}{\mathrm{\epsilon }}\right)\frac{1}{1+{\left(k\mathrm{\Lambda }\right)}^2}$, (2.3)

если положить в них ε = ε_w, ε_o = ε^*, Λ = λ. Однако для моделирования диэлектрических свойств воды величина корреляционной длины Λ должна составлять несколько Ангстрем (например, Λ = 5Ả в работе [20]), т.е. выбранное в работе [16] значение 24Ả для параметра Λ в формуле (2.1) слишком велико, чтобы используемое в [16] ядро (2.1) соответствовало принятой в NE одномодовой диэлектрической функции (2.2) для расчета электростатических полей внутри водной фазы.

Далее, величина ε_о в формулах (2.2) и (2.3) должна соответствовать диэлектрической постоянной воды в зоне прозрачности, разделяющей дебаевские и инфракрасные моды [21], величина которой составляет ε_о = ε₂ = 4.9 [20, 21], тогда как формула (2.1) была использована в работе [16] для гораздо меньшего значения параметра ε^* = ε₁ = 1.8.

Ввиду отсутствия в [16] аналитической формулы для ΔW, слишком большого радиуса полости, а также неудачно выбранных значений параметров (λ и ε_o) в ядре (2.1), результаты расчетов энергии сольватации иона в полости, полученных в [16], невозможно использовать для анализа стабилизации иона в водной полости радиуса 5Å, который рассматривался в работе [2] и анализируется в настоящей работе.

В работе [17] рассматривался точечный заряд, находящийся в центре полости, заполненной средой с одномодовой диэлектрической функцией, причем ее параметры считались тождественными их величинам для бесконечного водного раствора. В отличие от работы [16] была выведена аналитическая формула (Δw_DA, см. ниже) для расчета энергии переноса иона из объема раствора в центр сферической полости в рамках “диэлектрического приближения” (DA) для нелокального диэлектрического отклика, которое было введено в работах [22, 23]. При расчетах по этой формуле работы [17], как и в работе [16], значение ε_o в (2.3) считалось равным высокочастотной диэлектрической проницаемости воды ε_o = ε₁ = 1.8, тогда как этот параметр следует выбирать равным ε₂ = 4.9, как это сделано в настоящей работе.

Диэлектрическое приближение для нелокального диэлектрического отклика (см. ниже) очень широко использовалось при NE-описании пространственно-ограниченных полярных сред [17, 22–37]. Следует, однако, отметить его принципиальный недостаток: даже для систем с высокой симметрией (плоской, цилиндрической или сферической) результаты для распределения потенциала или энергий сольватации удается получать в виде аналитических выражений только для одномодовой аппроксимации диэлектрической функции среды (2.2) или (2.3).

Из-за предположения в работе [17] о том, что параметры диэлектрического отклика воды внутри полости тождественны их величинам для воды вокруг иона в растворе, в выражении для изменения энергии сольватации ΔW этой работы отсутствует NE-аналог второго слагаемого в правой части формулы (1.2), обусловленного специфическими свойствами воды внутри полости. С другой стороны, для расчета [17], где ион симулировался точечным зарядом, предположение о совпадении диэлектрических функций воды (как минимум, совпадение значений диэлектрической постоянной ε₂) в объеме раствора и внутри полости является обязательным условием для получения конечной величины энергии переноса иона.

В последние годы был разработан альтернативный подход [38], который позволяет получать аналитические результаты для систем, где полярная среда занимает ограниченный объем, при произвольной модели фурье-образа ее диэлектрического отклика во всем пространстве ε(k), что позволяет, в частности, проводить расчеты для трехмодовой модели диэлектрической функции среды. Этот подход в работах [39, 40] получил название “обратное диэлектрическое приближение” (IDA), которое мы будем использовать далее. Краткий обзор этого метода можно найти в Приложении 1.

Цель данной работы – вывести в рамках NE аналитические выражения для изменения энергии сольватации иона в центре полости по сравнению с объемом раствора (ΔW ) как для одно-, так и трехмодовой модели диэлектрической функции воды внутри полости ε_cav(r₁, r₂) (которая может совпадать с этой характеристикой в объеме раствора ε(r₁, r₂) или отличаться от нее), причем для иона конечного размера, и проанализировать, при каких ее параметрах величина свободной энергии переноса иона из объема раствора внутрь полости ΔG понижается до уровня, приводящего к значительному увеличению вероятности нахождения иона внутри подобной полости.

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

Методы нелокальной электростатики для пространственно ограниченных сред. Сравнение общих соотношений для подходов, основанных на диэлектрическом приближении [23], и обратного диэлектрического приближения [38]

Решение уравнений NE для среды, занимающей все пространство. Отметим разницу между нелокальной и локальной электростатикой, которую необходимо учитывать при расчете энергии сольватации иона внутри полости (подробнее см. в [38]). При локальном подходе среда рассматривается как локальный и однородный (за исключением объемов, занятых растворенными веществами или другими средами) диэлектрик, характеризуемый своей диэлектрической постоянной ε_s. Это означает, что электрическое поле E(r) (точнее, напряженность электрического поля) в точке r и его индукция (electric displacement) D(r) в той же точке пропорциональны между собой: D(r)= ε_sE(r). Такое описание вполне оправдано для полей (создаваемых внешними по отношению к среде источниками поля), которые мало изменяются на молекулярных размерах. Однако в электрохимических системах это условие часто нарушается – например, для полей вокруг иона малого размера [20, 21, 38, 41, 42]. Поэтому поляризация P в точке r₁ определяется электрическим полем Е не только в той же точке r₁, но и ее окрестности, что математически выражается в виде нелокальной связи между величинами Р и Е (пределы интегрирования в (3.1) – все пространство, включая внутренность иона):

$P_{\alpha }\left({\mathit{r}}_{\mathbf{1}}\right)=\int {\chi }_{\mathrm{\alpha \beta }}\left({\mathit{r}}_1\mathbf{,}{\mathit{r}}_2\right)E_{\mathrm{\beta }}\left({\mathit{r}}_2\right)\text{\hspace{0.17em}d}{\mathit{r}}_{\mathbf{2}}$, (3.1)

где χ_αβ(r₁, r₂) – тензор нелокальной (диэлектрической) восприимчивости, который отличен от нуля, только если расстояние между точками r₁ и r₂ не превышает (максимальной) длины корреляции системы (Λ или λ₃) [20, 21, 38, 41, 42], а по повторяющимся индексам β ведется суммирование от 1 до 3. Используя соотношение: D(r) = E(r) + 4π Р(r), получаем нелокальное выражение для индукции через распределение электрического поля (пределы интегрирования в (3.2) – все пространство, включая внутренность иона):

$D_{\alpha }\left({\mathit{r}}_1\right)=\int {\epsilon }_{\mathrm{\alpha \beta }}\left({\mathit{r}}_1,{\mathit{r}}_2\right)E_{\mathrm{\beta }}\left({\mathit{r}}_2\right)\text{d}{\mathit{r}}_2$, (3.2)

где (тензорная и нелокальная) диэлектрическая функция системы ε_αβ(r₁, r₂) дается формулой:

${\epsilon }_{\alpha \beta }\left({\mathit{r}}_1,{\mathit{r}}_2\right)={\delta }_{\alpha \beta }\delta \left({\mathit{r}}_1-{\mathit{r}}_2\right)+4\pi {\chi }_{\alpha \beta }\left({\mathit{r}}_1,{\mathit{r}}_2\right)$. (3.3)

Для нахождения распределения потенциала диэлектрическая функция системы должна быть задана во всем пространстве, т.е. в областях, занятых как полярной средой (водой), так и локальным диэлектриком. Можно ввести обратную (тоже тензорную и нелокальную) диэлектрическую функцию системы ε^–1_αβ(r₁, r₂), связанную интегральным соотношением с функцией ε_αβ(r₁, r₂), как это было сделано в [20, 21, 38, 41, 42] (пределы интегрирования в (3.4) и (4) – все пространство, включая внутренность иона):

$\int {\epsilon }_{\mathrm{\alpha \beta }}^{-1}\left({\mathit{r}}_1,{\mathit{r}}_2\right){\epsilon }_{\mathrm{\beta \gamma }}\left({\mathit{r}}_2,{\mathit{r}}_3\right)\text{\hspace{0.33em}d}{\mathit{r}}_2={\delta }_{\mathrm{\alpha \gamma }}\delta \left({\mathit{r}}_1-{\mathit{r}}_3\right)$. (3.4)

Тогда можно выразить электрическое поле через индукцию:

$E_{\alpha }\left({\mathit{r}}_{\mathbf{1}}\right)=\int {\epsilon }_{\mathrm{\alpha \beta }}^{-1}\left({\mathit{r}}_1,{\mathit{r}}_2\right)D_{\mathrm{\beta }}\left({\mathit{r}}_{\text{2}}\right)\text{\hspace{0.33em}d}{\mathit{r}}_2$. (4)

Соотношение (3.2) или (4) дополняется стандартными уравнениями электростатики:

$\text{div\hspace{0.17em}}\mathit{D}\left(\mathit{r}\right)=4\mathrm{\pi }\mathrm{\rho }\left(\mathit{r}\right)\text{,\hspace{1em}}\mathit{E}\left(\mathit{r}\right)=-\text{grad\hspace{0.17em}}\mathrm{\phi }\left(\mathit{r}\right)$, (5)

где ρ(r) – плотность внешних зарядов, φ(r) – потенциал электрического поля. Для нахождения всех этих распределений потенциала должна быть задана диэлектрическая функция системы (во всем пространстве) и распределение плотности зарядов ρ(r).

Решение задачи может быть найдено в аналитическом виде для любой диэлектрической функции среды и любого распределения внешнего заряда при условиях однородности среды и отсутствия межфазных границ: среда занимает все пространство, а внешние заряды распределены внутри среды без нарушения ее диэлектрических свойств [20, 21, 41, 42].

Решение уравнений NE для среды, занимающей лишь определенную пространственную область, в рамках “диэлектрического приближения” [23]. Расчет электрического поля E(r) и его потенциала φ(r) является гораздо более сложным, если полярная среда занимает лишь определенную пространственную область V, которая для рассматриваемой в настоящей работе системы заключена между двумя концентрическими сферами радиусов a и R (рис. 1а). Такая задача для систем, в которых необходимо учитывать эффект вырезания полости иона из растворителя, раньше решалась в основном в рамках “диэлектрического приближения”, предложенного в работах [22, 23]. Во-первых, как и во всех работах по этой проблеме, предполагается, что отсутствует корреляция между распределениями поляризации (3.1) в пространственных областях, занятых различными средами; в частности внутри областей, занятых локальным изотропным диэлектриком р с диэлектрической постоянной ε_р, выполняются локальные соотношения:

ε_αβ(r₁, r₂) = ε_р δ_αβ δ(r₁ – r₂),

ε^–1_αβ(r₁, r₂) = (ε_р)^–1 δ_αβ δ(r₁ – r₂),

тогда как внутри области V, занятой нелокальной полярной средой (водой), выполняются интегральные соотношения (3.1), (3.2), (3.4) и (4), где обе пространственные координаты r₁ и r₂ лежат в области V, т.е. интегрирование выполняется только по области V. Во-вторых, предполагается, что нелокальная диэлектрическая функция среды в формуле (3.2) имеет такой же вид, что и для этой среды в неограниченном пространстве: ε_αβ(cav)(r₁, r₂) = ε_αβ(cav)(r₁–r₂). В частности, для рассматриваемой системы (cм. рис. 1а), где полярная среда занимает область V в форме сферического пояса между поверхностью иона в центре полости (r = a) и границей полость/внешний диэлектрик (r = R), формула (3.2) принимает вид (интегрирование в (6.1) проводится по объему V):

_{$\begin{array}{c}{D}_{\alpha }\left({\mathit{r}}_1\right)=\underset{\text{V}}{\int }{\mathrm{\epsilon }}_{\mathrm{\alpha \beta }\left(\text{cav}\right)}\left({\mathit{r}}_1-{\mathit{r}}_2\right)E_{\mathrm{\beta }}\left({\mathit{r}}_2\right)\text{d}{\mathit{r}}_2,\\ a<r_1<R, a<r_2,\end{array}$} (6.1)

где (как и ниже) для всех радиус-векторов используется обозначение: r = | r |.

Ее комбинация с уравнениями (5) приводит к интегральному уравнению для электрического поля или интегро-дифференциальному уравнению для его потенциала. Однако нахождение их решений в аналитическом виде (даже для самых простых геометрий системы) возможно по существу только для одномодовой модели диэлектрической функции растворителя (2.2) – (2.3) [23], тогда как обобщение этой теории на более сложные формы диэлектрической функции затруднительно.

Решение уравнений NE для среды, занимающей лишь определенный объем, в рамках “обратного диэлектрического приближения” [38]. По указанной выше причине в работах [38–40] был предложен другой приближенный способ расчета электрических полей в системах, где полярная среда занимает ограниченную часть пространства. Подробнее об этом подходе можно прочитать в разделе “ПРИЛОЖЕНИЕ 1” в конце данной статьи. Вместо подхода [23], основанного на приближении для нелокальной диэлектрической функции (см. выше), предполагается [38–40], что пространственно-временная корреляция пульсаций поляризации нелокальной среды имеет такой же вид, как для этой среды в неограниченном пространстве. Это означает, что для рассматриваемой системы (cм. рис. 1а), где полярная среда занимает область V в форме сферического пояса между поверхностью иона в центре полости (r = a) и границей полость/внешний диэлектрик (r = R), формула (4) принимает вид (интегрирование в (6.2) проводится по объему V):

$\begin{array}{c}{E}_{\alpha }\left({\mathit{r}}_1\right)=\underset{V}{\int }{\mathrm{\epsilon }}_{\mathrm{\alpha \beta }\left(\text{cav}\right)}^{-1}\left({\mathit{r}}_1-{\mathit{r}}_2\right)D_{\mathrm{\beta }}\left({\mathit{r}}_1\right)\text{\hspace{0.17em}d}{\mathit{r}}_2,\\ a<r_1<R, a<r_2<R,\end{array}$ (6.2)

где индукция, создаваемая зарядом иона, имеет вид:

$D_{\mathrm{\beta }}\left(\mathit{r}\right)=\frac{e{r}_{\mathrm{\beta }}}{r^3},r>a$, (7)

а функция ε^–1_αβ(cav)(r₁ – r₂) связана по формуле (8) с фурье-образом диэлектрической функции воды внутри полости ε_cav(k) (пределы интегрирования во второй строчке (8) – все k пространство):

$\begin{array}{c}{\mathrm{\epsilon }}_{\mathrm{\alpha \beta }\left(\text{cav}\right)}^{-1}\left({\mathit{r}}_1-{\mathit{r}}_2\right)\equiv {\delta }_{\mathrm{\alpha \beta }}{\mathrm{\epsilon }}_{\text{cav}}^{-1}\left(\left|{\mathit{r}}_1-{\mathit{r}}_2\right|\right)=\\ ={\delta }_{\mathrm{\alpha \beta }}{\left(2\mathrm{\pi }\right)}^{-3}\int {\left[{\mathrm{\epsilon }}_{\text{cav}}\left(k\right)\right]}^{-1}\exp \left[i\mathit{k}\left({\mathit{r}}_1-{\mathit{r}}_2\right)\right]\text{\hspace{0.17em}d}\mathit{k}=\\ ={\delta }_{\mathrm{\alpha \beta }}{\left(2{\mathrm{\pi }}^2\left|{\mathit{r}}_1-{\mathit{r}}_2\right|\right)}^{-1}\underset{0}{\overset{\infty }{\int }}{\left[{\mathrm{\epsilon }}_{\text{cav}}\left(k\right)\right]}^{-1}k\sin \left[k\left|{\mathit{r}}_1-{\mathit{r}}_2\right|\right]\text{\hspace{0.17em}d}k,\\ a<r_1<R, a<r_2<R\end{array}$ (8)

Формула, аналогичная формуле (8), выведена в [20, 21] для связи функции ε_αβ(r₁ – r₂) в “диэлектрическом приближении” с ε(k).

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

$\begin{array}{c}{E}_{\alpha }\left(\mathit{r}\right)=\left(r_{\alpha }/r\right)E\left(r\right),E\left(r\right)=-\text{d}\mathrm{\phi }\left(r\right)/\text{d}r,\\ \mathrm{\phi }\left(r\right)=\mathrm{\phi }\left(R\right)+\underset{r}{\overset{R}{\int }}E\left(r\right)dr.\end{array}$ (9)

Формулы (6.2)–(9) показывают важнейшее преимущество подхода обратного диэлектрического приближения (IDA) [38] по сравнению с более ранним методом диэлектрического приближения (DA) [23]: для систем с высокой симметрией (например, сферически симметричных, как в данной работе) новый подход позволяет получать аналитические решения для распределения электрического поля и его потенциала при произвольном функциональном виде диэлектрической функции среды ε(k).

Общая формула для распределения потенциала вокруг иона в полости для произвольной модели ε(k), вывод которой основан на подходе обратного диэлектрического приближения (IDA)[38]

В данном разделе выводятся формулы для распределения потенциала φ(r) вокруг иона, находящегося в центре полости (рис. 1а), на основе соотношений (6.2)–(9).

Во внешней области (при r >R) находится локальный диэлектрик, так что

$\begin{array}{c}E\left(r_1\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}{r}_1^2},\mathrm{\phi }\left(r_1\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}{r}_1},r_1>R,\\ \mathrm{\phi }\left(R\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}R}.\end{array}$ (10.1)

Для нахождения выражений для распределения поля и потенциала внутри сферического слоя (a < r < R), заполненного полярной средой, нужно упростить выражения (6.2)–(9). Для этого воспользуемся преобразованиями, проведенными в работе [38] для распределения поля вокруг иона радиуса а, находящегося в неограниченной нелокальной среде. Выведенная формула (17) работы [38] для распределения поля сохраняет тот же вид, за исключением ограничения области интегрирования по координате r₂ в соответствии с формулой:

$E\left(r_1\right)=\frac{1}{r_1}\left[\underset{a}{\overset{R}{\int }}{\mathrm{\epsilon }}_{l=1}^{-1}\left(r_1,r_2\right)\frac{e}{r_2^2}{r}_{2}d{r}_2\right],a\le r_1,r_2\le R$. (10.2)

Функция ε^–1_l=1(r₁, r₂) под интегралом в формуле (10.2) представляет коэффициент разложения обратной диэлектрической функции ε^–1_αβ(r₁, r₂) по сферическим функциям Y_l0 для нижнего индекса l = 1, см. формулу (16) работы [38]. Для нее сохраняется вывод формулы (20) работы [38]:

_{$\begin{array}{c}{\mathrm{\epsilon }}_l^{-1}\left({\mathit{r}}_1,{\mathit{r}}_2\right)=\left(2/\mathrm{\pi }\right)r_1{r}_2\times \\ \times \underset{0}{\overset{\infty }{\int }}{\left[\mathrm{\epsilon }\left(k\right)\right]}^{-1}{j}_l\left(k{r}_1\right)j_l\left(k{r}_2\right)k^{2}dk,\end{array}$} (10.3)

где j_l(х) – сферическая функция Бесселя с индексом l = 1: j₁(x) = x^–1 (cos x – x^–1 sin x). Подстановка выражения (10.3) в формулу (10.2) и последующее интегрирование по формуле:

_{$\underset{a}{\overset{R}{\int }}{j}_l\left(kr\right)dr={\left(k^{2}R\right)}^{-1}\sin \left(kR\right)-{\left(k^{2}a\right)}^{-1}\sin \left(ka\right)$} (10.4)

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

$\begin{array}{c}E\left(r\right)=-\frac{2e}{\mathrm{\pi }}\underset{0}{\overset{\infty }{\int }}\frac{dk}{{\mathrm{\epsilon }}_{\mathrm{cav}}\left(k\right)}\left[\frac{\sin ka}{a}-\frac{\sin kR}{R}\right]j_1\left(kr\right),\\ a\le r\le R.\end{array}$ (10.5)

Подстановка этого выражения в формулу (9) для потенциала φ(r) и интегрирование по координате r по формуле (10.4) дает окончательное выражение для распределения потенциала:

$\begin{array}{c}{\phi }_{cav}\left(r\right)=\frac{e}{{\epsilon }_{p}R}+\\ +\frac{2e}{\pi }\underset{0}{\overset{\infty }{\int }}\frac{dk}{k^2{\epsilon }_{cav}\left(k\right)}\left\{\frac{\sin ka}{a}-\frac{\sin kR}{R}\right\}\left\{\frac{\sin kr}{r}-\frac{\sin kR}{R}\right\},\end{array}$ (11)

в частности, для значения потенциала на поверхности иона:

${\mathrm{\phi }}_{\mathrm{cav}}\left(a\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}R}+\frac{2e}{\mathrm{\pi }}\underset{0}{\overset{\infty }{\int }}\frac{dk}{{\mathrm{\epsilon }}_{\mathrm{cav}}\left(k\right)}{\left\{\frac{\sin ka}{ka}-\frac{\sin kR}{kR}\right\}}^2$. (12)

При R → ∞ формула (11) переходит в формулу (24) работы [38] для распределения потенциала вокруг иона радиуса a в бесконечном растворителе:

$\mathrm{\phi }\left(r\right)=\frac{2e}{\mathrm{\pi }}\underset{0}{\overset{\infty }{\int }}\frac{dk}{k^2\mathrm{\epsilon }\left(k\right)}\left\{\frac{\sin ka}{a}\right\}\left\{\frac{\sin kr}{r}\right\}$.

Одномодовая модель диэлектрической функции полярной среды. Сравнение результатов расчета изменения энергий сольватации иона на основе подходов диэлектрического приближения (DA) [23] и обратного диэлектрического приближения (IDA) [38]. Зависимость энергии сольватации от параметров диэлектрической функции воды в объеме раствора и внутри полости для подхода обратного диэлектрического приближения (IDA) [38]

Формулы для распределения потенциала вокруг иона внутри полости и изменения энергии сольватации иона в полости ΔW, для одномодовой модели ε(k), вывод которой основан на диэлектрическом приближении (DA). В [17] для вывода формулы распределения потенциала вокруг иона в полости, окруженной локальным диэлектриком, использовалось уравнение (6.1). Система уравнений (6.1)–(5) приводит к нелокальному интегро-дифференциальному уравнению Пуассона для электростатического потенциала φ(r), и имеет вид:

$\mathrm{div}\underset{\mathrm{V}}{\int }\mathrm{\epsilon }\left({\mathit{r}}_1-{\mathit{r}}_2\right)\nabla \mathrm{\phi }\left({\mathit{r}}_2\right)d{\mathit{r}}_2=-4\mathrm{\pi \rho }\left({\mathit{r}}_1\right)$. (13.1)

Векторы r₁ и r₂ в (13.1) находятся внутри полости радиуса R: a < r₁, r₂ < R, где a – радиус борновской сферы, моделирующий ион, по поверхности которой распределен заряд с объемной плотностью ρ(r).

В работе [17] было сделано дополнительное приближение: уравнение (13.1) рассматривалось для точечного иона. Для его решения было использовано ядро (2.1), соответствующее одномодовой модели диэлектрической функции (2.2). Тогда уравнение (13.1) сводилось к дифференциальному уравнению (см. уравнение (10) в [17]), решение которого для распределения потенциала φ(r) вокруг точечного иона внутри полости имеет вид:

_{$\begin{array}{c}{\mathrm{\phi }}_{\mathrm{cav}}\left(r,R\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{cav}}r}\times \\ \times \left[1+\left[\left({\mathrm{\epsilon }}_{\mathrm{cav}}/{\mathrm{\epsilon }}_{\mathrm{o}}\right)-1\right]\exp \left(-r/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)\right]+C\left({\mathrm{\epsilon }}_{\mathrm{p}},R\right),\end{array}$} (13.2)

_{$C\left({\mathrm{\epsilon }}_{\mathrm{p}},R\right)=\frac{e}{R}\left[1/{\mathrm{\epsilon }}_{\mathrm{p}}-\left[1/{\mathrm{\epsilon }}_{\mathrm{cav}}+\left(1/{\mathrm{\epsilon }}_{\mathrm{o}}-1/{\mathrm{\epsilon }}_{\mathrm{cav}}\right)\exp \left(-R/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)\right]\right].$} (13.3)

Распределение потенциала φ(r) вокруг точечного иона в полости бесконечного радиуса, но с параметрами такими же, как в полости конечного радиуса, имеет вид:

_{$\begin{array}{c}{\mathrm{\phi }}_{\mathrm{cav}}\left(r,R=\infty \right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{cav}}r}\times \\ \times \left[1+\left[\left({\mathrm{\epsilon }}_{\mathrm{cav}}/{\mathrm{\epsilon }}_{\mathrm{o}}\right)-1\right]\exp \left(-r/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)\right].\end{array}$} (13.4)

Отметим, что в работе [17] значение параметр ε_o было выбрано равным 1.8, а в данной работе при построении графиков принимается: ε_o = 4.9.

Поскольку в работе [17] при решении нелокального уравнения Пуассона в полости рассматривался ион точечных размеров, то для нахождения изменения энергии сольватации иона ΔW в полости надо осуществить предельный переход:

_{$\begin{array}{c}\Delta w({\epsilon }_p,R)=\left(e/2\right)\times \\ \times {\lim }_{r\to 0}\left[{\phi }_{cav}\left(r,R=\infty \right)-{\phi }_{cav}\left(r,R\right)\right].\end{array}$} (13.5)

Вывод формулы для ΔW осуществляется подстановкой соотношений (13.2)–(13.4) для потенциала в полости в формулу (13.5).

При сопоставлении с результатами публикации [17] следует отметить, что формула для Δw_DA работы [17] была выведена для одномодовой модели диэлектрической функции воды ε(k) в неограниченном пространстве (см. формулу (2.2)), тогда как отличная от нее диэлектрическая функция в полости ε_1M(cav)(k) может иметь другие параметры. В работе [17] в формуле (2.2) значение ε_o принималось равным 1.8. Здесь мы используем эту формулу для расчета изменения энергии сольватации при переходе иона радиуса a из раствора в центр водной полости при условии, что ε_o в диэлектрической функции равна ε_o = ε₂ = 4.9, так как такая модель, согласно [20, 21], лучше аппроксимирует расчетные энергии сольватации ионов соответствующим экспериментальным значениям. С учетом этой поправки, выведенную в [17] формулу для ∆w_DA-иона в центре полости радиуса R можно представить в виде:

$\Delta w_{\mathrm{DA}}\left(R\right)=\frac{e^2}{2R}\left\{\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{p}}}\right)+\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right)\exp \left(-R/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)\right\}$. (14)

При ее выводе предполагалось, что ион имеет беконечно малый размер, поэтому в формуле (14) нет зависимости от радиуса иона. Однако для случая, когда параметры воды в полости отличны от их значений в растворе, для расчета ΔW_DA необходимо к Δw_DA прибавить дополнительный член, как это сделано в формуле (15):

$\begin{array}{c}\Delta W_{\mathrm{DA}}\left(R,a\right)=\Delta w_{\mathrm{DA}}\left(R\right)+\\ +\left[w_{1\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)-w_{1\mathrm{M}}\left(a\right)\right]\end{array}$, (15)

где w_1Mcav – энергия сольватации иона радиуса a в бесконечном растворителе, но с параметрами ε_cav и Λ_cav такими же, как в полости:

$w_{1\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)=\frac{e^2}{2a}\left(1-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}+\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right)f\left(2a/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)\right)$, (16)

а w_1M – в свободном растворителе с неизмененными параметрами:

$w_{1\mathrm{M}}\left(a\right)=\frac{e^2}{2a}\left(1-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}+\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{\mathrm{\epsilon }}\right)f\left(2a/\mathrm{\Lambda }\right)\right)$. (17)

Формула (17) имеет тот же функциональный вид, что и формула (16), но как диэлектрическая проницаемость в ε (ε = 80), так и корреляционная длина Λ могут отличаться от их значений в полости. Энергии ΔW_DA приписан нижний индекс DA (диэлектрическая аппроксимация), чтобы отличать ее от формулы для ΔW, вывод которой, приведенный ниже, основан на новом НЭ-подходе работы [38].

Формулы для распределения потенциала вокруг иона в полости и изменения энергии сольватации иона в полости ΔW, для одномодовой модели ε(k), вывод которой основан на теории обратного диэлектрического приближения (IDA) [38]. Используя общую формулу (11), получим выражение для распределения потенциала вокруг иона в полости для одномодовой модели ε(k). Диэлектрическая функции для одномодовой модели воды для раствора ε_1M(k) была определена ранее, см. формулу (2.3). Диэлектрическая функция в полости ε_1Mcav(k) (18) имеет такой же функциональный вид, как ε_1M(k) в свободном растворе, но значения параметров воды в полости: ε_cav (длинноволновая диэлектрическая постоянная) и Λ_cav (корреляционная длина) могут отличаться от их величин для раствора ε и Λ:

$\frac{1}{{\mathrm{\epsilon }}_{1\mathrm{M}\left(\mathrm{cav}\right)}\left(k\right)}=\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right)\frac{1}{1+{\left(k{\mathrm{\Lambda }}_{\mathrm{cav}}\right)}^2}$, (18)

Подставляя в формулу (11) диэлектрическую функцию воды в полости (18) и вычисляя интеграл, получим для распределения потенциала вокруг иона в полости в рамках одномодовой модели диэлектрической функции формулу:

$\begin{array}{c}{\mathrm{\phi }}_{1\mathrm{M}\left(\mathrm{cav}\right)}\left(R,r\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}R}+\frac{e}{{\mathrm{\epsilon }}_{\mathrm{o}}}\left(\frac{1}{r}-\frac{1}{R}\right)-\\ -e\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right)A_{\mathrm{cav}}\left(R,r\right),\end{array}$ (19.1)

где введена вспомогательная функция A_cav(R, r) по формуле:

$A_{\mathrm{cav}}\left(R,r\right)=\frac{g\left({r/\mathrm{\Lambda }}_{\mathrm{cav}},a/{\mathrm{\Lambda }}_{cav}\right)}{r}-\frac{g\left(R/{\mathrm{\Lambda }}_{\mathrm{cav}},r/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)+g\left(R/{\mathrm{\Lambda }}_{\mathrm{cav}},a/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)-f\left(2R/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)}{R},$ (19.2)

где функция g двух аргументов x и y определяется как:

$g\left(x,y\right)=1-\exp \left(-x\right)\frac{\sinh \left(y\right)}{y}$, (19.3)

а функция f аргумента у (в формуле (19.2) этот аргумент равен y = 2R/Λ) имеет вид:

$f\left(y\right)=1-\left[1-\exp \left(-y\right)\right]/y$. (19.4)

Потенциал на границе борновской сферы иона получается подстановкой r = a в (19.1)–(19.2):

_{$\begin{array}{c}{\mathrm{\phi }}_{1\mathrm{Mcav}}\left(R,a\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}R}+\frac{e}{{\mathrm{\epsilon }}_{\mathrm{o}}}\left(\frac{1}{a}-\frac{1}{R}\right)-\\ -e\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right){\mathrm{\Phi }}_{\mathrm{cav}}\left(R,a\right),\end{array}$} (20.1)

${\mathrm{\Phi }}_{\mathrm{cav}}\left(R,a\right)=\frac{f\left(2a/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)}{a}+\frac{{\mathrm{B}}_{\mathrm{cav}}\left(R,a\right)}{R}$, (20.2)

$B_{\mathrm{cav}}\left(R,a\right)=f\left(2R/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)-2g\left(R/{\mathrm{\Lambda }}_{\mathrm{cav}},a/{\mathrm{\Lambda }}_{\mathrm{cav}}\right)$. (20.3)

Потенциал на границе борновской сферы в неограниченном растворителе получается из (20.1) предельным переходом R → ∞ и заменой значений параметров ε_cav и Λ_cav в полости на их значения в свободном растворе (ε = 80 и Λ = 5 Å):

${\phi }_{1\mathrm{M}}\left(a\right)=+\frac{e}{a{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{e}{a}\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{\mathrm{\epsilon }}\right)f\left(2a/\mathrm{\Lambda }\right)$. (20.4)

Приведем теперь формулы для ΔW, выведенные в настоящей работе при использовании теории работы [38]. Для расчета энергии сольватации W иона радиуса а в неограниченном растворителе имеется известная формула, связывающая W с потенциалом на границе борновской сферы иона:

W(a) = (e/2) [e/a – φ(a)]. (21)

Для иона в полости потенциал зависит от ее радиуса согласно выведенной выше формуле (12), поэтому для этого случая энергия сольватации иона тоже зависит от радиуса полости, а формула (21) примет вид:

$W_{\mathrm{cav}}\left(R,a\right)=\left(e/2\right)\left[e/a-{\mathrm{\phi }}_{\mathrm{cav}}\left(R,a\right)\right]$. (22)

Подставляя в (22) потенциал на границе иона из (20.1), получим для энергии сольватации иона W_1Mcav в полости в рамках одномодовой модели диэлектрической функции формулу:

_{$W_{1M\left(\mathrm{cav}\right)}\left(R,a\right)=\frac{e^2}{2R}\left\{\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{p}}}\right)+\left(\frac{1}{{\mathrm{\epsilon }}_{\mathrm{o}}}-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{cav}}}\right)B_{\mathrm{cav}}\left(R,a\right)\right\}+w_{1\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right),$} (23)

где последний член w_1M(cav) в (23), не зависящий от радиуса полости, был введен выше формулой (16), а параметр B_cav определяется формулой (20.3). Формула для изменения энергии сольватации при переходе иона радиуса a из раствора в полость для одномодовой модели ε(k) имеет вид:

$\Delta W_{1\mathrm{M}\left(\mathrm{cav}\right)}=W_{1\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)-w_{1\mathrm{M}}\left(a\right)$, (24)

в которой последний член w_1M рассчитывается по известной формуле, полученной в работе [20] для энергии сольватации иона в неограниченном растворителе для одномодовой модели диэлектрической функции и был определен выше формулой (17).

Сопоставление результатов расчета изменения энергий сольватации иона в полости на основе диэлектрического приближения (DA) [23] и обратного диэлектрического приближения (IDA) [38] при использовании одномодовой модели диэлектрической функции. Сравнение расчетов для модели IDA (23)–(24) (сплошные кривые) и для модели DA (14)–(16) (кружки и прямоугольники) при использовании одной и той же одномодовой диэлектрической функции представлены на рис. 2. Красная сплошная кривая 2 и синие точки (прямоугольники) рассчитаны при ε_cav = 4.9, в то время как черная сплошная кривая 1 и зеленые точки (кружки) – для ε_cav = 80. Считалось, что диэлектрическая функция в полости определяется по формуле (18), а в неограниченном растворителе формулой (2.3), причем параметры ε_o и корреляционная длина как в свободном растворе (Λ), так и в полости (Λ_cav) были выбраны одинаковыми: ε_o = 4.9, Λ = Λ_cav = 5Å.

 

Рис 2. Изменение свободной энергии Гиббса (ΔG = –ΔW, в единицах kT) при переходе иона K⁺ в центр полости в зависимости от диэлектрической проницаемости белкового окружения водной полости. Все расчеты выполнены по нелокальной электростатике (NE) и одномодовой модели диэлектрической функции воды. Сплошные кривые 1 и 2 рассчитаны по NE-формулам (23)–(24), выведенным в соответствии с обратным диэлектрическим приближением (IDA) [38]. Кружки и прямоугольники – расчет по NE- формулам (14)-(16), полученным на основе подхода диэлектрического приближения (DA) [23]. Для кривой 2 и прямоугольников диэлектрическая проницаемость в полости равна ε_cav = 4.9, для кривой 1 и кружков ε_cav = 80. Для всех случаев параметр ε_o = 4.9, корреляционные длины в свободном растворе (Λ) и в полости (Λ_cav) равны друг другу: Λ = Λ_cav = 5 Å; радиус полости R = 5Å.

 

Результаты расчета по обоим подходам дают численно практически одни и те же величины для изменения энергии сольватации иона, поэтому можно сделать вывод о надежности предсказаний различных подходов, основанных как на теории [23], так и теории [38] к решению уравнений нелокальной электростатики для иона в водной полости.

Вывод формулы для изменения энергии сольватации ΔW иона при его переходе в центр сферической полости для трехмодовой модели диэлектрической функции растворителя в рамках обратного диэлектрического приближения (IDA)

Формула (23)–(24) для ΔW справедлива только для одномодовой модели диэлектрической функции, которая позволяет дать лишь полуколичественное описание энергии сольватациии. Для иона в бесконечном растворителе одномодовая модель использовалась в [20]. Однако для оценки изменения стандартной энергии Гиббса при переходе иона между двумя неограниченными растворителями в [42] приходилось применять трехмодовую модель. Поэтому для детального анализа процесса стабилизации катионов в полости канала следует проводить расчет изменения энергии сольватации для трехмодовой модели диэлектрической функции. Такой расчет выполнен в данном разделе в рамках обратного диэлектрического приближения [38–40].

Диэлектрическая функция в свободном растворе и в полости. Рассчитывается изменение энергии сольватации ΔW (a, следовательно, и изменение химического потенциала иона или его свободной энергии ΔG =–ΔW) при переходе иона радиуса а из неограниченной нелокальной среды с трехмодовой диэлектрической функцией (3M) ε_3M(k) в центр сферической водной полости (cav), окруженной локальным диэлектриком (диэлектрическая постоянная: ε_p), см. рис. 1а. Внутри полости нелокальная среда характеризуется трехмодовой диэлектрической функцией ε_cav(k). Параметры в полости и в наружном растворе могут различаться.

Приведем формулу ε(k) для трехмодовой модели диэлектрической функции [21, 41] в свободном растворе ε_3M:

$\frac{1}{{\epsilon }_{3\mathrm{M}}\left(k\right)}=1-\sum _{i=1}^3\frac{C_i}{1+{\left(k{\mathrm{\lambda }}_i\right)}^2}$, (25)

где константы C_i выражаются через диэлектрические проницаемости ε_i:

$C_1=1-1/{\mathrm{\epsilon }}_1,C_2=1/{\mathrm{\epsilon }}_1-1/{\mathrm{\epsilon }}_2,C_3=1/{\mathrm{\epsilon }}_2-1/{\mathrm{\epsilon }}_3$. (26)

Нижний индекс “1” относится к электронной моде, “2” – к инфракрасной моде, а “3” – к дебаевской (ориентационной моде). Значения диэлектрических проницаемостей соответствующих мод: ε₁ = 1.8, ε₂ = 4.9, ε₃ = 80, а значения корреляционных длин в свободном растворе: λ₁ = 0.5Å, λ₂ = 1Å, λ₃ = 5Å.

Соответствующая формула для ε_3M(cav)(k) для трехмодовой модели диэлектрической функции в водной полости имеет аналогичный вид:

_{$\frac{1}{{\mathrm{\epsilon }}_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(k\right)}=1-\sum _{i=1}^3\frac{C_{i\left(\mathrm{cav}\right)}}{1+{\left(k{\lambda }_{i\left(\mathrm{cav}\right)}\right)}^2}$,} (27)

где коэффициенты C_{i (cav)}, диэлектрические проницаемости ε_{i (cav)} и значения корреляционных длин λ_{i (cav)} в полости выбраны как:

_{$C_{1\left(\mathrm{cav}\right)}=C_1,C_{2\left(\mathrm{cav}\right)}=C_2,C_{3\left(\mathrm{cav}\right)}=1/{\mathrm{\epsilon }}_2-1/{\mathrm{\epsilon }}_{3\left(\mathrm{cav}\right)}.$} (28)

В соотношениях (28) диэлектрические проницаемости соответствующих мод определяются как: ε_1(cav) = ε₁ = 1.8, ε_2(cav) = ε₂ = 4.9, 4.9 < ε_3(cav) < 80, а корреляционные длины в формуле (27) как: λ_1(cav) = λ₁ = 0.5Å, λ_2(cav) = λ₂ = 1Å, для дебаевской моды корреляционная длина в полости может принимать значения в пределах 2Å < λ_3(cav) < λ₃ = 5Å. Таким образом, предполагается, что диэлектрические отклики в неограниченном растворе и внутри полости различаются только величинами длинноволновой диэлектрической постоянной (ε₃ или ε_3(cav)) и соответствующей корреляционной длины (λ₃ или λ_3(cav)).

Распределение потенциала вокруг иона в полости для трехмодовой модели диэлектрической функции (обратное диэлектрическое приближение [38]). Подставляя в (12) диэлектрическую функцию воды в полости (27) и вычисляя интеграл, получим для распределения потенциала вокруг иона в полости в рамках трехмодовой модели диэлектрической функции формулу:

${\phi }_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,r\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}R}+e\left(\frac{1}{r}-\frac{1}{R}\right)-e\sum _{i=1}^3{C_i}_{\left(\mathrm{cav}\right)}{A_i}_{\left(\mathrm{cav}\right)}\left(R,r\right),$ (29)

в которой три функции A_{i (cav)}(R, r) определяются как:

$A_{i\left(\mathrm{cav}\right)}\left(R,r\right)=\frac{g\left(r/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)},a/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)}\right)}{r}-\frac{g\left(R/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)},r/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)}\right)+g\left(R/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)},a/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)}\right)-f\left(2R/{\mathrm{\lambda }}_{\mathrm{i}\left(\mathrm{cav}\right)}\right)}{R}$. (30)

Вспомогательные функции g(x, y) и f (y) в (30) были определены формулами (19.3)–(19.4).

Потенциал на границе борновской сферы иона получается из формулы (29):

$\begin{array}{c}{\mathrm{\phi }}_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)=\frac{e}{{\mathrm{\epsilon }}_{\mathrm{p}}R}+e\left(\frac{1}{a}-\frac{1}{R}\right)-\\ -e\sum _{i=1}^3{C}_{i\left(\mathrm{cav}\right)}{\mathrm{\Phi }}_{i\left(\mathrm{cav}\right)}\left(R,a\right),\end{array}$ (31)

В которой три вспомогательныя функции Ф_{i (cav)}(R, a) определяются как,

${\mathrm{\Phi }}_{i\left(\mathrm{cav}\right)}\left(R,a\right)=\frac{f\left(2a/{\lambda }_{i\left(\mathrm{cav}\right)}\right)}{a}+\frac{B_{i\left(\mathrm{cav}\right)}\left(R,a\right)}{R}$, (32)

а три функции B_{i (cav)}(R, a) определяются через введенные формулами (19.3)–(19.4) функции g(x, y) и f (y):

${B_i}_{\left(\mathrm{cav}\right)}\left(R,a\right)=f\left(2R/{{\mathrm{\lambda }}_{\mathrm{i}}}_{\left(\mathrm{cav}\right)}\right)-2g\left(R/{{\mathrm{\lambda }}_i}_{\left(\mathrm{cav}\right)},a/{{\mathrm{\lambda }}_i}_{\left(\mathrm{cav}\right)}\right)$. (33)

Изменение энергии сольватации иона в полости ΔW для трехмодовой модели диэлектрической функции (обратное диэлектрическое приближение [38]). В данном разделе приводится вывод формулы для расчета изменения энергии сольватации иона ΔW (а, следовательно, и изменения свободной энергии Гиббса иона ΔG = –ΔW) при его переходе из раствора в центр сферической полости радиуса R, окруженной локальным диэлектриком (ε_p), см. рис. 1а. Использован новый подход в теории ограниченных в пространстве нелокальных диэлектрических сред, развитый в [38]. Диэлектрический отклик воды как в неограниченном растворе, так и внутри полости, описывается трехмодовой моделью, параметры которой в этих средах могут различаться: формула (25) для ε_3M(k) и формула (27) для ε_3Mсav(k).

Для расчета энергии сольватации иона радиуса а внутри полости воспользуемся формулой (22). Подставляя в нее формулу (31) для потенциала на границе борновской сферы, получим для энергии сольватации в рамках трехмодовой модели диэлектрической функции формулу:

$W_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)=\delta W_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)+w_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)$, (34.1)

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

$\begin{array}{c}\delta W_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)=\\ =\frac{e^2}{2R}\left\{\left(1-\frac{1}{{\mathrm{\epsilon }}_{\mathrm{p}}}\right)+\sum _{i=1}^3{C}_{i\left(\mathrm{cav}\right)}{B}_{i\left(\mathrm{cav}\right)}\left(R,a\right)\right\},\end{array}$ (34.2)

а последний член в (34.1) зависит от радиуса иона a, но не зависит от радиуса полости R:

$w_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)=\left(e^2/2a\right)\sum _{i=1}^3{C}_{i\left(\mathrm{cav}\right)}f\left(2a/{\lambda }_{i\left(\mathrm{cav}\right)}\right)$. (35)

Коэффициенты B_i(cav) в (34) были определены формулой (33). Корреляционная длина λ_3(cav) внутри полости может отличаться от ее значения в свободном растворе, поэтому в (34)–(35) λ₃ заменена на λ_3(cav). Энергия сольватации иона радиуса а в свободном растворе, согласно нелокально-электростатической теории [21, 41], рассчитывается как:

$w_{3\mathrm{M}}\left(a\right)=\left(e^2/2a\right)\sum _1^3{C}_{i}f\left(2a/{\mathrm{\lambda }}_i\right)$. (36)

Окончательная формула для расчета изменения энергии сольватации при переходе иона из раствора в полость для трехмодовой модели диэлектрической функции имеет вид:

$\Delta W_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)=\delta W_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(R,a\right)+\delta w_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)$, (37)

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

$\begin{array}{c}\delta w_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)=w_{3\mathrm{M}\left(\mathrm{cav}\right)}\left(a\right)-w_{3\mathrm{M}}\left(a\right)=\\ =\left(e^2/2a\right)\left[C_{3\left(\mathrm{cav}\right)}f\left(2a/{\mathrm{\lambda }}_{3\left(\mathrm{cav}\right)}\right)-C_{3}f\left(2a/{\mathrm{\lambda }}_3\right)\right].\end{array}$ (38)

Нелокально-электростатические расчеты изменения свободной энергии при переходе моновалентного иона в центр водной полости из раствора

Сравнение результатов расчетов для одномодовой модели NE и для локально-диэлектрической модели. На рис. 3 приведено сравнение изменения свободной энергии при переходе иона K⁺ в центр полости канала для одномодовой модели НЭ (корреляционная длина в полости: 5Å) и по формуле локальной электростатики (1).

 

Рис. 3. Изменение свободной энергии Гиббса (ΔG = –ΔW, в единицах kT) при переходе иона K⁺ в центр полости в зависимости от диэлектрической проницаемости белкового окружения водной полости. Сплошные кривые 1 и 2 – расчет по формулам (23)–(24) одномодовой модели (1M) нелокальной электростатики (NE), штриховые кривые 3 и 4 – расчет (ΔG = ΔG_LE) по формулам (1.1)–(1.2) локальной электростатики (LE). Для кривых 1 и 3 значение ε_cav = ε = 80, для кривых 2 и 4 ε_cav = ε_o = 4.9. Радиус полости R = 5Å, другие параметры приведены в табл. 1.

 

Таблица 1. Параметры, используемые при построении графиков на рис. 3. Для расчета по модели локальной электростатики (LE) ε_sol = ε.

Номер кривой на рис. 3	Свободный раствор или полость	LE или 1M	εo	ε или εcav	εp	Λ или Λ cav, Å
3 штриховая	Раствор	LE	—	80	—	—
	Полость	LE	—	80	1< εp<5	—
1 сплошная	Раствор	1M	4.9	80	—	5
	Полость	1M	4.9	80	1< εp<5	5
4 штриховая	Раствор	LE	—	80	—	—
	Полость	LE	—	4.9	1< εp<5	—
2 сплошная	Раствор	1M	4.9	80	—	5
	Полость	1M	4.9	4.9	1< εp<5	5

 

Из рис. 3 следует, что в рамках локальной электростатики (штриховые лини) изменение свободной энергии резко возрастает при уменьшении диэлектрической проницаемости в полости канала ε_cav от значения 80 (красная штриховая кривая 3) до 4.9 (черная штриховая 4), что заведомо делает невозможным объяснение стабилизации иона в полости.

Сравнение результатов расчетов (красные кривые) по локальной (штриховая линия 3) и одномодовой модели НЭ (сплошная линия 1) при диэлектрической проницаемости в полости канала ε_cav, равной 80, показывает небольшое снижение ΔG при переходе от локально-диэлектрического расчета к NE-расчету с корреляционной длиной λ₃ = 5Ả .

Сравнение результатов расчетов по одномодовым нелокально-электростатическим моделям (сплошные кривые 1 и 2) показывает небольшой рост изменение свободной энергии при переходе от диэлектрической проницаемости в полости канала ε_cav = 80 (красная сплошная кривая 1) к ε_cav = 4.9 (черная сплошная кривая 2) и неизменных остальных параметрах в полости.

Сравнение результатов расчетов по 1-модовой и 3-модовой моделям диэлектрической функции, а также для локальной модели. Сравнение НЭ-расчетов по формулам (23)–(24) для 1-модовой модели и по формулам (34.1)–(37) для 3-модовой моделям представлено на рис. 4. Расчеты проводились для двух значений корреляционной длины λ_3(cav) в полости.

 

Рис. 4. Нелокально-электростатические (NE) расчеты изменения свободной энергии Гиббса (ΔG = –ΔW, в единицах kT) при переходе иона K⁺ в центр полости в зависимости от диэлектрической проницаемости белкового окружения водной полости. Кривые 2 и 4 – расчет по формулам (23)–(24) одномодовой модели (1M) NE. Кривые 1 и 3 – расчет по формулам (34.1)–(37) трехмодовой модели (3M) NE. Для сплошных кривых 1 и 2 корреляционная длина воды в полости равна λ_3(cav) = 5Å, для штриховых 3 и 4 значение λ_3(cav) = 2Å. Радиус полости R = 5Å, другие параметры приведены в табл. 2.

 

Таблица 2. Параметры, используемые при построении графиков на рис. 4. Для одномодовой модели (1М) параметр ε_o указан в колонке ε₂, параметры ε и ε_cav в колонке ε₃, а корреляционные длины в свободном растворе (Λ) и в полости (Λ_cav) помещены в колонку λ₃

№ кривой на рис. 4	Раствор или полость	1М или 3М	ε1 или ε1(cav)	ε2 или ε2(cav)	ε3 или ε3(cav)	εp	λ1, Å или λ1(cav), Å	λ2, Å или λ2(cav), Å	λ3, Å или λ3(cav), Å
2 cплошная	Раствор	1М	—	4.9	80	—	—	—	5
	Полость	1М	—	4.9	80	1< εp<5	—	—	5
4 штриховая	Раствор	1М	—	4.9	80	—	—	—	5
	Полость	1М	—	4.9	80	1< εp<5	—	—	2
1 сплошная	Раствор	3М	1.8	4.9	80	—	0.5	1	5
	Полость	3М	1.8	4.9	80	1< εp<5	0.5	1	5
3 штриховая	Раствор	3М	1.8	4.9	80	—	0.5	1	5
	Полость	3М	1.8	4.9	80	1< εp<5	0.5	1	2

 

Из рис. 4 следует, что при ε_p = 5 и снижении корреляционной длины в полости канала наблюдается значительное уменьшение ΔG в полости, как это видно из сравнения штриховых (для λ_3(cav) = 2Å) и сплошных кривых (для λ_3(cav) = 5Å). При этом расчет по трехмодовой модели дает большие величины ΔG (которые намного превышают kT даже при ε_p = 5, см. рис. 5), чем расчет по одномодовой модели. Наименьшее значение ΔG достигается для расчета по одномодовой модели при наименьшей корреляционной длине в полости, равной 2Å (кривая 4).

 

Рис. 5. Изменение свободной энергии (ΔG = –ΔW, в единицах kT) при переходе иона K⁺ в центр полости в зависимости от корреляционной длины λ_3(cav) в водной полости. Кривые 1, 1' – расчет по формулам одномодовой модели (1M) нелокальной электростатики (NE). Кривые 2, 2' – расчет по формулам трехмодовой модели (3M) NE. Горизонтальные прямые линии 3, 3' – расчет по локальной электростатике (LE). Для сплошных кривых 1, 2, 3 радиус полости R = 5Å, для штриховых кривых 1', 2', 3' радиус полости R = 6Å, другие параметры приведены в табл. 3.

 

Таблица 3. Параметры, используемые при построении графиков на рис. 5. Для одномодовой модели (1М) параметр ε_o указан в колонке ε₂, параметры ε и ε_cav в колонке ε₃, а корреляционные длины в свободном растворе (Λ) и в полости (Λ_cav) помещены в колонку λ₃. Для локальной модели (LE) параметры ε_sol и ε_cav помещены в колонку ε₃

№ кривой на рис. 5	Раствор или полость	1М или 3М или LE	ε1 или ε1(cav)	ε2 или ε2(cav)	ε3 или ε3(cav)	εp	λ1, Å или λ1(cav), Å	λ2, Å или λ2(cav), Å	λ3, Å или λ3(cav), Å
1, 1'	Раствор	1М	—	4.9	80	—	—	—	5
	Полость 1: R=5Å 1': R=6Å	1М	—	4.9	80	εp = 5	—	—	2 <λ3< 5
2, 2'	Раствор	3М	1.8	4.9	80	—	0.5	1	5
	Полость 2: R=5Å 2': R=6Å	3М	1.8	4.9	80	εp = 5	0.5	1	2 <λ3< 5
3, 3'	Раствор	ЛЭ	—	—	80	—	—	—	—
	Полость 3: R=5Å 3': R=6Å	ЛЭ	—	—	80	εp = 5	—	—	—

 

Теперь приведем сравнение (рис. 5) результатов расчетов по формулам (23)–(24) одномодовой (ΔG_1М) и по формулам (34.1)–(37) трехмодовой (ΔG_3М) моделей при изменении корреляционной длины λ_3(cav) внутри полости со значением ΔG_LE, рассчитанным по формуле локальной электростатики (1.1)–(1.2) при двух значениях радиуса водной полости.

Для полости радиуса 5Å расчет по формуле локальной электростатики (1) дает ΔG_LE = 10.5 kT (кривая 3), а расчет по трехмодовой модели НЭ при λ_3(cav) = 2Å приводит к ΔG_3M = 5.1 kT (кривая 2), т.е. к уменьшению на 5.4 kT. Увеличение радиуса полости до 6Å при той же величине λ_3(cav) (2Å) приводит к понижению энергии ΔG_3M примерно еще на 2.2 kT, т.е. до 3.0 kT. Для одномодовой модели величина свободной энергии ΔG_1M при этих значениях параметров оказывается еще меньше: 0.86 kT.

Таким образом, НЭ-расчеты, представленные на рис. 5, показывают значительное уменьшение свободной энергии переноса иона ΔG при уменьшении корреляционной длины воды λ_3(cav) в полости канала до значения 2Å. Расчеты по НЭ демонстрируют большую чувствительность ΔG к изменению радиуса полости. Так при его увеличении всего на один ангстрем величина ΔG уменьшается на 2 kT при λ_3(cav) = 2Å для расчета по трехмодовой модели, как это видно из рис. 5. Эти два обстоятельства открывают возможность объяснить стабилизацию иона в полости, если ввести в рассмотрение дополнительное взаимодействие иона с зарядами белковых α-спиралей, направленных в центр полости, как это сделано в работе [2].

ЗАКЛЮЧЕНИЕ

Наше исследование посвящено проблеме преодоления катионом K⁺ энергетического барьера при входе этого иона внутрь KcsA калиевого канала из его водного раствора, обусловленного значительным понижением его энергии сольватации в канале W_cav за счет областей с низкой диэлектрической проницаемостью вокруг него. Рассмотрение этого явления в основополагающей работе Роукса и МакКиннона (Roux и MacKinnon) [2], проведенное в рамках описания воды внутри полости радиуса 5Å как полярной среды с той же локальной диэлектрической проницаемостью, что и в объеме раствора (ε_cav = ε_sol), дало для высоты барьера слишком большую величину, чтобы канал мог функционировать для ионного транспорта между раствором и внутриклеточным пространством. Для преодоления этого противоречия между результатом этого расчета и экспериментальными свидетельствами авторы работы [2] были вынуждены существенно усложнить модель системы, введя гипотезу о снижении энергетического барьера для входа иона за счет промотирующей роли двух других ионов К⁺ внутри канала.

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

Ввиду многочисленных указаний на возможность резкого понижения эффективной диэлектрической проницаемости воды внутри полостей или вблизи поверхностей твердых тел нами был, прежде всего, исследован вариант расчета ΔW в рамках той же локальной электростатической модели, но с меньшим значением диэлектрической постоянной внутри полости, т.е. при условии ε_cav < ε_sol. Однако такой расчет привел к противоположному эффекту, т.е. к резкому увеличению высоты барьера. Этот результат обусловлен тем, что уменьшение величины диэлектрической постоянной внутри полости ε_cav приводит к уменьшению энергии сольватации иона в полости.

Поэтому последующий поиск снижения энергетического барьера проводился прежде всего с учетом известного результата, что использованная в работе [2] локально-электростатическая модель М. Борна для энергии сольватации дает для нее существенно завышенный результат в объеме раствора по сравнению с ее экспериментальными значениями. Возможность успешной интерпретации этих экспериментальных данных для энергий сольватации дается нелокально-электростатической (NE) моделью описания диэлектрических свойств воды, которая была предложена Р.Р. Догонадзе и А.А. Корнышевым [41] и была успешно использована для описания как свойств ионов в растворах, так и электрических полей в пространственно-ограниченных полярных средах.

Основной трудностью теоретического анализа является использование NE-уравнений для расчета характеристик систем, где полярная среда (например, вода) занимает лишь пространственно-ограниченную область V, за пределами которой диэлектрические свойства системы резко отличаются. Такая проблема “вырезания объема” возникает уже при вычислении энергии сольватации отдельного иона, поскольку вода не проникает внутрь области, занятой ионом. Аналогично, при расчете электрических полей в системах, где вода занимает полость, необходимо учитывать ее отсутствие за пределами полости.

Для решения подобных задач нами был предложен в 2018 г. подход, основанный на обратном диэлектрическом приближении (IDA), согласно которому напряженность электрического поля E в любой точке среды с НЭ-свойствами r выражается через интеграл по области V, занятой средой, от произведения индукции D(r' ) и обратной нелокальной диэлектрической функции полярной среды ε^–1(r – r' ).

Важнейшим преимуществом подхода IDA является возможность для систем с высокой симметрией (в частности, сферически симметричных) получить решение для основных характеристик системы (напряженности и потенциала электрического поля, энергии сольватации иона) в аналитической форме – в виде интегралов для произвольного вида диэлектрической функции полярной среды и даже в виде комбинации известных функций для фурье-образа диэлектрической функции ε(k) дробно-рационального типа, например, для одно- и трехмодовой моделей. В широко используемом в литературе альтернативном подходе, основанном на “диэлектрическом приближении” (DA), который был предложен нами в 1977 году, для нахождения этих же характеристик необходимо решение интегрального уравнения, что в аналитическом виде удается получить решение только для одномодовой модели диэлектрической функции среды.

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

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

По сравнению с расчетом по формуле локальной электростатики (1) NE-расчеты показывают значительное уменьшение свободной энергии переноса иона ΔG при уменьшении корреляционной длины воды λ_3(cav) в полости канала до значения 2Å. Расчеты по НЭ демонстрируют большую чувствительность ΔG к изменению радиуса полости. Так при его увеличении всего на один ангстрем величина ΔG уменьшается на 2kT при λ_3(cav) = 2Å для расчета по трехмодовой модели. Эти два обстоятельства открывают возможность объяснить стабилизацию иона в полости, особенно если ввести в рассмотрение дополнительное взаимодействие иона с зарядами белковых α-спиралей, направленных в центр полости, как это сделано в работе [2].

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

Приведем пример, как полученные нами результаты можно будет использовать для оценки селективности транспорта моновалентных ионов при их прохождении полостей, рассмотренных в данной работе. В формуле для потока J ионов типа i должен появиться множитель u_iexp(–ΔG_i), где u_i – подвижность, а ΔG_i изменение свободной энергии перехода иона из неограниченного раствора в полость. Отношение таких множителей для моновалентных ионов равно селективности транспорта ионов разного типа. Нелокально-электростатические расчеты величин ΔG_i представлены в данной работе. Таким образом, селективность транспорта различных моновалентных ионов зависит от параметров системы, изученной в данной работе, – соотношения диэлектрических проницаемостей полости и ее окружения, а также радиуса полости.

ПРИЛОЖЕНИЕ 1

Методы нелокальной электростатики для пространственно ограниченных сред. Соотношения для подхода, основанного на обратном диэлектрическом приближении (IDA) [38]

Метод обратного диэлектрического приближения был разработан М.А. Воротынцевым (M.A. Vorotyntsev) с соавторами в [38–40]. В связи с тем, что он имеет решающее значение для получения результатов данной работы, целесообразно кратко обсудить его основные моменты.

При определенных условиях, установленных в работе [38] (см. также Приложение 2 в книге [21]), индукция D(r) тождественно совпадает с полем в вакууме G_vac(r), которое создается той же системой внешних зарядов ρ(r), а именно:

1. Для локального однородного диэлектрика: поверхности всех непроводящих полостей должны быть эквипотенциальными поверхностями поля G_vac. Тогда D(r) ≡ G_vac(r).
2. Для локального неоднородного диэлектрика: если ε(r) меняется на границе скачком, то должно выполняться условие для локального однородного диэлектрика (поверхности всех непроводящих полостей должны быть эквипотенциальными поверхностями поля G_vac). Кроме того, должно выполняться условие: G_vac(r) параллельно grad ε(r). Тогда D(r) ≡ G_vac(r).
3. Для диэлектрика с пространственной дисперсией: если ε_αβ(r₁, r₂) и границы раздела обладают той же симметрией, что и G_vac, то D(r) ≡ G_vac(r).

В частности, это условие выполняется в системах, где границы раздела сред и их диэлектрические функции обладают той же симметрией, что и поле в вакууме G_vac(r), определяемое распределением внешнего заряда. К таким системам относятся: 1) однородная и изотропная среда, занимающая все пространство; 2) плоская граница раздела z = 0, разделяющая среды либо с локальным откликом, либо обладающие диэлектрическими функциями, которые являются однородными и изотропными вдоль поверхности раздела, где плотность внешних зарядов зависит только от z; 3) сферически симметричная граница (или несколько сферически симметричных границ с тем же центром) r = a, разделяющая среды либо с локальным откликом, либо обладающие диэлектрическими функциями с той же сферической симметрией, где плотность внешних зарядов зависит только от r; 4) аналогичная система с цилиндрической симметрией. В случаях 2–4 при контакте двух нелокально-диэлектрических сред предполагается, что отсутствует корреляция флуктуаций поляризации в точках по разные стороны от межфазной границы, т.е. для таких значений пространственных координат диэлектрическая функция в формуле (3.2) основного текста (интегрирование в ней проводится по всему пространству):

$D_{\alpha }\left({\mathit{r}}_1\right)=\int {\mathrm{\epsilon }}_{\mathrm{\alpha \beta }}\left({\mathit{r}}_1,{\mathit{r}}_2\right)E_{\mathrm{\beta }}\left({\mathit{r}}_2\right)d{\mathit{r}}_2$

равна нулю, и это свойство автоматически выполняется и в следующих формулах (интегрирование в них тоже проводится по всему пространству):

$P_{\alpha }\left({\mathit{r}}_1\right)=\int {\chi }_{\mathrm{\alpha \beta }}^{}\left({\mathit{r}}_1,{\mathit{r}}_2\right)E_{\mathrm{\beta }}\left({\mathit{r}}_2\right)d{\mathit{r}}_2$,

$E_{\alpha }\left({\mathit{r}}_1\right)=\int {\mathrm{\epsilon }}_{\mathrm{\alpha \beta }}^{-1}\left({\mathit{r}}_1,{\mathit{r}}_2\right)D_{\mathrm{\beta }}\left({\mathit{r}}_2\right)d{\mathit{r}}_2$,

$P_{\alpha }\left({\mathit{r}}_1\right)=\int {\chi }_{\mathrm{\alpha \beta }}^{\left(D\right)}\left({\mathit{r}}_1,{\mathit{r}}_2\right)D_{\mathrm{\beta }}\left({\mathit{r}}_2\right)d{\mathit{r}}_2$,

где восприимчивость среды χ_αβ^(D)(r₁, r₂) определяется соотношением:

${\chi }_{\mathrm{\alpha \beta }}^{\left(D\right)}\left({\mathit{r}}_1-{\mathit{r}}_2\right)={\left(4\mathrm{\pi }\right)}^{-1}\left[{\delta }_{\mathrm{\alpha \beta }}\delta \left({\mathit{r}}_1-{\mathit{r}}_2\right)-{\mathrm{\epsilon }}_{\mathrm{\alpha \beta }}^{-1}\left({\mathit{r}}_1-{\mathit{r}}_2\right)\right]$. (A1)

Тогда соотношение (A.2) дает окончательный ответ для распределения электрического поля внутри пространственной области V, занятой средой с обратной диэлектрической функцией ε_αβ^–1(r₁, r₂):

$E_{\alpha }\left({\mathit{r}}_1\right)=\underset{\mathrm{V}}{\int }{\mathrm{\epsilon }}_{\mathrm{\alpha \beta }}^{-1}\left({\mathit{r}}_1,{\mathit{r}}_2\right){\left(G_{\mathrm{vac}}\right)}_{\mathrm{\beta }}\left({\mathit{r}}_2\right)d{\mathit{r}}_2$, (A2)

где поле в вакууме, естественно, зависит только от распределения внешних зарядов, но не зависит ни от формы границы раздела, ни от диэлектрических свойств среды. Потенциал поля φ(r) определяется согласно формуле (5), т.е. интегрированием выражения для электрического поля.

Аналогично, соотношение (A3) дает выражение для распределения поляризации среды внутри той же области V:

$P_{\alpha }\left({\mathit{r}}_1\right)=\underset{\mathrm{V}}{\int }{\chi }_{\mathrm{\alpha \beta }}^{\left(D\right)}\left({\mathit{r}}_1,{\mathit{r}}_2\right){\left(G_{\mathrm{vac}}\right)}_{\mathrm{\beta }}\left({\mathit{r}}_2\right)d{\mathit{r}}_2$. (A3)

Заметим, что для любого распределения G_vac(r), удовлетворяющего условию симметрии системы, существует соответствующее распределение плотности внешних зарядов. Поэтому формулу (A3) можно рассматривать как линейный отклик поляризации диэлектрической среды на произвольное внешнее поле. Тогда флуктуационно-диссипационная теорема выражает ядро χ_αβ^(D)(r₁, r₂) в формуле (A3) через коррелятор флуктуаций поляризации в отсутствие внешних полей <P_α(r₁, t₁) P_β(r₂, t₂)>, который отражает пространственно-временную структуру полярной среды [21, 43].

Если предположить, что эта структура среды остается неизменной вплоть до ее границ, то восприимчивость среды χ_αβ^(D)(r₁, r₂), однородной и изотропной вдали от ее границ, сохраняет этот вид во всей пространственной области, занятой средой, т.е. зависит только от расстояния между аргументами функции: χ_αβ^(D)(| r₁, r₂ |). Согласно соотношению (A1) такое свойство имеет и обратная диэлектрическая функция ε_αβ^–1(r₁– r₂), которая в этом случае связана с фурье-образом диэлектрической функции ε(k) [20, 21] (интегрирование во втором равенстве (A4) по всему k пространству):

$\begin{array}{c}{\mathrm{\epsilon }}_{\mathrm{\alpha \beta }}^{-1}\left({\mathit{r}}_1-{\mathit{r}}_2\right)\equiv {\delta }_{\mathrm{\alpha \beta }}{\mathrm{\epsilon }}_{}^{-1}\left(\left|{\mathit{r}}_1-{\mathit{r}}_2\right|\right)=\\ ={\delta }_{\mathrm{\alpha \beta }}{\left(2\mathrm{\pi }\right)}^{-3}\int {\left[\mathrm{\epsilon }\left(k\right)\right]}^1\exp \left[i\mathit{k}\left({\mathit{r}}_1-{\mathit{r}}_2\right)\right]d\mathit{k}=\\ ={\delta }_{\mathrm{\alpha \beta }}{\left(2{\mathrm{\pi }}^2\left|{\mathit{r}}_1-{\mathit{r}}_2\right|\right)}^{-1}\underset{0}{\overset{\infty }{\int }}\left[\mathrm{\epsilon }{\left(k\right)}^{-1}k\sin \left[k\left|{\mathit{r}}_1-{\mathit{r}}_2\right|\right]\right]dk.\end{array}$ (A4)

Таким образом, задание геометрии системы (в частности, пространственной области V, занятой полярной средой), фурье-образа ее диэлектрической функции в отсутствие границ ε(k) и распределение электрического поля внешних зарядов в отсутствие диэлектрических сред G(r) внутри этой пространственной области V позволяют найти распределение электрического поля в присутствии диэлектрических сред Е(r) по формулам (А2) и (A4).

Практически для систем, перечисленных выше (для которых индукция D(r) совпадает с полем в вакууме), решение сводится к однократным интегралам по соответствующей переменной. Поэтому это вычисление может быть легко проведено численно для любого функционального вида диэлектрической функции, т.е. ε(k). Эта особенность данного подхода дает ему огромное преимущество по сравнению с широко используемым альтернативным методом, основанным на “диэлектрическом приближении” (DA), в котором аналогичное предположение делается относительно диэлектрической функции среды ε_αβ(r₁, r₂) в соотношении (3), что приводит к необходимости решать интегральное уравнение для нахождения обратной диэлектрической функции в ограниченной области, где формула (А4) уже не применима.

ПРИЛОЖЕНИЕ 2

Список обозначений

Сокращения в тексте:

NE нелокальная электростатика

LE локальная электростатика

DA диэлектрическая аппроксимация

IDA обратная диэлектрическая аппроксимация

ΔW (ΔW < 0) изменение энергии сольватации иона при его переходе из свободного раствора в центр водной полости (ΔG = –ΔW: работа по переносу иона из раствора внутрь полости)

1M одномодовая модель диэлектрической функции воды

3M трехмодовая модель диэлектрической функции воды

Сокращения в индексах

1M одномодовая модель диэлектрической функции воды

3M трехмодовая модель диэлектрической функции воды

cav полость (от английского “cavity”)

DA диэлектрическая аппроксимация

LE локальная электростатика

Латинские

a радиус иона, рис. 1а

B_cav(R, a) функция радиуса полости и радиуса иона, определяемая формулой (20.3) и используемая в формуле (23)

B_{i (cav)}(R, a) функции радиуса полости и радиуса иона, определяемые формулой (33) при i = 1, 2, 3 и используемые в формуле (34)

C₁, C₂, C₃ константы в формуле (27) для трехмодовой модели диэлектрической функции воды, связанные по формулам (26) с диэлектрическими константами воды, соответствующими трем модам

D(r) индукция электрического поля в точке r (electric displacement). Нелокальное выражение для индукции через распределение электрического поля для среды, занимающей все пространство, дается формулой (3.2)

e заряд моновалентного иона, равный абсолютной величине заряда электрона

E(r) напряженность электрического поля в точке r. Для среды, занимающей все пространство, электрическое поле связано с индукцией интегральным соотношением (4)

f (y) вспомогательная функция одного аргумента, определяемая формулой (19.4)

g(x, y) вспомогательная функция двух аргументов, определяемая формулой (19.3)

ΔG = –ΔW изменение свободной энергии иона при его переходе из свободного раствора в центр сферической полости, формула (1.1)

G_vac(r) поле в вакууме, которое создается системой внешних зарядов ρ(r). Используется в Приложении 1, см. формулы (A2) и (A3)

j₁(x) сферические функции Бесселя; их определение, см. формулу (10.3)

k, k волновой вектор и его модуль

kT по оси абсцисс на рисунках отложено изменение свободной энергии, деленное на kT, где k – постоянная Больцмана, T – абсолютная температура

P(r) поляризация в точке r (математически выражается в виде нелокальной связи между величинами Р и Е: формула (3.1))

R радиус сферической водной полости, рис. 1а

r, r радиус-вектор, идущий из центра сферической водной полости в некоторую точку внутри полости; r = │r│

V область, заключенная между двумя концентрическими сферами радиусов a и R, см. рис. 1а

W_sol энергия сольватации иона в свободном растворе

W_cav энергия сольватации иона в центре сферической водной полости

ΔW изменение энергии сольватации иона при его переходе из свободного раствора в центр водной полости, формула (1.2)

W_1M(cav)(R, a) энергии сольватации иона радиуса а в центре сферической водной полости радиуса R в рамках одномодовой модели диэлектрической функции; формула (23), выведенная на основе метода обратного диэлектрического приближения работы [38]

w_1M(cav)(a) энергия сольватации иона радиуса a в неограниченном растворителе, но с параметрами такими же, как в полости малого размера, в рамках одномодовой модели диэлектрической функции; формула (16)

ΔW_1M(cav) изменение энергии сольватации при переходе иона из свободного раствора в полость для одномодовой модели диэлектрической функции; формула (24), выведенная на основе метода обратного диэлектрического приближения [38]

w_1M(a) энергия сольватации иона радиуса a в свободном растворе в рамках одномодовой модели диэлектрической функции; формула (17)

Δw_DA(R) изменение энергии сольватации при переходе иона бесконечно малых размеров из свободного раствора в полость с неизменными параметрами для одномодовой модели диэлектрической функции, выведенной в рамках диэлектрического приближения; формула (14)

ΔW_DA изменение энергии сольватации при переходе иона радиуса a из свободного раствора в полость для одномодовой модели диэлектрической функции, выведенной в рамках диэлектрического приближения; формула (15)

W_3M(cav)(R, a) энергия сольватации иона радиуса а в центре сферической водной полости радиуса R для трехмодовой модели диэлектрической функции; формула (34), выведенная на основе метода обратного диэлектрического приближения работы [38]

w_3M(сav)(a) энергия сольватации иона радиуса a в неограниченном растворителе, но с параметрами такими же, как в полости малого размера, в рамках трехмодовой модели диэлектрической функции; формула (35)

w_3M(a) энергия сольватации иона радиуса a в свободном растворе в рамках трехмодовой модели диэлектрической функции; формула (36)

ΔW_3M(cav)(R, a) изменение энергии сольватации при переходе иона из свободного раствора в центр сферической полости для трехмодовой модели ε(k); формула (37), выведенная на основе метода обратного диэлектрического приближения работы [38]

Греческие

ε диэлектрическая постоянная воды (ε = ε₃ = 80)

ε_o диэлектрическая постоянная воды в зоне прозрачности, разделяющей дебаевские и инфракрасные моды (ε_о = ε₂ = 4.9)

ε_sol диэлектрическая постоянная воды в объеме свободного раствора (ε_sol = ε = 80)

ε_cav диэлектрическая постоянная воды внутри полости (ε_cav = ε_3(cav))

ε_p диэлектрическая проницаемость среды, окружающей сферическую полость

ε₁ высокочастотная (“оптическая”) диэлектрическая проницаемость воды в свободном растворе (ε₁ = 1.8)

ε₂ диэлектрическая проницаемость воды, соответствующая колебательной моде в свободном растворе (ε₂ = ε_о = 4.9)

ε₃ диэлектрическая проницаемость воды, соответствующая дебаевской (ориентационной) моде в свободном растворе (ε₃ = ε = 80)

ε_1(cav) высокочастотная (“оптическая”) диэлектрическая проницаемость воды в полости (ε_1(cav) = ε₁ = 1.8)

ε_2(cav) диэлектрическая проницаемость воды, соответствующая колебательной моде в полости (ε_2(cav) = ε₂ = 4.9)

ε_3(cav) диэлектрическая проницаемость воды, соответствующая дебаевской (ориентационной) моде в полости 4.9 < ε_3(cav) < 80

ε_αβ(r₁, r₂) (тензорная и нелокальная) диэлектрическая функция системы; формула (3.3)

ε_αβ(cav)(r₁, r₂) (тензорная и нелокальная) диэлектрическая функция системы в полости (cм. формулу (6.1))

ε_αβ^–1(r₁, r₂) обратная (тензорная и нелокальная) диэлектрическая функция системы, связанная интегральным соотношением (3.4) с функцией ε_αβ(r₁, r₂)

ε_αβ^–1_(cav)(r₁ – r₂) обратная (тензорная и нелокальная) диэлектрическая функция в полости, связанная с фурье-образом диэлектрической функции воды внутри полости ε_cav(k) (формула (8))

ε(k) диэлектрическая функция неограниченной нелокальной среды

ε_cav(k) диэлектрическая функция водной фазы внутри полости

ε_1M(k) диэлектрическая функция неограниченной нелокальной среды для одномодовой модели; формулы (2.2)–(2.3)

ε_1M(cav)(k) диэлектрическая функция водной фазы внутри полости для одномодовой модели; формула (18)

ε_3M(k) диэлектрическая функция неограниченной нелокальной среды для трехмодовой модели; формула (25)

ε_3M(cav)(k) диэлектрическая функция водной фазы внутри полости для трехмодовой модели; формула (27)

ε^–1_l=1(r₁, r₂) коэффициент разложения обратной диэлектрической функции ε_αβ^–1(r₁, r₂) по сферическим функциям Y_l0 для нижнего индекса l = 1

λ₁ корреляционная длина воды, соответствующая электронной моде в свободном растворе (λ₁ = 0.5Å)

λ₂ корреляционная длина воды, соответствующая колебательной моде в свободном растворе (λ₂ = 1Å)

λ₃ корреляционная длина воды, соответствующая дебаевской (ориентационной) моде в свободном растворе (λ₃ = 5Å)

λ_1(cav) корреляционная длина воды, соответствующая электронной моде в полости (λ_1(cav) = 0.5Å)

λ_2(cav) корреляционная длина воды, соответствующая колебательной моде в полости (λ_2(cav) = 1Å)

λ_3(cav) корреляционная длина, соответствующая дебаевской (ориентационной) моде в полости (2Å < λ_3(cav) < 5Å)

Λ максимальная из трех корреляционных длин воды в свободном растворе (Λ = λ₃ = 5Å)

Λ_cav максимальная из трех корреляционных длин воды в полости (Λ_cav = λ_3(cav))

ρ(r) плотность внешних зарядов

φ(r) распределение потенциала электрического поля

φ_1M(cav)(R, r) распределение потенциала вокруг иона в полости в рамках одномодовой модели диэлектрической функции, рассчитываемое по формуле (19.1) при a < r < R

φ_1M(cav)(R, a) потенциал на границе борновской сферы иона для одномодовой модели диэлектрической функции, рассчитываемый по формуле (20.1)

φ_3M(cav)(R, r) распределение потенциала вокруг иона в полости в рамках трехмодовой модели диэлектрической функции, рассчитываемое по формуле (29) при a < r < R

φ_3M(cav)(R, a) потенциал на границе борновской сферы иона для трехмодовой модели диэлектрической функции, рассчитываемый по формуле (31)

χ_αβ(r₁, r₂) тензор нелокальной (диэлектрической) восприимчивости

ФИНАНСИРОВАНИЕ РАБОТЫ

Номер госзадания 124022600337-0 (Рубашкин А.А., Вигонт В.А.); номер госзадания 122011300058 (Воротынцев М.А.).

КОНФЛИКТ ИНТЕРЕСОВ

Авторы заявляют, что у них нет конфликта интересов.

About the authors

A. A. Rubashkin

Institute of Cytology RAS

Author for correspondence.
Email: andrey.rubashkin@gmail.com
Russian Federation, St. Petersburg

V. A. Vigont

Institute of Cytology RAS

Email: andrey.rubashkin@gmail.com
Russian Federation, St. Petersburg

M. A. Vorotyntsev

Institute of Physical Chemistry and Electrochemistry named after A.N. Frumkin RAS

Email: mivo2010@yandex.com
Russian Federation, Moscow

References

Supplementary files