CC BY
8ISSN 0868-5886 НАУЧНОЕ ПРИБОРОСТРОЕНИЕ, 2016, том 26, № 4, c. 64-76
- МАТЕМАТИЧЕСКИЕ МЕТОДЫ г _
И МОДЕЛИРОВАНИЕ В ПРИБОРОСТРОЕНИИ
УДК 537.534.7 © А. П. Щербаков
ПРОЦЕДУРА МОДЕЛИРОВАНИЯ СТОЛКНОВЕНИЙ ИОНОВ С МОЛЕКУЛАМИ В НЕОДНОРОДНЫХ ПЕРЕМЕННЫХ ВО ВРЕМЕНИ ЭЛЕКТРОГАЗОДИНАМИЧЕСКИХ ПОЛЯХ
Разработаны и протестированы процедуры компьютерного моделирования движения ионов в неоднородных переменных во времени электрогазодинамических полях. Процедуры основаны на комбинированной модели взаимодействия ионов с молекулами газа, описывающей рассеяние ионов в широкой области энергий и учитывающей единым образом взаимодействие иона с молекулой как в области отталкивания, так и в области притяжения. Разработанные модели и соответствующие вычислительные процедуры позволяют вводить неизотропные модели рассеяния. Предложенная монте-карловская процедура моделирования длины или времени свободного пробега основана на дискретизации по времени, что естественным образом сопрягается с дискретизацией уравнений движения.
Кл. сл.: сечение столкновения, поляризационный потенциал, потенциал упругих шаров, подвижность, диффузия, электрогазодинамическое поле
ВВЕДЕНИЕ
Транспортирующие газонаполненные ионные тракты на основе радиочастотных полей с мульти-польной пространственной структурой нашли широкое применение в масс-спектрометрии. Они используются в качестве промежуточной системы, соединяющей область высокого давления ионного источника (Electrospray, MALDI) с высоковакуумной областью масс-анализатора [1, 2]. Помимо высоких транспортирующих качеств такого рода устройства эффективно сжимают фазовый объем пучка, что имеет большое значение при анализе потоков заряженных частиц с большим начальным разбросом по энергиям и углам. Последнее свойство дало основание назвать эти устройства транспортирующим ионным трактом со столкно-вительной фокусировкой (Collisional Focusing Ion Guide, CFIG) [2].
Методика транспортировки ионов радиочастотным мультиполем нашла применение и в ядерной физике для масс-анализа экзотических корот-коживущих ядер с помощью ионных ловушек. Получаемые в сепараторах продуктов ядерных реакций потоки быстрых заряженных частиц с большими эмиттансами тормозятся в газовой ячейке, а затем с помощью радиочастотных квад-руполей эти потоки трансформируются в пространственно сжатые низкоэнергетические ионные пучки, пригодные для анализа в ловушках [3].
Эффективный потенциал, формируемый неоднородным радиочастотным полем [4], ограничива-
ет ионы в поперечном сечении, столкновения с молекулами буферного газа приводят к диссипации начальной кинетической энергии ионов. Пучок охлаждается и стягивается к оси системы — фокусируется. Заметим, что эффективный потенциал не зависит от знака заряда иона. Это позволяет собирать и удерживать ионы разной полярности в общей области, что широко используется в газофазной химии при исследовании реакций ионов противоположной полярности, в которых по крайней мере один из реактантов является многозарядным — такого рода ионы легко образуются при электроспрей-ионизации биополимеров [5].
Для оптимизации радиочастотных мультипо-лей, используемых в качестве транспортирующих систем или ловушек, необходима модель динамики ионов в таких устройствах. Эта модель должна адекватно описывать траектории ионов в переменных во времени электрических полях заданной пространственной конфигурации, взаимодействие ионов с молекулами буферного газа с учетом их хаотического движения и потока газа в целом. Такая модель необходима также для интерпретации результатов масс-спектрометрического эксперимента, например для анализа причин масс-дискриминации, для их устранения или учета при обработке данных.
Адекватность модели взаимодействия означает ее пригодность и достаточную точность при описании столкновений ионов в требуемом диапазоне энергий. Этот диапазон простирается от тепловых энергий (~ 10-2 эВ) до нескольких сотен эВ. По-
этому модель должна учитывать как отталкива-тельную ветвь потенциала взаимодействия, ответственную за столкновения при достаточно больших энергиях, так и притягивающую ветвь, определяющую рассеяние при малых энергиях.
Динамика движения ионов в радиочастотных мультиполях рассматривалась в работах [2, 6]. Однако развитые в этих работах модели, как аналитические [2], так и численные [6], основаны на концепции взаимодействия абсолютно упругих шаров. Между тем это приближение при малых энергиях справедливо лишь для достаточно больших молекул, когда сечение "жесткой сердцевины" больше сечения поляризационного взаимодействия. Для ионов средних и малых масс при небольших значениях энергии столкновения оба сечения становятся сравнимыми. Поэтому при вычислении пространственно-энергетических распределений охлажденного пучка ионов модель упругих шаров должна быть дополнена моделью поляризационного взаимодействия.
Необходимость совместного учета как поляризационного, так и взаимодействия абсолютно упругих шаров отмечалась еще в работе [7]. В работе была развита процедура компьютерного моделирования, основанная на модели постоянной средней частоты столкновений (модель Максвелла, поляризационное взаимодействие) и на модели постоянной средней длины свободного пробега (приближение постоянного эффективного сечения, модель жесткой сферы). Комбинация этих двух моделей в рамках приближения фиктивных столкновений была предложена в работе [8]. Отметим, что приближение фиктивных столкновений позволяет учитывать в рамках единой процедуры несколько типов взаимодействия и различные типы реакций (например, помимо упругого рассеяния, фрагментацию ионов).
Приближение фиктивных столкновений оказалось эффективным при моделировании движения ионов в неоднородных электрогазодинамических полях. Развитые на его основе численные процедуры были использованы при компьютерном моделировании и исследовании газодинамических источников ионов [9, 10].
Наличие фактора времени и соответствующих параметров, его описывающих, — частота (или спектр частот) внешнего переменного во времени электрического поля — требуют анализа соотношений между этой частотой и временн м параметром, характеризующим процесс взаимодействия с газом, — частотой столкновений.
В настоящей работе предложена комбинированная модель, описывающая рассеяние ионов в широкой области энергий и учитывающая единым образом взаимодействие иона с молекулой как в области отталкивания, так и в области при-
тяжения. Модель основана на приближении эффективного радиуса взаимодействия в области отталкивания и на приближении Ланжевена в области поляризационного притяжения. Модель учитывает скоростную зависимость сечений взаимодействия, позволяет вводить неизотропные модели рассеяния и применима при расчете траекторий движения ионов в неоднородных газовых потоках. Динамика ионов представляется в виде монте-карловской процедуры, включающей в себя моделирование длины или времени свободного пробега, вычисление траектории иона в неоднородном радиочастотном электрическом поле и моделирование столкновения с молекулой буферного газа. Тестирование вычислительной процедуры производится на известных результатах задачи о диффузии и дрейфе ионов в газе под действием однородного электрического поля.
СЕЧЕНИЕ СТОЛКНОВЕНИЯ
В приближении эффективного радиуса столкновения [13] сечение рассеяния иона на молекуле
газа определяется расстоянием
наибольшего
сближения сталкивающихся частиц, которое находится из равенства
V(г0) = Е,
(1)
где V (г) — потенциал взаимодействия, Е = ¡и2 /2 — энергия относительного движения, Л = тМ/ (т + М) — приведенная масса, т — масса иона, М — масса молекулы, и — скорость относительного движения.
Если для отталкивательной ветви потенциала взаимодействия взять часто используемую аппроксимацию [12]
V (г) = А ехр (-г / а), (2)
то для расстояния наибольшего сближения полу-
и сечение <г аппрок-
чаем
г0 = а 1п(2А / (¡и2)),
симируется выражением
Л,2 Л
<г = <0 1П
У,
(3)
В этой формуле введены параметры <0 = па2 и и02 = 2А / л . В отличие от модели абсолютно жестких шаров сечение (3) зависит от скорости, хотя и слабо.
Приведем некоторые значения параметров из работы [12]. Для пары Аг—02: А = 14.18кэВ,
а = 0.240Л, < = 0.181 х10-16см2, и0 =
= 3.91 х 107см/с. Для пары Ar—N2: A = 28.76кэВ ,
■ = 0.289 А,
< = 0.262 х10-16см2
= 5.79 х107см/с.
Ветвь потенциальной кривой, соответствующая притяжению иона и молекулы с нулевым диполь-ным моментом, определяется поляризационным взаимодействием — взаимодействием заряда с наведенным дипольным моментом:
V (r) = -
а q
2r4
(4)
<p =
2nq
а
U
(5)
2.3 х 10-9 z
<p =■
v
(6)
МОДЕЛИРОВАНИЕ СТОЛКНОВЕНИЯ
Скорость иона массы т после его столкновения с молекулой массы М определяется из соотношения [14]
M
и =-
vk + -
m
-U +-
M
(7)
т +М т +М т + М
Здесь и = и1 - им, и = |и| — скорость относительного движения, и1м — скорости до столкновения
соответственно иона и молекулы.
Плотность вероятности распределения случайного единичного вектора рассеяния к соответствует дифференциальному сечению рассеяния
и определяется потенциалом взаимодействия сталкивающихся частиц.
Для моделирования неравномерных распределений необходимо перейти в систему координат, одна из осей которой направлена вдоль вектора относительной скорости. Пусть в системе координат ХУХ направление вектора относительной скорости задается полярным а и азимутальным /3 углами (см. рис. 1). Тогда его компоненты суть
vX = vsinacos f, vY = vsinasin ff,
v7 = vcosa.
(8)
где ар — поляризуемость нейтрального атома или
молекулы, q — заряд иона. Этот потенциал допускает орбитирование [13].
Сечение рассеяния на потенциале (4) обычно аппроксимируют сечением поляризационного захвата на орбиту (орбитирования), которое выражается формулой Ланжевена [13]
Замечательным свойством сечения поляризационного захвата является постоянство константы скорости реакции crpv = const, которым объясняется, например, независимость подвижности ионов в газе от температуры при малых значениях величины электрического поля [11].
Если в последнем равенстве выразить сечение < в см2, заряд иона q в единицах элементарного
заряда z, скорость v в см/с, поляризуемость ар
в А3, приведенную массу u в атомных единицах массы, то получим
Для того чтобы ось X была направлена вдоль вектора и, необходимо сделать два поворота: первый — на угол /3 вокруг оси X, второй — на угол а вокруг оси X' — нового положения оси X (см. рис. 1).
Матрица первого поворота
O z (Р) =
cos . sin Р 01
- sin Р cos Р 0
0 0 1,
(9)
* У
Рис. 1. Кинематика столкновения. XXX — исходная система координат; и — направление вектора относительной скорости до столкновения, задается полярным а и азимутальным / углами; X У Х' — система координат, в которой полярная ось направлена вдоль вектора относительной скорости
V
и
Матрица второго поворота
Or (a) =
(cosa 0 - sinai 0 1 0 sin a 0 cosa
(10)
Пересчет координат вектора при переходе от системы XXX к системе ХТ'X' дается произведением матриц
O7 (a) Oz (p) =
cosacosP cosasinp -sina
-sinp cosP 0 sinacosP sinasin p cosa
. (11)
Обратный переход от системы ХТX' к системе XXX производится с помощью транспонированной матрицы
Í,
Oz (-P) Or (-a) =
cosacosP -sinp sinacosP cosasinp cosP sinasinp
-sina
0
cosa
. (12)
dff _ a2
dQ= T'
(14)
где а — радиус взаимодействия (полное сечение рассеяния — па2), dQ = 2nsm3d3 . Откуда нормированная плотность распределения полярного угла рассеяния 3 на промежутке [0,п)
p(3) = (1/2) sin $.
(15)
Розыгрыш случайной величины % , распределенной с плотностью вероятности р(х), производится согласно лемме [15], лежащей в основе практической реализации метода статистического моделирования:
R =
Ъ
J p(x)dx,
(16)
где R — случайная величина, равномерно распределенная в интервале (0,1). Используя (15), получаем
cos$ = 1 - 2R,
р = 2л R.
(17)
Таким образом, согласно соотношениям (13), розыгрыш случайных компонент единичного вектора рассеяния в случае изотропного рассеяния производится по формулам
кг = 1 - 2R1,
кх, =yJ1- К cos (2лR2), ky, =yj 1 - k2, sin (2лR2),
Компоненты единичного вектора рассеяния к в преобразованной системе координат ХТ X' определяются полярным 3 и азимутальным р углами:
кх, = 8ш3с08(, кг = 8т38тр, Кг = С083. (13)
Для центральных потенциалов взаимодействия между сталкивающимися частицами азимутальный угол рассеяния р распределен равномерно на промежутке [0,2п) . Плотность распределения полярного угла рассеяния 3 определяется потенциалом взаимодействия.
Для изотропного рассеяния, имеющего место, например, при столкновении абсолютно упругих шаров, дифференциальное сечение рассеяния:
(18)
где Rl и R2 — две независимые реализации равномерно распределенной на (0,1) случайной величины.
В общем случае для произвольного центрального потенциала при моделировании полярного угла рассеяния удобнее сначала разыгрывать случайный прицельный параметр Ь с плотностью вероятности р(Ь) ~ Ь , а затем вычислять искомый угол с помощью классической функции отклонения 3 = 3(Ь) [16], которая в ряде важных случаев имеет простой вид. Например, для степенного потенциала V (г) = с ¡гп в области малых углов
т=
с4ж Г((п +1)/2)
Ebn
Г(п/ 2)
(19)
Сгенерировав случайные углы рассеяния 3 и р, мы определяем компоненты единичного вектора рассеяния к в системе координат ХТ'X' (формулы (13) или (18)), которые с помощью преобразования (12) пересчитываются в исходную систему координат XXX. Выполнив это преобразование, из выражения (7) получаем для декартовых компонент вектора скорости иона после столкновения:
(20)
, т
»х =-77 » х +
т + М
М
+-— (Кх» 0083 - Ку,и8Ш3 + к2,1>х + »мх),
т + М
,т
»т =-77 » +
т + М
М
+--(кх8Ш3 + Ку,и0083 + КХ+ иму ),
т + М
,т
» =-77 »12 +
т + М
М
+-77 (-Кх »хт + кх ^ +»мх ).
т + М
Здесь учтены соотношения (8) для вектора относительной скорости и введено обозначение ихт = »81па .
Таким образом, зная компоненты векторов скорости иона и молекулы до столкновения, с помощью соотношения (20) можно вычислить компоненты вектора скорости иона после столкновения. При этом скорость молекулы до столкновения разыгрывается как случайная величина, декартовы компоненты которой в нашей модели распределены по нормальному закону с постоянным смещением и, описывающим скорость потока газа в целом:
(»м),=Л
3 + из,
3 = X, у, г,
(21)
(»м) , = 91.2 х102Л/—N, ] = X, у, г.
Моменты времени 0 < t1 < t2 < ... или соответственно координаты точек 0 < 51 < s2 < ... столкновений образуют поток событий [17]. События, которые происходят по альтернативным сценариям (рассеяние на отталкивательной ветви потенциала взаимодействия и рассеяние на притягивающей ветви), можно трактовать как суперпозицию соответствующих потоков событий, и для моделирования длины или времени свободного пробега использовать соответствующие процедуры, приведенные, например, в [17]. В настоящей работе используется другой подход, основанный на дискретизации процесса движения во времени, естественным образом вытекающий из процедуры дискретизации по времени уравнений движения.
Вероятность того, что время t между столкновениями лежит в интервале (0,т), определяет функцию распределения случайной величины х — времени свободного пробега между столкновениями:
Р(0 < t <х) = F(х) = 1 - ехр { )dt },
(23)
где N — три независимые реализации случайной
нормально распределенной величины с нулевым средним и единичной дисперсией. Такая процедура соответствует максвелловскому распределению молекул газа по скоростям при температуре Т.
Полагая, что М задана в атомных единицах массы, Т — в градусах Кельвина, а скорость — в см/с, получаем
(22)
Следующим важным этапом моделирования столкновения является розыгрыш времени или длины свободного пробега. Случайная величина — число столкновений иона с молекулами газа за определенный промежуток времени или на определенном отрезке длины траектории (в зависимости от используемой параметризации) представляет собой пуассоновский процесс, вообще говоря, неоднородный для неоднородного ЭГД-поля.
где V = па(»)» — средняя частота столкновений; п — концентрация молекул или атомов среды; и — скорость относительного движения сталкивающихся частиц; а (и) — сечение столкновения, зависящее от относительной скорости. Интегрирование в соотношении (23) проводится вдоль траектории иона, вдоль которой изменяются, вообще говоря, концентрация молекул и скорость иона » .
Вероятность того, что расстояние s между столкновениями лежит в интервале (0,Л), определяет функцию распределения случайной величины Л — длины свободного пробега между столкновениями:
Р(0 < s < Л) = F(Л) = 1 - ехр {-|щ}, (24)
где Л = 1/(па (и)) — средняя длина свободного пробега. Интегрирование в соотношении (24) также проводится вдоль траектории иона.
При вычислении средней частоты столкновений или средней длины свободного пробега в качестве сечения столкновения а (и) в соотношениях (23) и (24) необходимо брать наибольшее по величине сечение аг или ар .
Соотношение (16) может быть переписано в виде
F(х) = R или F(Л) = R , (25)
в котором функции распределения даются формулами (23) и (24).
Для генерирования случайных величин т и Л необходимо разрешить уравнения (25) относительно случайных величин т и Л. Очевидно, это возможно в случаях, когда постоянны среднее время между столкновениями т = 1/у (постоянна средняя частота столкновений) или средняя длина свободного пробега:
т = const или X = const. И мы получаем в первом случае
т = -т ln R
и во втором
X = -Xln R .
(26)
(27)
(28)
As << П |Vn| ,
эту концентрацию можно считать постоянной.
Если, кроме того, А << Л , из соотношения (24) следует, что вероятность столкновения на этом промежутке
P(As) = As/ X.
(30)
Если оперировать не длиной промежутка, а временем движения по нему А и если выполняется условие Ах << т, то из соотношения (23) следует, что вероятность столкновения на этом промежутке
Р(А) = А/ т. (31)
В развитой нами процедуре обе модели — и модель Максвелла, и модель абсолютно жестких шаров — объединены в единый алгоритм. Расчет траектории иона производится на малом промежутке как по времени, так и по расстоянию. На данном промежутке берется локальное значение концентрации молекул газа и вычисляются значения среднего времени т и средней длины пробега Л между столкновениями. Далее из двух значений вероятности столкновения (30) или (31) выбирается наибольшая
Pcoll = max(P(As), P(At)),
(32)
В случае постоянной концентрации молекул газа первый случай (в литературе он получил название "модель Максвелла") реализуется при условии (г(и)и = const. Этому условию удовлетворяет сечение поляризационного захвата (5), которое адекватно описывает рассеяние ионов на нейтральных, не обладающих дипольным моментом молекулах при малых энергиях (от тепловых до нескольких эВ).
В случае постоянной концентрации молекул газа второй случай — постоянна средняя длина свободного пробега — реализуется при условии ст(и) = const, т. е. при условии независимости сечения взаимодействия от энергии иона (в литературе он получил название "модель абсолютно упругих шаров"). Это условие удовлетворительно выполняется, когда столкновение определяется близкими отталкивательными силами (отталкива-тельной ветвью потенциала взаимодействия).
В реальной ситуации концентрация молекул газа может иметь большие градиенты. Однако если рассматривать движение иона на достаточно малых промежутках As, на которых выполняется условие
и разыгрывается столкновение с вероятностью попадания равномерно распределенной на (0,1) случайной величины R в промежуток (0, Рсо11).
РАСЧЕТ ТРАЕКТОРИИ
В соответствии с принятой моделью в промежутке между столкновениями ион движется только под действием электрического поля. Уравнение движения в терминах безразмерного времени д = 0Х, 0=2п/, где / — частота внешнего электрического поля:
d^r
d?2
= A(r,?).
(33)
В частности, для квадрупольного радиочастотного потенциала в плоскости (х, у) и постоянного во времени однородного электрического поля, направленного вдоль оси X :
Ф(г,X) = (и0 - К^О))(х2 - у2)/ г2 -Е,г, (34)
= Ъ - ^
x | 4 2
= | ^ - q^ 42
cos? Iy,
(35)
(29) A =
q
mQ2
-E.
где
a2 = 8
qUо mQ2 r02
q2 = 4
mQ2r02
(36)
есть стандартные обозначения параметров в квад-рупольной масс-спектрометрии. Переменные в уравнении (33) для потенциала (34) разделяются.
Динамика ионов в переменном во времени квадрупольном потенциале в бесстолкновитель-ном режиме полностью определяется параметрами а2 и q2 (36), а не отдельными значениями величин т, О, г0, и0 и У0. Эти законы подобия нару-
x
шаются, если в рассмотрение вводятся столкновения с молекулами газа. В этом случае появляются еще два масштаба: средняя длина свободного пробега Л и средняя частота столкновений V, что приводит к дополнительным параметрам задачи
1 = Л/г
Г=V//.
(37)
г, = г + vh + А Л /2,
п+1 п п п '
Ап+1 = А(Гп+1, ?п+1),
и„+1 = ип + (А„ + А„+1)/2.
(38)
h □ О/у .
(39)
Ф = = 1 q у —
Фе(Г тО2 г04 4 Ъ ° r02,
г2 = х2 + у2,
Первый из этих параметров характерен для задач, в которых преобладают столкновения при достаточно больших энергиях, когда рассеяние определяется потенциалом "жесткой сердцевины" и определенное значение имеет средняя длина свободного пробега. Второй параметр характерен для задач, в которых преобладают столкновения при малых энергиях, когда рассеяние определяется поляризационным потенциалом и определенное значение имеет частота столкновений.
Уравнение движения (33) в работе интегрировалось методом Верле в скоростной форме [18]: для (п+ 1)-го шага по времени дп+1 = (п + 1)h имеем:
который является потенциалом гармонического осциллятора. Частота колебаний F в поле этого потенциала пропорциональна частоте внешнего поля /
F =
242
/.
(41)
Модуль динамической составляющей поперечного поля для квадрупольного потенциала не зависит от азимутального угла и дается простым соотношением
г
Е = 2У —
^ 2 •
(42)
Параметром q2 определяется амплитуда быстрых осцилляций
qED 1
аО =—^ =_ q2 г. О тО2 2 2
(43)
Шаг интегрирования выбирается из двух условий.
1) Достаточное число точек дискретизации на временном периоде — h << 2п .
2) Мала вероятность столкновения на выбранном шаге P(Дt) □ 1. Откуда следует
За время х между столкновениями ион сдвинется в поперечной плоскости на
Л=
Е.х'1
т2
Отношение амплитуды быстрых осцилляций к величине этого сдвига
Поскольку типичные значения О □ □ 106рад • с-1, а, например, при давлении буферного газа 1 Торр V □ 107сч, то это условие является довольно жестким. Кроме того, при выборе шага интегрирования необходимо учесть условие (29), накладываемое газодинамическими параметрами задачи.
Алгоритм Верле обладает вторым порядком точности при минимальных затратах машинного времени, что очень важно при статистическом моделировании физических процессов.
Обратимся теперь к рассмотрению движения в квадрупольном чисто динамическом поле без статической компоненты (и0 = 0). При условии q2 □ 1 (адиабатическое приближение) движение заряженной частицы в поперечной плоскости может быть представлено в виде быстрых осцилля-ций с частотой / , наложенных на синусоидальные колебания в поле эффективного потенциала [4]
ап_ 1 *2_ 1
Л 2п2 /2 2п2
г2-
(44)
Таким образом, при V < / имеем Л > аО — одно столкновение иона с молекулами буферного газа происходит в среднем за несколько быстрых осцилляций. Наоборот, при V > / Л < аО — за период быстрых осцилляций происходит несколько столкновений. Интересна переходная область, когда V « /, и можно ожидать появления особенностей, возможно, резонансных в свойствах рассматриваемой системы. Данный вопрос требует дополнительного исследования. Заметим, что движение ионов в ЭГД-полях носит дробный характер. Оно состоит в чередовании промежутков времени свободного (без столкновений) движения в электрическом поле и промежутков времени взаимодействия (столкновений) с молекулами, которые на много порядков меньше первых, т. е. могут рассматриваться точечными во временном масштабе. Именно этим обстоятельством можно
0
0
Результаты моделирования Линейная аппроксимация
Результаты моделирования Квадратичная аппроксимация
Е^ В/см
Рис. 2. Дрейфовая скорость ий ионов в случае поляризационного взаимодействия для т = 28, М = 28 в зависимости от величины напряженности продольного поля Е
Е^ В/см
Рис. 3. Средняя энергия W ионов в случае поляризационного взаимодействия для т = 28, М = 28 в зависимости от величины напряженности продольного поля Е
объяснить дополнительный разогрев ионов в радиочастотном поле, который связан со случайным характером распределения моментов времени столкновений ионов с молекулами в пределах периода изменения внешнего поля.
ТЕСТИРОВАНИЕ
Тестирование процедуры моделирования движения ионов в газе под действием электрического поля проводилось как с точки зрения столкнови-тельной динамики, так и с точки зрения динамики в переменном электрическом поле.
В первом случае моделировался дрейф ионов в газе молекул N (М = 28 а.е., поляризуемость ар = 1.76 А3) под действием продольного однородного и постоянного во времени электрического поля: Е5 = (0,0,Е2), Е0 = 0. Варьировалась величина продольной составляющей поля Ег и вычислялись дрейфовая скорость ионов ил =(о2) и соответственно коэффициент подвижности k =
. 1 2 = °л/Ег, средняя энергия ионов Ж =— т(и ) ,
а также коэффициенты поперечной и продольной диффузии
D±=< х2 + у2) / ( 4Лй ), Da =<( г-< г))2) / ( Ъй).
(45)
Здесь и далее символ <•••) означает среднее по ансамблю частиц, полученное в результате компьютерного моделирования в фиксированный момент времени дрейфа . Отметим, что результаты моделирования с хорошей точностью удовлетворяют условию эргодичности: средняя по времени z-компонента скорости отдельного иона (скорость дрейфа) равна средней скорости по ансамблю
У ).
Результаты представлены на рис. 2, 3 и в таблице для поляризационного потенциала и на рис. 4, 5 для потенциала абсолютно упругих шаров.
В случае чисто поляризационного потенциала ряд вычисленных значений дрейфовой скорости с высокой точностью (погрешность менее 1%) аппроксимируется прямой (см. рис. 2), и, следовательно, коэффициент подвижности не зависит от напряженности электрического поля. Его значение с большой точностью соответствует известному выражению для подвижности в случае поляризационного потенциала [7]
k =
ЛV
(46)
где, как и выше, q — заряд иона, л — приведенная масса, V — средняя частота столкновений. Наилучшая аппроксимация достигается при таких значениях шага дискретизации уравнений движения по времени, при которых вероятность столкновения на данном временном шаге менее 0.1.
Значения параметров задачи, полученные в результате компьютерного эксперимента по моделированию движения ионов с массами т = 1, 28 и 133 а.е.м. в газе молекул N2 при р = 1 мм рт. ст. и Т = = 300 К (концентрация молекул п = 3.22*1016 см3)
Характеристика т = 1 т = 28 т = 133
Ег, В/см 2 10 20 2 10 20 2 10 20
102 X см2/с 2.50 9.75 31.2 0.566 0.805 1.47 0.480 0.620 0.985
Д, 102 х см2/с 2.63 9.72 28.3 0.640 0.867 1.49 0.507 0.633 0.963
В>7, 102 х см2/с 2.52 9.06 29.5 0.573 0.680 0.933 0.458 0.499 0.604
В, 102 х см7с 2.56 9.07 27.5 0.628 0.701 0.974 0.485 0.517 0.619
к, 103 х см2/(В ■ с) 8.87 9.03 8.96 2.32 2.39 2.39 1.85 1.88 1.87
Результаты моделирования
Ег, В/см
Рис. 4. Дрейфовая скорость ил ионов в случае потенциала взаимодействии упругих шаров для т = = 28, М = 28 в зависимости от величины напряженности продольного поля Е
Ш
О 4 -
Результаты моделирования -Линейная аппроксимация
10 12 14
Ег, В/см
Рис. 5. Средняя энергия W ионов в случае потенциала взаимодействии упругих шаров для т = 28, М = 28 в зависимости от величины напряженности продольного поля Е
Значения средней энергии W в зависимости от напряженности электрического поля для чисто поляризационного потенциала взаимодействия с высокой точностью (погрешность менее 1 %) допускают квадратичную аппроксимацию (см. рис. 3) и соответствуют соотношению Ванье для средней энергии дрейфующих ионов [11]
1 3 1
Ж № =- т(и2 >=- кТ + - (т + М )и2.
2 2 2
(47)
В таблице представлены рассчитанные по результатам компьютерного моделирования в соответствии с формулами (45) коэффициенты поперечной и продольной диффузии для ряда масс ионов и значений напряженности электрического поля. Для сравнения там же приведены значения
коэффициентов диффузии, рассчитанных в соответствии с соотношениями Ванье [11]:
± \ х / ^ \ ху / ?
а 2а
= ^ (<^2 >-^2 ) = ^ -ц,)2 >.
(48)
Здесь а = qEJm — ускорение иона под действием электрического поля. Наши расчеты соответствуют этим выражениям в пределах относительной погрешности, меньшей 10 %.
В таблице приведены также полученные в результате компьютерного моделирования значения коэффициентов подвижности. Экспериментальные значения коэффициентов подвижности и диффузии в пределе нулевого поля [11], пересчитанные
для р = 1ммрт.ст., составляют: 1) для N2 в N — k = 0.152х 104см2/(В• с) и В = 39.7см2/с; 2) для Cs+ в N2 — k = 0.185 х104см2/(В• с) и В = = 47.7см2/с.
Для тяжелых ионов соответствие с экспериментом отличное. Для ионов N2 , дрейфующих в N, в реальном эксперименте присутствует резонансная перезарядка. В нашем компьютерном эксперименте резонансная перезарядка не моделировалась, что, естественно, приводит к более высоким значениям рассчитанных коэффициентов подвижности и диффузии по сравнению с экспериментальными данными. Однако наши значения близки к имеющимся в литературе (см., например, [11]) экспериментальным данным для близких по массе ионам, которые не испытывают резонансной перезарядки при столкновении с молекулой N. Например, подвижность № + в N при р = 1ммрт.ст. составляет 0.22 х104см2/(В • с).
Зависимости дрейфовой скорости и средней энергии ионов от величины электрического поля в случае потенциала взаимодействия абсолютно упругих шаров представлены на рис. 4 и 5. Средняя энергия дрейфующих ионов линейно зависит от Ег, дрейфовая скорость растет пропорционально ^Е^ , соответственно коэффициент подвижности убывает обратно пропорционально ^/Е^ в полном соответствии с теоретическими представлениями [11].
Комбинированная модель взаимодействия — отталкивание на малых расстояниях в соответствии с потенциалом (2) в приближении эффективного радиуса взаимодействия (3) и притяжение на больших расстояниях в соответствии с поляризационным потенциалом (4) в приближении Лан-жевена (5) — тестировалась на примере моделирования дрейфа ионов с массами т = 28, 133 и 1000 а.е.м. в газе молекул Для этих ионов константа скорости поляризационного захвата стри при столкновении с молекулой N составляет соответственно 8.26 х10-10, 6.43 х10-10 и 5.82х10-10см3/с. Для давления р = 1 Торр (концентрация молекул п = 3.2 х 1016см-3) формула (46) дает следующие значения подвижности: k = 2.3 х 103, 1.8 х 103 и 1.6 х 103 см2 / (В • с) соответственно. Наши расчеты дали k = = 2.7 х 103см2 / (В • с) для т = 28 при Е2 < < 100В/см; k = 2.1 х 103см2 /(В• с) для т = 133 при Е2 < 70 В / см и k = 1.7 х 103см2 / (В • с) для т = 1000 при Ег < 40В/см. При больших значе-
ниях напряженности электрического поля, как показывает контроль характера столкновений, предусмотренный в программе, доля поляризационных столкновений начинает убывать, растет доля столкновений, обусловленных отталкивательной ветвью потенциала взаимодействия. Наблюдается заметное отклонение вычисленных значений подвижности от значений, даваемых формулой (46), — подвижность начинает убывать обратно пропорционально уЕ .
Таким образом, предложенная процедура дает значения подвижности с относительной погрешностью, меньшей 20 %. Для чисто поляризационного взаимодействия в области малых напряжен-ностей электрического поля погрешность меньше 1 %.
С точки зрения динамики иона в высокочастотном электрическом поле тестирование предложенной процедуры и ее программной реализации проводилось на основании исследования устойчивости движения иона в квадрупольном радиочастотном поле.
На рис. 6 представлены траектории ионов, стартующих из одной точки в поперечной плоскости квадрупольного радиочастотного поля. Ионы дрейфуют под действием продольного постоянного поля в буферном газе. Отчетливо наблюдается фокусировка ионов — стягивание к оси квадру-польного поля.
ЗАКЛЮЧЕНИЕ
Разработаны и протестированы процедуры компьютерного моделирования движения ионов в неоднородных переменных во времени электрогазодинамических полях. Процедуры основаны на комбинированной модели взаимодействия ионов с молекулами газа, описывающей рассеяние ионов в широкой области энергий и учитывающей единым образом взаимодействие иона с молекулой как в области отталкивания в приближении эффективного радиуса, так и в области притяжения в приближении Ланжевена для поляризационного потенциала. Тестирование процедур проводилось на основе вычисления коэффициентов подвижности и диффузии. Полученные в результате компьютерного моделирования значения коэффициентов подвижности при малых напряженностях электрического поля удовлетворяют соотношению Ванье для чисто поляризационного потенциала с относительной погрешностью, меньшей 1 %, значения коэффициентов диффузии — с относительной погрешностью, меньшей 10 %. Рассчитанные значения удовлетворительно совпадают с имеющимися экспериментальными данными.
Рис. 6. Траектории ионов в квадру-польном поле (и0 = 0, ¥0 = 598 В) для т = 133, М = 28, п = = 0.322х1016 см3 (р = 0.1 Торр, V = = 0.207х107 с-1).
Приложено продольное поле Ег = = 20 В/см. Область дрейфа 0 < г < < 10 см, г0 = 2 см, точка старта г = = 1 см. Варьируется параметр q2 (36)
Правильно передается переход от поляризационного потенциала (независимость коэффициентов подвижности от напряженности электрического поля) к потенциалу отталкивания (убывание
обратно пропорционально ).
Разработанные модели и соответствующие вычислительные процедуры позволяют вводить неизотропные модели рассеяния (анизотропный потенциал взаимодействия) и применимы при расчете траекторий движения ионов в неоднородных газовых потоках и в переменных во времени неоднородных электрических полях.
СПИСОК ЛИТЕРАТУРЫ
Douglas D.J., French J.B. Collisional focusing effects in radio frequency quadrupoles // J. Am. Soc. Mass Spectram. 1992. Vol. 3. P. 398-408. Doi: 10.1016/1044-0305(92)87067-9.
Tolmachev A.V., Chernushevich I.V., Dodonov A.F., Standing K.G. A collisional focusing ion guide for coupling an atmospheric pressure ion source to a mass spectrometer // Nucl. Instr. Meth. Phys. Res. B. 1997. Vol. 124. P. 112-119. Doi: 10.1016/S0168-583X(97)00068-2. Schwarz S., Bollen G., Lawton D., Lofy P. et al. The low-energy-beam and ion-trap facility at NSCL/MSU // Nucl. Instr. Meth. Phys. Res. B. 2003. Vol. 204. P. 507-511. Doi: 10.1016/S0168-583X(97)00068-2.
4. Gerlich D. Inhomogeneous RF fields: a versatile tool for the study of processes with slow ions // Adv. Chem. Phys. 1998. Vol. 82. P. 1-176. Doi: 10.1002/ 9780470141397.ch1.
5. Reid G.E., Wells J.M., Badman E.R., McLuckey S.A. Performance of a quadrupole ion trap mass spectrometer adapted for ion/ion reaction studies // Int. J. Mass Spec-trom. 2003. Vol. 222. P. 243-258. Doi: 10.1016/S1387-3806(02)00984-3.
6. Tolmachev A.V., Udseth H.R., Smith R.D. Modeling the ion density distribution in collisional cooling RF multipole ion guides // Int. J. Mass Spectrom. 2003. Vol. 222. P. 155-174. Doi: 10.1016/S1387-3806(02)00960-0.
7. Щербаков А.П. Численное моделирование транспортировки ионных пучков в электрогазодинамических полях // Научное приборостроение. Л.: Наука, 1988. С. 46-55.
8. Бородинов А.Г., Щербаков А.П. О моделировании процесса формирования ионного пучка в газодинамических источниках ионов // Научное приборостроение. Л.: Наука, 1990. С. 10-14.
9. Бородинов А.Г., Веренчиков А.Н., Щербаков А.П. Исследование транспортировки ионов в газодинамических системах // ЖТФ. 1991. Т. 61, № 6. С. 1-7.
10. Бородинов А.Г., Огородников А.К., Щербаков А.П. Исследование транспортировки ионов в газодинамических источниках и оптимизация системы формирования и фокусировки пучка // Научное приборостроение. 1993. Т. 3, № 2. С. 79-85.
11. Мак-Даниель И., Мэзон Э. Подвижность и диффузия ионов в газах. М.: Мир, 1976. 422 с.
12. Леонас В.Б. Исследование короткодействующих межмолекулярных сил // УФН. 1972. Т. 107, № 1. С. 29-56.
13. Мак-Даниель И. Процессы столкновений в ионизованных газах. М.: Мир, 1967. 832 с.
14. Ландау Л.Д., Лифшиц Е.М. Механика. М.: Наука, 1973. 208 с.
15. Ермаков С.М. Метод Монте-Карло и смежные вопросы. Изд. 2-е. М.: ГРФМЛ "Наука", 1975. 472 с.
16. Щербаков А.П. Методы компьютерного моделирования процессов атомного рассеяния в задачах научного приборостроения // Научное приборостроение. 2003. Т. 13, № 1. С. 14-23. URL: http://213.170.69.26/mag/ 2003/full1/Art2.pdf.
17. Феллер В. Введение в теорию вероятностей и ее приложения. В 2-х томах. Т. 2. Пер. с англ. М.: Мир, 1984. 738 с.
18. Хеерман Д.В. Методы компьютерного эксперимента в теоретической физике. М.: Наука, 1990. 176 с.
Институт аналитического приборостроения РАН, г. Санкт-Петербург
Контакты: Щербаков Анатолий Петрович, anpshch@yandex. ru
Материал поступил в редакцию: 14.10.2016
ISSN 0868-5886
NAUCHNOE PRIBOROSTROENIE, 2016, Vol. 26, No. 4, pp. 64-76
SIMULATION OF ION-MOLECULE COLLISIONS IN INHOMOGENEOUS TIME-DEPENDED ELECTROGASDYNAMIC FIELDS
A. P. Scherbakov
Institute for Analytical Instrumentation of RAS, Saint-Petersburg, Russia
Computer based procedures have been developed and tested to simulate the ion motion in inhomogeneous time-depended electrogasdynamic fields. Procedures are based on combined model of ion-molecule collision. This combined approach takes into account both repulsive and attractive region of the interaction potential and correctly describes the ion scattering over a wide range of collision energies. Models and computational procedures developed enable to introduce the anisotropic modes of scattering. Monte Carlo simulation of free path length ore free time is based on time discretization technique, which is well combined with discretization of motion equations.
Keywords: collision cross section, mobility, diffusion, polarization potential, hard sphere potential, electrogasdynamic field
REFERENCES
1. Douglas D.J., French J.B. Collisional focusing effects in radio frequency quadrupoles. J. Am. Soc. Mass Spec-trom., 1992, vol. 3, pp. 398-408. Doi: 10.1016/1044-0305(92)87067-9.
2. Tolmachev A.V., Chernushevich I.V., Dodonov A.F., Standing K.G. A collisional focusing ion guide for coupling an atmospheric pressure ion source to a mass spectrometer. Nucl. Instr. Meth. Phys. Res. B., 1997, vol. 124, pp. 112-119. Doi: 10.1016/S0168-583X(97)00068-2.
3. Schwarz S., Bollen G., Lawton D., Lofy P. Morrissey D.J., Ottarson J., Ringle R., Schury P., Sun T., Varentsov V., Weissman L. The low-energy-beam and ion-trap facility at NSCL/MSU. Nucl. Instr. Meth. Phys. Res. B., 2003, vol. 204, pp. 507-511. Doi: 10.1016/S0168-583X(97)00068-2.
4. Gerlich D. Inhomogeneous RF fields: a versatile tool for the study of processes with slow ions. Adv. Chem. Phys., 1998, vol. 82, pp. 1-176. Doi: 10.1002/ 9780470141397.ch1.
5. Reid G.E., Wells J.M., Badman E.R., McLuckey S.A. Performance of a quadrupole ion trap mass spectrometer adapted for ion/ion reaction studies. Int. J. Mass Spec-trom, 2003, vol. 222, pp. 243-258. Doi: 10.1016/S1387-3806(02)00984-3.
6. Tolmachev A.V., Udseth H.R., Smith R.D. Modeling the ion density distribution in collisional cooling RF multipole ion guides. Int. J. Mass Spectrom., 2003, vol. 222, pp. 155-174. Doi: 10.1016/S1387-3806(02)00960-0.
7. Scherbakov A.P. [Numerical modeling of transportation of ionic bunches in electrogasdynamic fields]. Nauchnoe Priborostroenie [Scientific Instrumentation], Leningrad, Nauka Publ., 1988, pp. 46-55. (In Russ.).
8. Borodinov A.G., Scherbakov A.P. [About modeling of process of formation of an ionic bunch in gasdynamic sources of ions]. Nauchnoe Priborostroenie [Scientific Instrumentation], Leningrad, Nauka Publ., 1990, pp. 10-
14. (In Russ.).
9. Borodinov A.G., Verenchikov A.N., Scherbakov A.P. [Research of transportation of ions in gasdynamic systems]. Zhurnal tekhnicheskoj fiziki [Journal of technical physics], 1991, vol. 61, no. 6, pp. 1-7. (In Russ.).
10. Borodinov A.G., Ogorodnikov A.K., Scherbakov A.P. [Research of transportation of ions in gasdynamic sources and optimization of system of formation and focusing of a bunch]. Nauchnoe Priborostroenie [Scientific Instrumentation], 1993, vol. 3, no. 2, pp. 79-85. (In Russ.).
11. Mason E., McDaniel E.W. Transport Properties of Ions in Gases. Wiley, 1976 (Russ. ed.: Mack-Daniel' I., Mezon E. Podvizhnost' i diffuziya ionov v gazah. Moscow, Mir Publ., 1976. 422 p.).
12. Leonas V.B. [Research of short-range intermolecular forces]. UFN [Advances in Physical Sciences], 1972, vol. 107, no. 1, pp. 29-56. (In Russ.).
13. McDaniel E. Collision Phenomena in Ionized Gases. New York, Wiley, 1964 (Russ. ed.: Mack-Daniel' I. Processy stolknovenij v ionizovannyh gazah. Moscow, Mir Publ., 1967. 832 p.).
14. Landau L.D., Lifshic E.M. Mekhanika [Mechanic]. Moscow, Nauka Publ., 1973. 208 p. (In Russ.).
15. Ermakov S.M. Metod Monte-Karlo i smezhnye voprosy [Monte Carlo method and related issues]. Second Edition. Moscow, Nauka Publ., 1975. 472 p. (In Russ.).
16. Scherbakov A.P. [Computer simulation of atom scattering as applied to problems of scientific instrumentation]. Nauchnoe Priborostroenie [Scientific Instrumentation], 2003, vol. 13, no. 1, pp. 14-23. URL: http:// 213.170.69.26/mag/2003/full1/Art2.pdf. (In Russ.).
17. Feller W. An Introduction to Probability Theory and Its Applications. Vol. 2. Second Edition. Wiley, 1983. 704 p. (Russ. ed. : Feller V. Vvedenie v teoriyu veroyatnostej i ee prilozheniya. Vol. 2. Translation from English. Moscow, Mir Publ., 1984. 738 p.).
18. Heerman D.W. Computer simulations methods in theoret-
ical physics. Springer-Verlag, 1986 (Russ. ed.: Heerman D.V. Metody komp'yuternogo ehksperimenta v teoreti-cheskoj fizike. Moscow, Nauka Publ., 1990. 176 p.).
Contacts: Scherbakov Anatoliy Petrovich, anpshch@yandex. ru
Article received in edition: 14.10.2016
