Письма в ЖЭТФ, том 112, вып. 9, с. 632 - 636
© 2020 г. 10 ноября
Комбинированная схема восстановления функции
распределения частиц по размерам с использованием данных
малоуглового рассеяния
В.В.Волков1), П.В.Конарев, А.Е.Крюкова
Институт кристаллографии им. А. В. Шубникова,
Федеральный научно-исследовательский центр “Кристаллография и фотоника” РАН, 119333 Москва, Россия
Поступила в редакцию 7 октября 2020 г.
После переработки 7 октября 2020 г.
Принята к публикации 7 октября 2020 г.
Для анализа полидисперсных систем в наноразмерном диапазоне существует ряд алгоритмов, поз-
воляющих определять профиль распределения частиц по размерам из данных малоуглового рентгенов-
ского и нейтронного рассеяния. Соответствующая задача наименьших квадратов плохо обусловлена,
что может приводить к сильной зависимости решения от стартового приближения и параметров ал-
горитмов поиска. Процедура, основанная на решении системы линейных уравнений с регуляризацией,
дает стабильное решение, но не свободное от артефактов. Решения, предлагаемые другими алгоритма-
ми, которые ищут функцию распределения как в аналитическом виде, так и в виде непараметрической
гистограммы, могут зависеть от параметров поиска и начального приближения. В работе предложена
комбинированная схема использования этих алгоритмов, позволяющая повысить устойчивость решения
и, как следствие, его надежность.
DOI: 10.31857/S1234567820210107
Введение. Свойства многих современных ма-
(2) прямой поиск распределения частиц по разме-
териалов определяются структурными особенностя-
рам в виде гистограммы методами нелинейных наи-
ми в наноразмерном диапазоне, поэтому для их ис-
меньших квадратов (программа VOLDIS, неопубл.)
следования часто применяют методы малоуглового
и случайного поиска McSAS [4];
рентгеновского и нейтронного рассеяния (МУРР и
(3) постулирование вида распределения в анали-
МУРН, соответственно), которые позволяют каче-
тическом виде (например, нормального распреде-
ственно и количественно исследовать строение ве-
ления или распределения Шульца [5]) и проведе-
щества без специальной подготовки образцов в диа-
ние приближения данных методами нелинейных наи-
пазоне разрешений от долей нанометра до несколь-
меньших квадратов (программа MIXTURE (и ее мо-
ких микрон. Анализ таких систем представляет со-
дифицированный вариант POLYMIX) [6] и програм-
бой плохо обусловленную задачу, что может приво-
ма SASFIT [7]).
дить к большим разбросам решений. Так, число обу-
Каждый из этих алгоритмов обладает своими
словленности матрицы вторых производных мини-
преимуществами и недостатками, и использование
мизируемой целевой функции может достигать де-
только одного из них не всегда позволяет получать
сятков и сотен миллионов, что приводит к сильной
решения, свободные от артефактов. В данной работе
зависимости решения не только от ошибок в исход-
нами предложена новая схема совместного исполь-
ных данных, но и от реализации алгоритма поиска, и
зования алгоритмов, позволяющая повысить надеж-
от величин стартовых значений параметров модели.
ность восстановления функции распределения час-
Для анализа полидисперсных систем по малоуг-
тиц по размерам. Работа схемы проиллюстрирова-
ловым данным существует ряд алгоритмов, среди ко-
на на примерах двухкомпонентных и трехкомпонент-
торых можно выделить несколько основных:
ных систем полидисперсных сферических частиц с
(1) методы линейных наименьших квадратов, в
разнесенными и частично перекрывающимися рас-
том числе с регуляризацией решения по Тихонову
пределениями компонент.
(программa GNOM [1] из пакета ATSAS [2] и про-
Алгоритмы восстановления функции рас-
грамма GIFT [3]);
пределения частиц по размерам по данным
малоуглового рассеяния. Для вычисления функ-
1)e-mail: volkicras@mail.ru
ции распределения по размерам частиц в полидис-
632
Письма в ЖЭТФ том 112 вып. 9 - 10
2020
Комбинированная схема восстановления функции распределения частиц по размерам...
633
персных системах по данным интенсивности рассея-
I(s) =
vkIk(s),
(2)
ния I(s), программы GNOM [1] и GIFT [3] решают
k=1
интегральное уравнение
где vk - относительная объемная доля k-ого компо-
нента с интенсивностью рассеяния Ik(s), которая, в
свою очередь, зависит от вида (объемной) функции
I(s) =
DV (r)I0(s, r)dr,
(1)
(r) и заданного квадрата форм-
распределения DV
k
Rmin
фактора частицы I0
k
(s, r):
где DV (r) - искомая функция распределения, опре-
деленная в диапазоне размеров [Rmin, Rmax], I0(s, r)
Ik(s) =
DV
k
(r)I0
k
(s, r)dr,
(3)
- заданный форм-фактор частицы (в данной работе
R
min
мы использовали сферы), s = (4π/λ) sin θ - модуль
где r - размер частицы. Тогда результирующее рас-
вектор рассеяния, 2θ - угол рассеяния, λ - длина вол-
пределение есть линейная комбинация парциальных
ны. Данное интегральное уравнение, представлен-
функций:
ное в конечных разностях в виде системы линейных
уравнений с матрицей коэффициентов, составленной
DV (r) = vkDV
(r).
(4)
k
из интенсивностей рассеяния сферическими частица-
k=1
ми, рассчитанными на экспериментальной сетке уг-
Решение задачи состоит в определении параметров
ловых отсчетов интенсивности и искусственной сет-
аналитических распределений и их отностительных
ке радиусов от Rmin, (обычно от 1Å) до предполага-
вкладов в сумму (2) путем минимизации квадра-
емого максимального значения Rmax, представляет
тичного отклонения теоретической интенсивности от
собой плохо обусловленную задачу. Для получения
экспериментальной.
устойчивых решений, имеющих физический смысл,
В программе MIXTURE для поиска решения ис-
в программе GNOM применена регуляризация реше-
пользуется метод переменной метрики в вариан-
ния по Тихонову [8], т.е. минимизируется функцио-
те Бройдена-Флетчера-Голдфарба-Шанно [9] с про-
нал, в который входит как квадрат невязки между
стыми ограничениями на параметры, в програм-
экспериментальной и рассчитанной кривой рассея-
мах POLYMIX и VOLDIS реализован модифициро-
ния, так и квадрат нормы первой производной функ-
ванный вариант алгоритма Левенберга-Марквардта
ции распределения, умноженный на параметр регу-
[10], тогда как в программе McSAS применен отно-
ляризации и обеспечивающий гладкость контура [1].
сительно медленно сходящийся метод Монте-Карло.
Решение чувствительно к эффектам обрыва I(s) и
Использование стохастических методов, хотя и обла-
может содержать осцилляции, существенно искажа-
дающих свойствами поиска глобального экстремума,
ющие форму распределения. Параметр регуляриза-
в данном случае нам представляется необоснован-
ции в данном подходе фактически представляет со-
ным, так как нестабильность поиска непараметри-
бой степень сглаживания решения, позволяющий по-
ческой гистограммы распределения вызвана плохой
высить его устойчивость, но требует подбора области
обусловленностью задачи и ограниченностью точно-
определения [Rmin, Rmax] функции DV (r).
сти машинных вычислений, а не ее неоднозначно-
Два других алгоритма используют методы
стью.
нелинейной минимизации, причем в одном случае
Сравнение эффективности восстановления
функцию распределения непосредственно ищут
функции распределения для двух и трех-
в виде гистограммы (программы McSAS
[4] и
компонентных систем полидисперсных сфе-
VOLDIS), в другом случае функция распределения
рических частиц. В данной работе мы рассмотре-
задается в аналитическом виде (программа SASFIT
ли случаи двухкомпонентных и трехкомпонентных
[7]) или комбинации аналитических выражений
систем полидисперсных сферических частиц с разне-
(MIXTURE/POLYMIX [6]).
сенными и частично перекрывающимися распределе-
В последнем случае интенсивность рассеяния
ниями компонентов. Для этих систем были смодели-
представляют в виде линейной комбинации вкладов
рованы наборы данных малоуглового рассеяния как
рассеяния от K компонент, при этом каждая компо-
без шума, так и с добавлением шума, подчиняюще-
нента есть система частиц заданной формы, распре-
гося Пуассоновскому распределению. Структурные
деление которых описывается одним аналитическим
параметры систем приведены в табл. 1.
выражением с учетом, при необходимости, межча-
По незашумленным модельным данным все
стичной интерференции:
три алгоритма (программы GNOM, VOLDIS и
Письма в ЖЭТФ том 112 вып. 9 - 10
2020
634
В.В.Волков, П.В.Конарев, А.Е.Крюкова
Таблица 1. Структурные параметры двух- и трехкомпонентных систем полидисперсных сферических частиц. Ri - средний
радиус сфер i-компоненты, dRi - ширина распределения i-компонента, vi - объемная доля i-компоненты (i = 1 ÷ 3)
Номер
Компонент 1
Компонент 2
Компонент 3
системы
R1 = 5.0 нм
R2 = 12.0 нм
1
dR1 = 1.5 нм
dR2 = 3.0 нм
v1 = 0.70
v2 = 0.30
R1 = 6.0 нм
R2 = 11.0 нм
2
dR1 = 4.0 нм
dR2 = 2.0 нм
v1 = 0.10
v2 = 0.90
R1 = 4.0 нм
R2 = 15 нм
3
dR1 = 0.7 нм
dR2 = 1.5 нм
v1 = 0.20
v2 = 0.80
R1 = 4.0 нм
R2 = 8.0 нм
R3 = 15.0 нм
4
dR1 = 0.5 нм
dR2 = 1.5 нм
dR3 = 3.0 нм
v1 = 0.70
v2 = 0.20
v3 = 0.10
R1 = 3.0 нм
R2 = 8.0 нм
R3 = 22.0 нм
5
dR1 = 0.4 нм
dR2 = 2.0 нм
dR3 = 3.0 нм
v1 = 0.005
v2 = 0.078
v3 = 0.917
R1 = 37.0 нм
R2 = 65.0 нм
R3 = 130.0 нм
6
dR1 = 7.0 нм
dR2 = 4.0 нм
dR3 = 5.5 нм
v1 = 0.01
v2 = 0.89
v3 = 0.10
MIXTURE/POLYMIX) восстанавливают решение
точками) также содержат отличия от точного реше-
с удовлетворительной для практики точностью
ния (сплошные линии), причем особенно неустойчи-
(относительная ошибка не превышает 0.1-0.2 %, на
выми также оказываются области малых размеров
графиках отклонение не видно, и поэтому рисунки
(особенно хорошо это видно на примере рис. 2(5)).
мы не приводим). Естественно, по зашумленным
Этот эффект связан с тем, что парциальные кривые
кривым рассеяния, показанным на рис. 1-2, удается
рассеяния от частиц близких размеров очень похожи,
найти только приближенное распределение. При
и программа, увеличивая вклад частиц какого-либо
этом форма функции распределения, восстанавлива-
радиуса, компенсирует чрезмерный рост интенсивно-
емая программой GNOM (синие пунктирные кривые
сти рассеяния уменьшением вклада от соседних час-
на рисунках), оказываeтся наиболее отличающейся
тиц с близкими размерами. При малых же размерах
от точного решения (сплошные красные линии на
сходство форм кривых рассеяния (например, по кри-
рисунках). Если положения пиков распределений
терию корреляции) намного больше, чем в области
определяются достаточно точно (отклонения не
больших размеров. Это связано с тем, что рассея-
превышают 5 %), то ширины распределений оказы-
ние от малоразмерных частиц в экспериментальной
ваются намного больше заданных при построении
области углов представляет собой монотонные мед-
моделей. Кроме того, найденные распределения
ленно спадающие функции, тогда как при больших
имеют артефакты в виде осцилляций из-за эффекта
размерах парциальные кривые содержат особенно-
обрыва функции интенсивности рассеяния. Регуля-
сти в виде максимумов (по крайней мере, на про-
ризация, используемая в программе GNOM, хотя
изводных), имеющих разные угловые положения.
и позволяет находить достаточно гладкие решения
Программа POLYMIX (аналогично MIXTURE)
(требование гладкости и приводит к уширению
позволяет находить решения, практически совпада-
пиков) с малой амплитудой артефактов, тем не
ющие с точными во всех представленных случаях
менее может неверно передавать форму распреде-
(сплошные красные кривые на рисунках, визуаль-
ления в области начальных (малых) размеров из-за
но не отличимые от точных контуров) при старте
чрезмерного ограничения первой производной по
с приближений, близких к решениям, полученным
контуру.
программами GNOM и VOLDIS. Однако, решения
Результаты программы непараметрического по-
могут содержать и значительные отличия, если стар-
иска гистограммы распределения VOLDIS для за-
товое приближение далеко от истинного. Как было
шумленных данных (кривые на рисунках показаны
показано ранее в работах [11, 12], диапазон области
Письма в ЖЭТФ том 112 вып. 9 - 10
2020
Комбинированная схема восстановления функции распределения частиц по размерам...
635
Рис. 1. (Цветной онлайн) Столбец слева: малоугловые
Рис. 2. (Цветной онлайн) Данные от трехкомпонентной
данные от двухкомпонентной системы полидисперсных
системы полидисперсных сферических частиц, номера
сферических частиц, номера строк соответствуют мо-
строк соответствуют моделям 4-6 из табл. 1. Обозначе-
делям 1-3 из табл. 1. Строки: крестики - зашумленные
ния и цвета аналогичны рис. 1
модельные данные, синяя пунктирная кривая - реше-
ние GNOM, черные точки - решение VOLDIS; сшлош-
пользовать в качестве стартовых параметры, оценен-
ная красная кривая - решение POLYMIX, полученное
с использованием стартовых параметров, оцененных
ные по решениям GNOM и MIXTURE/POLYMIX.
по решениям GNOM и VOLDIS. Рассеяние от точной
Схема комбинированного использования
модели визуально совпадает с решением POLYMIX.
алгоритмов восстановления функции распре-
Столбец справа: контуры полученных распределений
деления. Основываясь на сравнении полученных
по радиусам частиц, обозначения кривых те же
результатов восстановления функции распределения
частиц по размерам разными алгоритмами, мы
предлагаем следующую схему их комбинированного
стартовых значений параметров распределений, с ко-
использования:
торых возможно восстановление правильного реше-
1) использовать программу VOLDIS для оценки
ния, зависит от многих факторов, в частности, от от-
гистограммы распределения, числа компонентов (т.е.
носительного вклада параметров компонентов смеси
вкладов, представленных унимодальными распреде-
в общую кривую рассеяния и от уровня и вида шума
лениями) и среднего радиуса частиц в каждом ком-
в данных. Поэтому естественной кажется мысль ис-
поненте. При этом особое внимание следует обратить
Письма в ЖЭТФ том 112 вып. 9 - 10
2020
636
В.В.Волков, П.В.Конарев, А.Е.Крюкова
на область углов рассеяния от 0 до начального экс-
при частичной финансовой поддержке Российско-
периментального. Теоретическая кривая интенсивно-
го фонда фундаментальных исследований (грант
сти должна монотонно экстраполировать экспери-
#19-32-90190).
ментальный контур, в противном случае найденное
решение будет содержать артефакты в виде значи-
1. D. I. Svergun, J. Appl. Crystallogr. 25, 495 (1992).
тельных вкладов больших размеров частиц. Для то-
2. D. Franke, M. V. Petoukhov, P. V. Konarev,
го чтобы этого избежать, необходимо провести по-
A. Pankovich, A. Tuukkanen, H.D. T. Mertens,
вторный поиск, искусственно уменьшая Rmax;
A. G. Kikhney, N.R. Hajizadeh, J. M. Franklin,
2) использовать программу GNOM для нахож-
C. M. Jeffries, and D. I. Svergun, J. Appl. Crystallogr.
50, 1212 (2017).
дения предварительной оценки функции распреде-
3. O. Glatter, J. Appl. Crystallogr. 12, 166 (1979).
ления, указывая в ней максимальный размер час-
4. I. Bressler, B. R. Pauw, and A. F. Thuenemann, J. Appl.
тиц Rmax, подобранный на предыдущем шаге. Цель
Crystallogr. 48, 962 (2015).
данного численного эксперимента является получе-
5. G. V. Schulz, Z. Phys. Chem. Abt. B 30, 379 (1935).
ние гладкого контура, уточнения числа и параметров
6. P. V. Konarev, V. V. Volkov, A. V. Sokolova,
компонентов;
M. H. J. Koch, and D. I. Svergun, J. Appl. Crystallogr.
3) уточнить параметры распределения с помо-
36, 1277 (2003).
щью программ POLYMIX/MIXTURE, используя в
7. I. Bressler, J. Kohlbrecher, and A. F. Thuenemann,
качестве стартового приближения параметры рас-
J. Appl. Crystallogr. 48, 1587 (2015).
пределений, оцененные на предыдущих шагах.
8. A. N. Tikhonov and V. Ya. Arsenin, Solution of ill Posed
Предложенная схема позволяет улучшить устой-
Problems, Wiley, N.Y. (1977), 258 p.
чивость восстановления функции распределения
9. P. E. Gill, W. Murray, and M. H. Wright, Practical
частиц по размерам с использованием малоугловых
Optimization, Academic Press, London (1981), 401 p.
данных. На третьем шаге схемы можно дополни-
10. J. J. More, The Levenberg-Marquardt Algorithm,
тельно проводить комбинированиe локальных и
Implementation and Theory. Lecture Notes in
глобальных минимизационных методов, что позво-
Mathematics, ed. by G. A. Watson, Springer-Verlag,
лит еще больше расширить диапазон сходимости
Berlin (1978), v. 630, p. 105.
к правильному решению и, таким образом, повы-
11. A. E. Kryukova, P. V. Konarev, and V. V. Volkov,
сить достоверность определения распределений по
Crystallogr. Rep. 63, 26 (2018).
размерам, как было показано в работе [13].
12. A. E. Kryukova, A. S. Kozlova, P. V. Konarev,
Работа выполнена при поддержке Министерства
V. V. Volkov, and V. E. Asadchikov, Crystallogr.
науки и высшего образования в рамках выполне-
Rep. 63, 531 (2018).
ния работ по Государственному заданию ФНИЦ
13. A. E. Kryukova, P. V. Konarev, V. V. Volkov, and
“Кристаллография и фотоника” РАН, а также
V. E. Asadchikov, J. Mol. Liq. 283, 221 (2019).
Письма в ЖЭТФ том 112 вып. 9 - 10
2020