CC BY
26
1. МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
1.1. ПАРАЛЛЕЛЬНОЕ МОДЕЛИРОВАНИЕ ЭЛЕКТРИЧЕСКИХ ПОЛЕЙ ЛОВУШКИ МАСС-СПЕКТРОМЕТРА ДЛЯ ПОВЫШЕНИЯ ТОЧНОСТИ
ИЗМЕРЕНИЯ МАСС ИОНОВ
Попов Александр Михайлович, профессор, доктор физ.-мат. наук, заведующий кафедрой, факультет Вычислительной Математики и Кибернетики. Московский государственный университет имени М.В. Ломоносова, e-mail: [email protected]
Аннотация: работа посвящена численному моделированию поведения ионов в ловушках масс-спектрометров ионно-циклотронного резонанса с Фурье-преобразованием. Определение отношения массы к заряду по циклотронной частоте иона в фиксированном магнитном поле является наиболее точным методом. Определение масс частиц, составляющих молекулу, наиболее важно для понимания структур протеинов, создания новых лекарств. В работе представлен трехмерный параллельный численный код для моделирования поведения ионов в реальной геометрии электродов ловушки и с учетом полей индуцируемых ионами на стенках ловушки.
Метод вычисления полей основан на параллельном итеративном решении интегральных пограничных уравнений теории потенциала.
Разработанная параллельная версия численного кода реализована на суперкомпьютере IBM Blue Gene/P. Новый код позволяет моделировать движение миллиона ионов в течение реального времени детектирования сигнала в эксперименте.
Ключевые слова: математическое моделирование, параллельные вычисления, масс-спектрометр, поведение ионных облаков.
PARALLEL SIMULATIONS OF ELECTRIC FIELDS IN MASS-SPECTROMETER TRAP FOR INCREASING OF IONS MASSES MEASUREMENTS ACCURACY
Popov Alexander Mikhailovich, doctor in physics and mathematics, professor of the Lomonosov Moscow State University, e-mail: [email protected]
Abstract: the work is directed to numerical simulations of ions motion in mass-spectrometer trap on the base of ion-cyclotron resonance with Fourier transform. Determination of mass to charge ratio using cyclotron ion frequency in magnetic field is the most precise method of particles masses included in interesting molecular. Measurements of particles masses included in molecular are important for determination of protein structure and drug design. New three-dimensional parallel numerical code is presented for simulation if ions behavior in real geometry of the trap and takes into account the fields induced on the electrodes of the trap. Method of electric fields calculations is based on parallel solution of integral equations of potential theory.
Developed parallel code version is realized on IBM Blue Gene/P supercomputer at Moscow State University. New code permits one to model of million ions motion in real time of masses detection in experiment.
Index terms: mathematical modeling, parallel computations. Mass-spectrometer, ions clouds behavior.
1. Введение.
Работа посвящена численному моделированию поведения ионов в ловушках масс-спектрометра ионно-циклотронного резонанса с Фурье-преобразованием [1]. Определение отношения массы к заряду по циклотронной частоте иона в фиксированном магнитном поле является наиболее точным методом. Одним из наиболее важных приложений исследований является определение молекулярной структуры протеинов с использованием масс-спектров.
Точность определения масс определяется, в том числе, и процессами, происходящими с ионами в ловушке масс-спектрометра. Не все процессы, влияющие на эволюцию ионов и точность их определения, учитываются при анализе
экспериментальных масс-спектров. Это обуславливает появление математических моделей поведения ионов в ловушке масс-спектрометра. Модели учитывают различные механизмы, которые возникают при взаимодействии ионов, и позволяют оценивать их влияние на точность измерения масс.
Первые попытки моделирования и разработки математических программ предложены Д. Далом при создании кода БМЮЫ [2] для персонального компьютера. Следующий шаг в развитии моделирования представлен в разработанном трехмерном коде частиц в ячейке Д. Митчеллом [3] для моделирования масс-спектрометра на основе ионно-циклотронного резонанса. Многое обусловлено точностью восстановления электрических полей реальной ловушки и полями наведенны-
ми ионами на стенках ловушки. Методы конечных разностей, использованные в полном коде PIC3D [2], [3], могут быть использованы только для замкнутых поверхностей электродов.
Новый программный комплекс, ориентированный на современные суперкомпьютеры для моделирования движения ионов, учитывающий реалистичную конфигурацию полей, представлен в работах [4], [5],[6].
В новых версиях трехмерного кода частиц в ячейке PIC3D используется конфигурация полей, вычисленная с помощью параллельного кода.
В работе [7] представлена новая параллельная гибридная версия кода PIC3D - программа, сочетающая CPU (Central Processing Units, центральные процессоры) и GPU (Graphic Processing Units, графические ускорители) вычислительных систем. Расчеты выполнялись на системах с различными GPU, в том числе, на суперкомпьютере «Ломоносов» НИВЦ МГУ, где имеется возможность использования гибридных вычислительных систем с CPU и GPU. Система позволяет реализовывать различные стратегии распределения параллельных процессов по CPU устройствам с параллельным выполнением части задач на GPU.
В данной работе представлен метод вычисления потенциала электрического поля ловушки, основанного на граничных интегральных уравнениях теории потенциала. По сравнению с другими методами, используемыми в масс-спектрометрических вычислениях, предложенный метод имеет целый ряд преимуществ. Во-первых, может быть рассмотрена любая конфигурация геометрии незамкнутых электродов для создания электрических полей. Поле в замкнутом проводнике при постоянном потенциале должно быть равно нулю и использование таких полей является, вообще говоря, неверным. Поля от незамкнутых электродов могут быть вычислены во всем пространстве с учетом щелей между электродами и наведенными ионами электрических полей во всем пространстве, влияющих на точность измерений. Во-вторых, поля могут быть вычислены на произвольных сетках внутри ловушки, позволяя иметь желаемую точность в различных областях. Рассмотрено введение системы специальных секций электродов с целью компенсации разницы между идеальными условиями ловушки и реалистичными конфигурациями полей. Верификация кода проведена для ловушки двухкамерного масс-спектрометра «O-trap», предложенной в [8].
Моделирование движения ионов и свойства спектра получены с помощью PIC3D. Представлено моделирование движения ионов и свойства спектра в скомпенсированном поле с использованием протокола SWIFT.
Разработанная параллельная версия численного кода, реализована на суперкомпьютере IBM Blue Gene/P, установленном на факультете ВМК МГУ. Параллельная версия кода позволяет моделировать движение миллиона ионов в течение реального времени детектирования сигнала в эксперименте.
2. Моделирование движения ионов с использованием кода PIC3D в ловушке масс-спектрометра
Для моделирования движения ионов в ловушке масс-спектрометра используется модель частиц-в-ячейке [9]. Уравнения движения ионов внутри ловушки под действием куло-новских сил и силы Лоренца записываются в виде:
drj dt
= Vi
dvt
= 4iEtrap (rt, t) + 4i [vrB]
+
i у q^iCn—j 4ne0j=i |ri- r
+ Ewall(ri)
j=1 ri j
+Fcoil(Vi)
+Eexct(ri) cos Y(t)
где электрическое поле E определяется по вычисленным потенциалам: потенциалу удерживающего поля ловушки, потенциалу высокочастотного поля, возбуждающего ионы, кулонов-скому потенциалу парного взаимодействия частиц, который необходимо рассчитывать на каждом шаге по времени и учитывать в форме частица-сетка, значению потенциала на электродах ловушки, который вычисляется на каждом шаге и учитывает влияние стенки на движение ионов. Потенциал коллективного взаимодействия ионов рассчитывается путем решения уравнения Пуассона на разностной трехмерной сетке.
В модели частиц в ячейке вводится пространственная трехмерная разностная сетка. В рассматриваемой области потенциал кулоновского взаимодействия определяется плотностью заряда:
р(к, I, т) = ^^ qlWkilim(xl,yt,zt)
Вычисления проводятся на сетке при помощи интерполяции зарядов частиц внутри ячейки, определяемой сеткой. Здесь Vol - объем ячейки, W - веса интерполяции. Тогда коллективный потенциал кулоновского поля находится из решения первой краевой задачи для уравнения Пуассона на пространственной сетке:
АФ = —1 p(rl,t)
Ф| Ь ound = V(r,t) и по найденному потенциалу вычисляется сила, действующая на все ионы внутри ячейки.
3. Идеальный 3D аксиальный квадрупольный потенциал в ловушке Пеннинга. Аналитическое решение.
В области
г < а; —L < z < L. идеальный потенциал представляется в следующей форме:
Vid = —
2 U,
V2 I ■ а ■
_г 2 / mm ^ mm4-|
а~[Z - ~2 - (~2 ~)]
¡2 I min 1 min "г" 2
X2 + V2
Via = - Va[z2--^--С]
Vo =
2 U o
С =
'-min ^min
i2 I min 1 min "г" 2
Боковые электроды имеют потенциал
V(0,z= ± lmln )= - Uo.
Для настоящей работы важны следующие обстоятельства. Для каждой из поверхностей электродов мы предполагаем, что на поверхности проводников должны поддерживать постоянные значения потенциалов. Однако, в реальности для проводящего электрода это означает, что он представляет собой замкнутую поверхность, на которой потенциал постоянен, а внутри поверхности имеет нулевой потенциал. В реальности поверхности имеют разрез и представляют собой части двухсторонних поверностей и внутри такого цилиндра потенциал ненулевой. Если разрез тонкий, то искажение поля на концах мало. Однако искажение существует и может оказывать влияние на движение иона. Так будет с любыми частями ловушки масс-спектрометра. Эффект будет влиять на точность определения масс ионов.
4. Определение потенциалов удерживающего поля от электродов произвольной формы. Интегральные пограничные уравнения для электрического поля ловушки
Задача заключается в нахождении потенциала У(х,у,г), удовлетворяющего уравнению Лапласа в области ионной ловушки:
АУ = 0
с граничными условиями первого рода на электродах, окружающих ловушку:
= '
где к = 1,2,3,...,К - номер поверхности электрода. Правая часть граничного условия определяет заданный постоянный потенциал на каждой поверхности. Ищем решение в виде потенциала двойного слоя
Ш(М) = -
II — (—
■)р(Р)йар
Шк (М) = -
II — (—
]] дп„ Ямр
(ъ—ЖРк )йар.
ме такие члены с интегрируемой особенностью должны иметь специальную аппроксимацию. Внешняя часть поверхности является отталкивающей, а внутренняя - притягивающей.
2тк(Р(к)) -
II — (—
-К (Р)ЛаР
= иг
+ I II ш;(^
При этом точка наблюдения тоже находится на этой поверхности. В правой части уравнения имеется два члена. Первый член есть приложенный потенциал на поверхности к, а второй - сумма потенциалов от всех поверхностей, которые индуцируют электрическое поле на данной поверхности. Уравнения должны решаться итерационно для каждой поверхности системы. Таким образом, система для дипольных моментов получается связанной. Идея итерационной процедуры состоит в фиксации дипольных моментов для всех поверхностей кроме й-ой. Поскольку правая часть уравнения известна на текущем итерационном шаге, она может быть вычислена. Как только решения для всех поверхностей на этом шаге будут найдены, мы сможем вычислить и новую правую часть.
Таким образом, мы переписываем систему уравнений в форме:
2пук(Р(Ю) -
^ дпр (И1к]
-)рк (Р)йар = Фк
с правой частью:
Здесь Р есть точка на поверхности и М - точка наблюдения трехмерной разностной сетки, покрывающей объем ловушки, Рмр - расстояние между точкой М и точкой Р, пр - нормаль к поверхности в точке Р, ур является неизвестной поверхностной плотностью дипольного момента, которая различна на каждой поверхности. Для любой нормали рр потенциал W(M) удовлетворяет уравнению Лапласа по переменной М. На поверхности данный потенциал имеет разрыв. В случае одной поверхности мы приходим к одному интегральному уравнению для неизвестной плотности момента р(Р):
и1 = Ш(Ро) + 2пу(Р0), В случае многих поверхностей решение ищется как сумма потенциалов двойного слоя, создаваемого каждой поверхностью:
ф>=и> +11^(р1
(Ю
2пу!(Р^) -I Ц — (Р)йа„ = и1
1=11
1р
1=1 ^ р ^
К
= 1
-т
V, (Р)йар + 2крк(Р}к-1) = ик
Дипольные моменты получаются из граничного условия для потенциала двойного слоя, когда мы помещаем точку М сначала на первую поверхность (первое уравнение), затем на вторую поверхность (второе уравнение), и на к-ю поверхность (к-е уравнение). Таким образом, мы получаем систему уравнений для дипольных моментов, которая на каждой поверхности определяется граничными условиями, а именно значениями потенциала на каждой поверхности и потенциалами, наведенными другими поверхностями. В каждом уравнении один из интегралов в сумме имеет сингулярность. В численной схе-
Эта система уравнений для дипольных моментов на каждой поверхности определена граничными значениями потенциалов на каждой поверхности, индуцированных другими поверхностями. В каждом уравнении один интеграл в сумме имеет особенность, когда £ = к и соответствующие координаты точек Р. Это интегрируемая особенность, но для численной схемы такие члены должны иметь специальную аппроксимацию. Внешняя сторона поверхности отражательная, а внутренняя часть притягивающая. Угол между внутренней нормалью и
направлением от точки Р на поверхности к фиксированной точке М. Поверхности имеют две стороны. После нахождения дипольных моментов, потенциал может быть найден как потенциал двойного слоя:
к
у(м) = --Щ1ыр ^^
где F учитывает квадратурные формулы, и есть решение на п-ой поверхности. Полная система алгебраических уравнений приобретает следующую блочную форму:
N1 N2 ( N1 N2 Л
2uu\m -УУ F™, ^ - V J V 'V vW ,>Н _ п
kmij ij
1=0 j=0
У)уУ ^kmijUiJ f = fkm
=2 =o =o
5. Аппроксимация 2Р интегралов на сети
Рассмотрим поверхность типа (а) - однозначное проектирование на плоскость (х,у). Введем сетку на области О:
х^ < х < хК; 2а <2 < .
х0 = Хь, = ; Х1 = Х0 + I ' ^х; I = . .,
1\х
¿0 = 2а; К = г" - ^; = 2о +} ' к2;) = 1.....Ы2
2Р интегралы представлены в форме:
¡а(хро,гро) = j ^р j Н (хро,гро; хр,гр)
ъа Х1
n0-l N1 N2
1 2
У {УУ F&jvtf] + 2п и°т - УУ FЙ> ^
п=1 1=0j=0
=o =o
Nn N1 N2
У УУ F^}=a
n=n0+1 i=0 j=0
Nn-1 ( N1 N2
n=1 I ¡=0 j=0 J I ¡=0 j=0
Сетка введена для каждой независимой переменной
Хр ^ Х1, 2р ^ Zj, Хро ^ Хк, 2ро ^ 2т.
Теперь мы введем интегрируемую функцию в каждой точке на сетке:
Н(Хк, 2т; Х0 2) ) = ¥а(Хк, ; Хь, 'Ч(ХЬ, 2) ), где и(х^, ) = р(х^, ) соответствуют дипольному моменту. Так, и есть разыскиваемое решение на сетке. Здесь мы используем квадратурную формулу трапеций для интегралов.
N¿-1 Ях-1 Ых-1
¡а^УУ Н(хи г]) кхк2 + 2 У (Н(хи г0) + Н(хи гЫг))кхкг
]=1 1=1 1=1
нх-1
+ 2 У (Н(х0'zj) + н(xNx, Zj))hxhz
где п = 1,...,Мп есть индексы поверхности, и есть полное число поверхностей. Видно, что система имеет блочную структуру, отражающую влияние каждой поверхности. Мы используем блочный итеративный метод Гаусса-Зейделя для решения системы. Пусть 5 есть индекс итераций. Затем мы получаем следующую итерационную процедуру:
-УУ&чг* 1=0 ¡=0
' N1 N2
У\УУ РктЦи(Т f = fkm п=2 I i=0 j=0 )
+ 4 Щ(х0, г0) + Н(хМх, г0) + н(х0,
После этого мы имеем интегральное уравнение в алгебраической форме:
2ликт - У У РВкт.. ' Щ] = [кт 1=0 ¡=0
Эти уравнения мы перепишем в форме стандартной алгебраической системы АУ = Ф, где У есть вектор решения, Ф есть линейный вектор известной правой части и А есть двумерная матрица, которая индуцирует геометрические характеристики электродов. Так как матрица плотная, мы используем Ш декомпозицию для решения этой системы уравнений
6. Решение системы интегральных уравнений для набора поверхностей
После дискретизации интегралов полная система уравнений для набора поверхностей может быть записана в форме алгебраических уравнений. Определим интегралы для п-ой поверхности после дискретизации в следующем виде:
(n0 1) N1 N2
n=1 i=0j=0
N1 N2
/V Х"1 Z7("o,s+1),,(no,s+1)-i
{LiL Fkmii UiJ }
i=0 j=0 Nn N1 N2
n=n0+1 i=0 j=0
Nn-1 N1 N2 ( N1 N2
У {УУр^и^^
n=1 i=0 j=0
kmij "ij } j 2_,У' kmij "ij
F(Nn,s+1)u(Nn,s+1) I
,i=0 j=0
+ 2nuNn,s+1 = fNn -t-2nukm — lkm.
N1 N2
УУ
i=0 j=0
p(n) u(n) kmij ij
Так, для некоторой поверхности с номером п0 мы имеем следующую систему уравнений, которую необходимо решить:
N1 N2
i=0 j=0
(V V r(no,s+lL,(n0,s+l)x
{L>h Fkmij Ulj }
no-l f N1
= &+
n=1 ^¡=0 j = 0
Nn N1 N2 + [
,s+1)u(n,s+1) I kmij ij 1
п=п0+1 ¿=0 ¡=0
Здесь правая часть известна. Первый член в правой части
есть известный потенциал, второй член находится на предыдущем шаге процедуры Зейделя, третий член известен с предыдущего шага итераций. Второй и третий член представляют влияние других поверхностей на распределение диполь-ного момента на данной поверхности. Каждое такое уравнение мы перепишем в форме стандартной алгебраической системы:
АУ = Ф,
где У есть вектор решения, Ф есть вектор известной правой части и А есть двумерная матрица, которая включает геометрические характеристики электродов.
7. Параллельный алгоритм
Идея параллельного алгоритма представлена на рис. 1.
Рис. 1 Структура параллельного алгоритма. Здесь «пб», «пс!», - номера поверхностей электродов, уПс|(пб,пс,....) - потенциал поля индуцированный поверхностью с номером «пС» на поверхности с номером «пб». Все потенциалы складываются в поле 113сС и параллельный алгоритм вычисления потенциалов на каждой поверхности повторяется.
На каждой итерации алгоритма электрический потенциал от каждой поверхности электрода может быть определен независимо. Это реализовано параллельными процессами, выполняющими вычисления потенциала на каждой поверхности ns ns = 1,2,...,т. На каждой итерации определяются поля Vind(ns,па, ., .), индуцируемые поверхностью ns, на всех других поверхностях па. Вычисленные поля отсылаются процессу с номером 0, используя операции передачи сообщений MPI. Вычисленное нулевым процессом тотальное индуцированное поле отсылается обратно в соответствующий вычислительный процесс с номером ns, который использует полученные данные для вычисления новых правых частей уравнения. На следующей итерации, модифицируемые граничные условия используются для вычисления соответствующей поверхности.
Вычисления, описанные выше, могут быть выполнены и с использованием последовательного алгоритма и только одного вычислительного процесса. Реализации последовательной и параллельной версий алгоритма выполнены на суперкомпьютере IBM Blue Gene/P, установленном на факультете ВМК МГУ. Расчеты были проведены для 19 по-
верхностей электродов с сеткой размером 23X23 на каждой поверхности. Пространственный объем расчетов содержит 32x32x32 узлов сетки. Только одна итерация алгоритма была сделана в обоих случаях.
Время выполнения одной итерации с использованием последовательного алгоритма на одном процессоре составило 18 минут 52 секунды. Параллельная версия алгоритма на 20-ти процессорах суперкомпьютера Blue Gene/P потребовала 44 секунды. Ускорение вычислений с использованием параллельной версии алгоритма составило порядка 20 раз по сравнению с последовательными вычислениями. Заметим, что только принципиальная параллелизация вычислений выполнена. В методе есть дополнительные места, где может быть использовано распараллеливание.
8. Зависимость точности получения масс-спектра молекулярных ионов от особенностей математической модели эволюции ионов и возможностей численной модели
Сначала обсудим особенности результатов расчетов для реальной геометрии электродов.
Приведем расчет удерживающего поля ловушки в масс-спектрометре О-трэп. Это поле не может быть задано анали-
тически. Электрическое поле присутствует в двух секциях -возбуждения и детектирования. Линии уровня поля показаны на рис. 2. Линии уровня приведены для удерживающего электрического потенциала для системы электродов для случая с отсутствием индуцированных полей (а) и случая, когда индуцированные поля берутся в расчет (б). Без учета присутствия ионов деформация поля на рис 2(а) отсутствует. Учет индуцированного электрического поля, создаваемого ионами в секции детектирования, оказывается существенным для эволюции ионов.
(б)
Рис. 2. Линии уровня удерживающего электрического потенциала
для системы электродов для случая с отсутствием индуцированньк полей (а) и когда индуцированные поля берутся в расчет (б).
9. Моделирование поля в секции возбуждающих электродов
В идеальной ловушке Пеннинга возбуждающее поле является полем плоского конденсатора. В реальных условиях важны не только краевые эффекты плоского случая, но и кривизна поверхности электродов.
Для моделирования используется представление поверхностей в виде набора секций. Для 6-ти секций конфигураия электродов показана на рис. 3.
Рис. 3 Конечная модель бесконечного конденсатора используемая в возбуждающей секции. Наличие подсекций позволяет регулировать необходимую конфигурацию поля.
J:
V ... - —1—___ /
Рис. 4. Вычисленные профили электрических полей в удерживающей и возбуждающей ячейке ловушки О-трэп. Для четырех секций показано постепенное приближение поля к необходимому. В центре ячейки
поле близко к ожидаемому. В первой, возбуждающей секции линии уровня поля выглядят следующим образом.
Рис. 5. Профили поля в возбуждающей секции ловушки О-трэп. Линии уровня электрического поля. 10. Численное моделирование масс-спектра
т/<\ (6)
Рис. 6. (а) Масс-спектр для молекулы содержащей 500 ионов, (б) иллюстрация точности воспроизведения спектральной линии, обусловленная кулоновским взаимодействием ионов, полученным методом Монте-Карло. Рисунок 7 показывает эволюцию облаков ионов в ловушке Пеннинга. В масс-спектре наблюдается слияние спектральных
пиков для трех близких масс. Слияние определяет точность воспроизведения масс.
б) Число ионов Ыр=7500. Рис. 7 Эволюция ионных облаков для разного количества ионов: (а) Ыр=25000, (Ь) Ыр=7500 .Три различных массы т/р = 99.67, 100.01, 100.31 присутствуют в кластере ионов. Слияние ионов для бл изких масс определяет точность метода. Сравнение расчетов с экспериментом.
б) Численное моделирование с помощью кода РОЭ. Рис. 8. Сравнение результатов (а) эксперимента, выполненного в Институте Атомной и Молекулярной физики в Амстердаме (АМОи [6]) и (б) математическим моделированием спектра данной молекулы с помощью кода РОЭ.
11. Выводы
В статье представлен параллельный численный код, разработанный на основе имевшейся программы PIC3D. Код позволяет проводить расчет для произвольной геометрии ловушки масс-спектрометра с различным числом электродов с учетом наведенных зарядов. Реализована полная схема вычислений масс-спектра: удержание ионов с учетом их взаимодействия и столкновений, возбуждение ионов, детектирование сигнала и его Фурье представление.
Показана возможность коррекции экспериментальных измерений с учетом реальной формы электродов, наведенного ионами электрического поля на электродах, формы режима возбуждения и детектирования. Выполнена оптимизация параллельного алгоритма при реализации кода на суперкомпьютере IBM Blue Gene/P, что позволяет проводить вычисления эволюции миллиона ионов в течение реального времени эксперимента.
Список литературы:
1. A.G. Marshal, C.G. Hendrickson, and G.S. Jackson. Fourier Transform Ion Cyclotron Resonance Mass Spectrometry: A Primer. - Mass Spectrom. Rev., 17, pp. 1-35, (1998)
2. D.A. Dahl.-SIMION for personal computer in reflection.- Int. J. Mass Spectrom 200, pp. 3-25, (2000)
3. D.W.Mitchell.- Realistic simulation of the ion cyclotron resonance mass spectrometer using a distributed three-dimensional particle-in-cell code.-J.Am. Soc. Mass Spectrom 10, pp.136-152, (1999)
4. E.N. Nikolaev, A.M. Popov, R.M.A. Heeren, et.al. Realistic modeling of ion motion in FT ICR cell/- Proceedings of the 5th North American FT-ICR MS Conference Key West, Florida, p.27, (2005)
5. E.N. Nikolaev, R.M.A. Heeren, A. M. Popov et. al. Current progress in supercomputer modeling of ion clouds behavior in FT ICR cells.- Proceedings of the 5th ASMS Conference on Mass Spectrometry Allied Topics .- Seattle Washington May 28-June 1, (2006)
6. E.N. Nikolaev, A.M. Popov, R.M.A. Heeren, et. al.; Rapid Commun. Mass Spectrom. 2007, 21, 3527.
7. Н.Н. Попова, Н.Г. Никишин. Численный код частиц в ячейке на основе гибридных CPU-GPU вычислительных систем. Программные продукты и системы. 2014. Вып. 1. 36-43.
8. Alexander S Misharin and Roman A. Zubarev.- Coaxial multi-electrode cell («O-trap») for high-sensivity detection at a multiple frequency in Fourier transform ion cyclotron raonance mass spectrometry: bmain design and modeling results.-Rapid Commun. Mass Spectrom. 20, pp.3223-3228 (2006).
9. Р. Хокни, Дж. Иствуд.- Численное моделирование методом частиц.- Москва «Мир» 1987, 638 С.
10. А.Н. Тихонов, А.А.Самарский. Уравнения математической физики. Изд. «Наука»,1977, 736 C.
11. Using MPI: Portable Parallel Programming with the Message-Passing Interface. Wiliam Gropp. Ewing Lusk .- The MIT Press, (1997).
