CC BY
11
УДК 543.51
П.А. Рюмин, И.А. Болдин, Д.М. Автономов, Е.Н. Николаев Институт энергетических проблем химической физики РАН
Метод емкостей в моделировании движения ионных ансамблей в ионных ловушках и системах транспорта ионов с электродами произвольной формы
Для вычисления электрического поля, действующего на индивидуальную заряженную частицу со стороны других ионов, электродов и зарядов изображений в ионных ловушках и устройствах транспорта ионов с электродами произвольной формы в масс-спектрометрических устройствах используется метод ёмкостей. Метод заключается в замене реальных электродов, на поверхностях которых наводятся заряды изображения зарядами внутри области, ёмкостями бесконечно малого размера, на которые помещаются заряды таким образом, чтобы сделать поверхности электродов эквипотенциальными. Было промоделировано движение облаков ионов внутри ловушки Кингдона (Орбитрэп) с использованием метода ёмкостей в рамках метода частицы в ячейке при количестве емкостей на поверхностях электродов 1600.
Ключевые слова: компьютерное моделирование, движение ионов, метод емкостей, орбитрэп, заряд изображения.
I. Введение
Одним из важных параметров масс-спектрометра, используемого при решении аналитических задач, является динамический диапазон концентраций анализируемых веществ и количество ионов, с которыми он может работать, не теряя своих характеристик, таких, как точность измерения масс. Динамический диапазон 106 означает, что прибор в состоянии разрешить и измерить массы единичных ионов и миллиона ионов, созданных одновременно в ионном источнике. При количестве ионов 106 на их движении сказывается поле от зарядов других ионов, и в некоторых устройствах, таких как ионные ловушки, его влияние становится доминирующим. При моделировании движения ионов в масс-спектрометре, обладающим большим динамическим диапазоном, необходимо учитывать влияние ион-ионного взаимодействия не только непосредственно, но и через заряды, наведенные в электродах различных ионно-оптических устройств масс-спектрометра. Наиболее широко используемая в масс-спектрометрическом сообществе для расчета электрических полей и траекторий движения ионов в них программа Simion [1] не предоставляет адекватного инструментария для такого типа анализа.
Ион-ионное взаимодействие может быть учтено путем расчета электрического поля, создаваемого каждым ионом в точках, где находятся все остальные ионы. Этот метод называется Particle-Particle (частица-частица). При этом взаимодействие с зарядами-изображениями вообще не учитывается, то есть метод целесообразно использовать, только если ионы не подходят близко к стенкам электродов. Другой недостаток — быстрый рост объемов вычислений с ростом числа частиц. Метод, позволяющий избежать этой проблемы и работать с большим количеством частиц, — это метод Particle-In-Cell (частица в ячейке).
Метод Particle-In-Cell (PIC) известен с 1955 г. [2] и широко применяется в физике плазмы [3]. Идея метода состоит в том, чтобы не рассчитывать ион-ионные взаимодействия по одному (как в методе Particle-Particle), а интерполировать заряды в узлы вычислительной сетки. Тогда количество зарядов получается равным числу узлов сетки. Затем уравнение Пуассона решается на этой сетке, для чего могут быть использованы различные алгоритмы, но самым быстрым при использовании параллельных вычислений является метод быстрого преобразования Фурье Fast Fourier Transform (FFT). Кроме того, FFT является прямым методом, то есть после выполнения определённого числа шагов он дает точный результат, в отличие от итерационных, где каждый шаг алгоритма уменьшает ошибку в определенное число раз. Проблема метода FFT состоит в том, что он применим только к устройствам прямоугольной геометрии, тогда как реальные
масс-спектрометрические устройства могут иметь электроды любой формы. Чтобы преодолеть это ограничение метода FFT, нами был применен алгоритм, использующий метод емкостной матрицы.
Рис. 2. Конфигурация электродов помещается в Рис. 1. Ловушка Кингдона (Orbitrap) в разрезе соответствующую прямоугольную вычислитель-
ную сетку, и электроды заменяются набором точек. На рисунке, для наглядности, заменен только внутренний электрод
II. Алгоритм вычислений
Заряженные частицы, находящиеся в окружении электродов, индуцируют заряды на поверхности этих электродов, в результате чего потенциал оказывается постоянным вдоль электродов. Этот процесс можно моделировать, распределяя виртуальные заряды рядом с поверхностями электродов таким образом, чтобы потенциал в точках каждого электрода оказался равным наперед заданному значению. Один из методов нахождения нужного распределения зарядов — это метод емкостной матрицы [4, 5]. Емкостная матрица (C) — это матрица в системе линейных уравнений Cp =q, где р — столбец, содержащий значения потенциала в определенном наборе точек, а q — величины зарядов, находящихся в этих точках (в действительности для лучшей устойчивости решения точки, в которых размещены заряды, немного смещены от точек, в которых измеряется потенциал). Рассмотрим метод емкостной матрицы на примере ловушки Кингдона (Orbitrap) [6] (рис. 1). Сначала нужно поместить поверхности электродов в прямоугольную вычислительную сетку (на этой сетке будут решаться уравнения Пуассона методом FFT) и заменить электроды емкостями (рис. 2). Каждый элемент поверхности (емкость) моделируется двумя близко лежащими точками — в одной будет вычисляться потенциал, в другую помещаются заряды. Для того чтобы найти емкостную матрицу, нужно разместить в соответствующую точку каждого элемента поверхности единичный заряд и вычислить методом FFT потенциалы, наведенные этим зарядом на каждую емкость. Значения потенциалов в соответствующих точках образуют столбец матри-
цы A (рис. 3). | A A есть обратная матрица C:
" qi " p1
qn pn
. Как видно из дальнейшего рассуждения, матрица
A
-і
A
" qi " ( \ р1 q1 1 \ Р1
qn = ( а-і ) pn qn =(C) pn
. (1)
Таким образом, после нахождения матрицы A надо ее обратить, получив матрицу C. После этого матрица C сохраняется в памяти компьютера и используется при расчете электрического поля в основной программе. Матрицу С необходимо вычислить один раз для каждой конфигурации электродов.
Полный алгоритм расчета движения ионов на каждом временном шаге выглядит так (шаги, помеченные звездочкой, подробно описаны в [7]):
1) интерполяция зарядов ионов в узлы вычислительной сетки*;
2) решение уравнения Пуассона методом FFT*;
3) экстраполяция значений потенциалов из узлов сетки в точки емкостей (алгоритм идентичен алгоритму интерполяции поля в точки, где находятся ионы) *;
4) определение столбца зарядов емкостей путем умножения матрицы C на столбец потенциалов, представляющий собой разность изначально заданного потенциала электрода, и потенциала, полученного в пункте 3;
5) интерполяция зарядов ионов и зарядов емкостей на узлы вычислительной сетки*;
6) решение уравнения Пуассона методом FFT с учетом как ионов, так и зарядов емкостей*;
7) вычисление напряженности электрического поля в узлах вычислительной сетки*;
8) экстраполяция напряженности поля на ионы*;
9) интегрирование уравнений движения*.
Помимо электрического поля, создаваемого самими ионами, в масс-спектрометрических
устройствах надо учитывать поле электродов: на электроды может подаваться переменное или постоянное напряжение, которое создает дополнительное электрическое поле в ловушке. При вычислении суммарного электрического поля можно использовать два разных подхода. Один из них был описан выше: значения потенциалов электродов подставляются в формулу (1) и вычисляются значения зарядов емкостей, которые создадут в точках электрода нужный потенциал. Это наиболее очевидный метод, но электрическое поле, создаваемое электродами, вычисляется этим методом недостаточно точно. Этот метод использовался нами ранее [8]. Определенного числа клеток сетки, используемых алгоритмом FFT, может быть достаточно для вычисления поля, создаваемого ионами, но недостаточно для вычисления поля электродов с удовлетворительной точностью. Вообще, не оправдано использование метода FFT для вычисления поля электродов, так как он предназначен для расчета полей, создаваемых облаками, состоящими из большого количества заряженных частиц.
Другой метод основывается на том, что в силу принципа суперпозиции можно рассчитать поле электродов и поле ионов отдельно, а затем сложить их, получив результирующее поле. Поле электродов при этом можно рассчитать один раз и прибавлять его к полю ионов, рассчитываемому на каждом временном шаге. Сначала однократно и с существенно большей точностью вычисляется поле электродов при отсутствии ионов и записывается в память компьютера. Это может быть сделано множеством методов, например методом емкостной матрицы с большим числом емкостей, или аналогичным, но более точным методом поверхностных зарядов (Surface Charge Method, SCM).
Рис. 4. Ловушка Кингдона, моделируемая с использованием метода емкостной матрицы
' 0 ' ’ ф\ ’
Чп
X —
_ 0 _ . Фп -
Рис. 3. На рисунке схематически изображено приближение электродов 2-мя наборами точек и вычисление матрицы А. Черные сплошные точки — места, где измеряется потенциал, незакрашенные точки — места, куда помещается поверхностный заряд
В некоторых случаях можно даже не вычислять поле электродов, если известно их аналитическое выражение, например ловушки Поля [9] и Кингдона[6] (если не учитывается искажение поля вследствие конечности размеров ловушки). Такой подход применим и в случае переменного поля — тогда надо по очереди подать потенциал 1 В на каждый электрод, а на остальные нулевой, тогда поле при любых потенциалах электродов может быть найдено как линейная комбинация вычисленных. После вычисления поля электродов начинается основной цикл программы, в котором на каждом шаге методом емкостной матрицы вычисляется поле ионов (при этом потенциал электродов считается нулевым), а к нему добавляется вычисленное заранее поле электродов. Этот метод использовался в настоящей работе.
Рис. 5. Ионные облака, временной сигнал и частот- Рис. 6. Ионные облака, временной сигнал и частот-
ный спектр 2000 ионов ный спектр 2 млн ионов
III. Результаты вычислений
Вышеописанный метод был применен нами для расчета движения ионных облаков в масс-ана-лизаторе масс-спектрометра типа ОгЪ^гар [6] (ловушка Кингдона) (рис. 4).
Вычисления проводились в кубе со стороной 40,6 мм. Форма электродов задаётся уравнением г2 = Г— 1П ((* — номер электрода, 1 соответствует внешнему электроду, 2 —
внутреннему), Я1 = 15 мм, й =6 мм, И1 = 20,3 мм. Поле электродов задавалось аналитически формулой р(т,г) = 2 2 — Г0 + 2, к = 1,412 • 107 В/м. Напряжение между электрода-
ми 2 кВ. Для исследования взаимодействия ионных облаков близких масс в ловушку запускались ионы двух масс: 100 Да и 100,5 Да. Начальная скорость ионов была подобрана таким образом, чтобы они двигались по окружности (в координатах ху) — 4,65 • 104 м/с (соответствует энергии 1,12 кэВ). Был использован временной шаг 3• 10-8 с. Количество емкостей — 1600. Для получения временного сигнала вычислялась разность суммы зарядов на емкостях одной половины ловушки и суммы зарядов противоположной половины. Рассчитанное время движения 4 мс.
Для исследования влияния эффектов объемного заряда на движение ионных облаков было рассчитано движение ионных облаков, состоящих из различного количества ионов. При 1000 ионов каждой массы движение облаков, временной сигнал и фурье-спектр выглядят совершенно невозмущенными ион-ионными взаимодействиями (рис. 5). Для моделирования большего количества ионов нами была применена модель квазиионов, то есть 106 ионов с массой 100 и зарядом 1 моделируются 103 ионов с массой 105 и зарядом 103. Адекватность модели была проверена путем сравнения динамики движения ионных облаков в случае 103 квазиионов с массой 105 зарядом 103 и в случае 104 квазиионов с массой 104 зарядом 102; расчет показал, что движение облаков в этих случаях идентично. На рис. 6 представлены результаты расчета с квазиионами массой 105 и зарядом 103, моделирующими 106 ионов каждой массы (рис. 6). При таком количестве ионов эффекты объемного заряда играют существенную роль. Ионное облако разделяется на три части:
две — это кольца с равномерным распределением ионов, а третья представляет собой небольшой сгусток. Этот сгусток состоит из частиц обеих масс и в течение эксперимента разделяется на два сгустка, намного медленнее, чем основная масса частиц разделяется на кольца. Сгусток двигается с частотой, отличной от частоты движения колец, образованных ионами обеих масс, и на фурье-спектре дает третий пик в дополнение к пикам, соответсвующим ионам двух измеряемых масс.
Расчет 4 мс движения ионов занимает около 10 часов на одном процессоре 1,1 Ггц серверной системы 1ВМ еБегуег рБепез 690. Во всем алгоритме с ростом количества емкостей наиболее вре-мязатратной частью является операция решения системы линейных алгебраических уравнений. В табл. 1 приведены характерные времена для различных процедур используемого алгоритма, требуемые для расчета 1 шага, при разном количестве емкостей, аппроксимирующих электрод.
Таблица 1
Количество емкостей, аппроксимирующих электроды 294 726 1536 3174 5766
Решение уравнения Пуассона методом FFT 0,04 с 0,04 с 0,04 с 0,04 с 0,04 с
Решение СЛАУ 0,005 с 0,015 с 0,045 с 0,3 с 0,9 с
Остальные процедуры 0,2 с 0,2 с 0,2 с 0,2 с 0,2 с
Суммарное время 0,25 с 0,26 с 0,29 с 0,54 с 1,14 с
Характерные времена вычисления различных процедур, выполняемых на каждом шаге, измеренные в секундах. Вычисления производились на 1 процессоре 1,1 Ггц серверной системы ІВМ е8егуег pSeries 690
IV. Выводы
Создан алгоритм, позволяющий применять метод Particle-In-Cell, использовавшийся ранее только для вычислений в ловушках кубической геометрии, для моделирования движения ионов в полях, создаваемых электродами произвольной формы. Метод применен для расчета движения ионных облаков в ловушке Кингдона.
Работа выполнена при поддержке ФЦП «Научные и научно-педагогические кадры инновационной России», грантов РФФИ 10-04-13306-РТ_оми, 09-04-00725-а, 09-03-92500-ИК_а и CRDF RUC1-2941-MO-09.
Литература
1. SIMION 7.0. Idaho National Engineering and Environmental Laboratory.
2. Harlow F.H. A Machine Calculation Method for Hydrodynamic Problems // Los Alamos Scientific Laboratory report LAMS-1956.
3. Dawson J.M. Particle simulation of plasmas // Reviews of Modern Physics. — 1983. — V. 55. — P. 403.
4. O. Buneman A compact non-iterative Poisson solver // Tech. Report 294, Stanford University Institute for Plasma Research, Stanford, CA, 1969.
5. Hockney R.W., Eastwood J.W. Computer simulation using particles. — New York: Adam Hilger, 1988.
6. Hu Q, Noll RJ, Li H., Makarov A., Hardman M., Cooks G. The Orbitrap: a new mass spectrometer // J Mass Spectrom. — 2005. — V. 40. — P. 430-43.
7. Nikolaev E.N., Heeren R.M.A., Popov A.M., Pozdneev A.V., Chingin K. Realistic modeling of ion clouds motion in FT ICR cell by use of a Particle in Cell approach // Rapid Communications in Mass Spectrometry. — 2007. — V. 21. — P. 3527-3546.
8. Nikolaev E., Heeren R., Popov A., Pozdneev A., Vladimirov G. The new possibilities in ion clouds dynamic simulation using supercomputers. Application to FTICR, Kingdon trap and accumulation quadrupole devices //55 th ASMS Conference (2007).
9. Paul W., Steinwedel H. Ein neues Massenspektrometer ohne Magnetfeld // Z fur Naturforsch A. — 1953. — V. 8. — P. 448-450.
Поступила в редакцию 04-04-2011.
