CC BY
13
УДК 621.385.833
ПРИМЕНИМОСТЬ БЕССЕТОЧНОГО МЕТОДА ГИДРОДИНАМИКИ СГЛАЖЕННЫХ ЧАСТИЦ ДЛЯ МОДЕЛИРОВАНИЯ ПРОЦЕССА ИЗОЛЯЦИИ ЗОНДОВ ЭЛЕКТРОХИМИЧЕСКОГО СТМ
ТЮРИКОВ А. В., ШЕЛКОВНИКОВ Е. Ю., ЖУЙКОВ Б. Л., КИЗНЕРЦЕВ С. Р. Институт механики Уральского отделения РАН, 426067, г. Ижевск, ул. Т. Барамзиной, 34
АННОТАЦИЯ. В статье рассмотрены вопросы применения бессеточного метода гидродинамики сглаженных частиц для моделирования процесса нанесения изолирующего покрытия апиезона Ь на острие зонда электрохимического сканирующего туннельного микроскопа. Описаны особенности метода и отличие от классических сеточных методов, показаны его преимущества при моделировании жидкостей со свободными границами и переменными вязкостью и поверхностным натяжением. Разработана специализированная программа для проведения компьютерного эксперимента, использующая для ускорения расчетов современный графический адаптер и технологии БтесЮошр^е.
КЛЮЧЕВЫЕ СЛОВА: электрохимический сканирующий туннельный микроскоп, зондирующее острие, изоляция, бессеточный метод сглаженных частиц.
ВВЕДЕНИЕ
В настоящее время вычислительная гидродинамика получила широкое распространение и применяется для решения самых разнообразных задач, связанных с изучением движения жидкостей в различных областях науки и техники. Классические методы гидродинамики, использующие математические уравнения для описания движения жидкости, в основном базируются на сеточных вычислениях [1]. Однако, в некоторых случаях, например, при решении задач моделирования процессов в системах с изменяющимися границами и свободной поверхностью, решение с помощью классических методов не является оптимальным. К таким задачам можно отнести нанесение изолирующего материала на острие зондов электрохимического сканирующего туннельного микроскопа (ЭСТМ), который является незаменимым средством исследования поверхностей образцов, находящихся в растворах электролитов, с разрешением близким к атомарному. Особенностью такого микроскопа является дополнительная обработка зондов, заключающаяся в их изоляции непроводящим вязким текучим материалом для предотвращения протекания с их поверхности паразитных фарадеевских токов, которые затрудняют детектирование туннельного тока и, как следствие, формирование качественных СТМ-топографий. При этом кончик острия зонда обязательно должен оставаться неизолированным и быть по возможности минимальным для осуществления протекания и регистрации туннельного тока, а также, одновременно, минимизирующим паразитные фарадеевские токи. Поэтому процесс изоляции зондов ЭСТМ является не менее важным этапом при их формировании, чем изготовление самих острий с нановыступами на конце [2 - 4] для получения качественных СТМ-топографий исследуемой поверхности.
МАТЕМАТИЧЕСКАЯ МОДЕЛЬ
Основными физическими характеристиками изотермической вязкой жидкости являются скорость V, плотность р и давление р, которые являются непрерывными величинами в пределах жидкости. Классическая формулировка движения для несжимаемого потока жидкости в течение времени I дается уравнениями Навье-Стокса [1]:
^+(V = Луу = о, (1)
где т - вязкость; /ех1 - сумма внешних сил, действующих на жидкость (в том числе, и гравитация). В классической эйлеровой гидродинамике жидкость представляется составленной из жидких ячеек, заключенных в регулярную сетку, которая обеспечивает численное решение уравнений (1) одним из хорошо известных методов (например, методом конечных разностей). Несмотря на то, что сеточный метод Эйлера дает хорошее описание некоторых характеристик жидкости (например, плотности и поля давления), его главным недостатком является наличие самой сетки. Жидкость при моделировании вынуждена существовать с сеткой и не может существовать вне области без сетки. Таким образом, есть существенная проблема описания естественного потока жидкости, границы которого постоянно меняются. Также, одним из существенных недостатков сеточных методов является значительное динамическое потребление памяти ЭВМ (особенно это касается моделей с адаптивными сетками, использование которых призвано компенсировать подвижные границы моделируемой системы). Задача моделирования процесса изоляции острий ЭСТМ связана именно с таким течением жидкости с изменяемой вязкостью (что связано с ее застыванием). Это позволяет утверждать, что применение сеточных методов для создания модели процесса изоляции острий не является оптимальным.
Формулировка метода гидродинамики сглаженных частиц (БРИ) (называемого также лагранжевой гидродинамикой [5 - 7]) происходит из области вычислительной астрофизики [6] и разработана для описания сжимаемых потоков. Метод широко использовался для моделирования астрофизических явлений, где сложные задачи могли быть описаны и поняты более интуитивно. БРИ является интерполяционным методом для расчета непрерывных параметров потока жидкости и их производных. При этом используются дискретные части жидкости, под которыми понимаются сглаженные частицы, имеющие конкретные параметры (массу, координаты, скорость и т.д.). Множество этих частиц также несет информацию о полях физических величин, связанных с решаемой проблемой, например, плотности, температуре, давлении и т.д. Описанные характеристики можно получить как взвешенные средние, вычисленные по соседним частицам. По сравнению с другими известными методами для численного приближения производных (например, методом конечных разностей, который требует, чтобы частицы были закреплены на регулярной сетке) БРИ позволяет вычислять производные непрерывных величин с применением аналитического дифференцирования. Предполагается, что каждая частица занимает часть исследуемого пространства, поэтому для расчета более точных взвешенных средних величин дискретные частицы должны располагаться как можно более плотно.
Интерполяция в БРИ основана на теории интегральных интерполянтов, использующих сглаживающие ядра. Так интегральная интерполяция величины функции А (г), определяемая по области О, вычисляется как:
А1 (г ) = | А (г ')Ж (г - г ',И)с1г', (2)
О
где г - радиус-вектор точки в О ; Ж - сглаживающее ядро радиуса И . Радиус (или ширина) ядра является масштабирующим фактором, контролирующим сглаживающие характеристики ядра. Дискретным эквивалентом (2) является выражение, применяемое при численных расчетах:
As (г) = Е АУЖ (г - г, И), (3)
]
где суммирование ведется по частицам в области О ; V - объем, занимаемый ] -й частицей; rj - ее радиус-вектор; А}- = А (г}-). Полагая, что V = mi / р{ (где и р{ - масса и плотность
соответствующей частицы), интерполяция величины А, а также ее градиент и лапласиан записываются как:
_ ш _ ш
А5 (г ) = 2 А} (Г - т], И); УА5 (г ) = 2 А] Р УЖ (г - т], И);
1 р 1 Г (4)
АА8 (г ) = 2 А; Г АЖ (г - г}, И).
1 р1
Более точное выражение для градиента (обладающее, к тому же, необходимым для некоторых взаимодействий между частицами свойством симметрии) можно записать в виде [6]:
( А А Л
УА, (г) = р2 -2+-Г (г -Т},И). (5)
1 IР} Р у
Для поля силы, не обладающей свойствами симметрии (такой, например, как вязкое трение) также применимо:
^ (г)=-2(1)туЖ(г -т],И). (6)
р 1
Сглаживающее ядро, согласно [7], должно обладать следующими свойствами:
| Ж (Г, И) с1г = 1; Ж (Г, И )> 0; 1пп Ж (г, И ) = #( г), (7)
п
где б(г) - дельта-функция Дирака.
Уравнения Навье-Стокса (1) (описывающие эйлерову гидродинамику) представляют собой условие баланса сил в ячейках конечно-разностной сетки и включают в себя
конвективное слагаемое (V -У) V, в котором нет необходимости в случае лагранжева
подхода, поскольку он напрямую связан с описанием движения частиц жидкости. Поэтому уравнение баланса сил для каждой частицы запишется в гораздо более простом виде:
р± = /у-Ур + /ехГ (8)
аХ
В данное выражение входят следующие удельные силы (плотности сил): сила вязкостного трения = /Ау , сила давления /ргеб,6, = Ур и внешняя сила /ех(, в случае с жидкостью со свободной поверхностью образованная гравитационной составляющей и силой поверхностного натяжения /ех1 = + . Необходимо заметить, что сила поверхностного натяжения действует лишь на частицы, находящиеся близко к поверхности со стороны их «поверхностных» соседей. Поскольку / = р§, а = -окп
(где § - ускорение свободного падения; о- поверхностное натяжение; к- гауссова кривизна; п нормаль к поверхности), то выражение (8) можно записать как:
р— = /АУ -Ур + -окп. (9)
ах
Для определения силы поверхностного натяжения применяется поле величины (в литературе «цветное» поле [6]), определяемой как:
|1, внутри частицы; I 0, иначе.
При такой формулировке нормаль n = Vcj натяжения можно представить в виде:
fsurf S
Vc
VcAc
Vc
а k = Vn. Силу поверхностного
(11)
Ядра, применяемые для интерполяции вышеописанных величин, подробно описаны в [6] и обозначаются как Ж^, Ж, Жу^с, а выражения для компонент, входящих в (9),
можно записать как:
m
р(r ) = EmWf (r -r,h); Vc(r) = E—VWf (r -r,h);
j j Pj
f press (r ) =-— E
f n \
pL + P
2 2 V — Pj J
m
Wpress (r - , h); fgrav (rt) = —g;
m
/ ч mL , ч ^ Vc(г)_ mL , ч
fic(r) =mE (Vj - v) —AWvlSc (r - rj, h); fsurf (r) = -^^E —L AWf (r -rj, h);
(12)
dvi dt
■ a,
Vc (r )| l Pj
fpress + fgrav + fvisc + fsurf
Pi
Последнее выражение в (12) включает в себя расчет ускорения для каждой частицы. Для определения ее положения в пространстве на очередном временном шаге необходимо использовать одну из численных схем интегрирования, например, метод Верле, основанный на вычислении следующего положения частицы без определения скорости. Из разложения радиус-вектора частицы в ряд Тейлора в моменты времени t + & и t — &
^ ^ ^ aAt2 b, At3
rt+At = rt + vtAt + -+ - t
2
6
+ O(At4);
_ _ _ atAt2 btAt3 _.. 4.
rt-At = rt - v(At + ^---V + O(At4)
(13)
2
6
очевидно, что
П+At = 2rt - rt-At + atAt2 + O(At4).
(14)
РЕЗУЛЬТАТЫ И ИХ ОБСУЖДЕНИЕ
Модель представляла собой конус (кончик платиновой иглы), описываемый углом а и высотой Н, на который нанесен слой жидкости (апиезона) толщиной й. Остывание жидкости моделировалось следующим образом: задавались начальная и конечная температуры жидкости и время начала и конца ее остывания. Температура в течение времени изменялась линейно. В зависимости от температуры задавались константы вязкости и поверхностного натяжения. Вязкость определялась в соответствии с [8], исходя из вязкостных характеристик апиезона Ь. Поверхностное натяжение апиезона Ь было изучено в более ранних работах [9], результаты которых использованы в текущих расчетах. В таблице указаны параметры, общие для всех экспериментов.
Таблица
Параметры компьютерных экспериментов
Параметр Обозначение, ед. изм. Значение
Плотность р, кг/м3 896
Масса частицы m, пг 71,49
Радиус ядра h, мкм 8
Начальная температура T K start' 348
Конечная температура Tend , K 298
Угол конуса a, градусы 30
Высота конуса H, мкм 600
Толщина слоя жидкости d, мкм 40
Компьютерные эксперименты были выполнены с применением специализированной разработанной программы, использующей для ускорения расчетов современный графический адаптер и технологии DirectCompute. Анализ показал, что при проведении экспериментов на стенде с GPU-ускорителем GeForce GTX 960 прирост производительности по сравнению с версией на центральном процессоре Intel Core i7 6700 составил ~2 порядка. Результаты экспериментов иллюстрируются рисунком, на котором представлены модельные изображения кончика зондирующего острия ЭСТМ с изолирующим материалом при различных временных режимах его застывания.
а) б) в) г)
а) - частицы жидкости в начальный момент времени; б) - время застывания 50 мс; в) - 100 мс; г) - 300 мс
Рис. Модельные изображения кончика зондирующего острия ЭСТМ с изолирующим материалом при различных временных режимах его застывания
На рис. а показано расположение частиц изолирующей жидкости в начальный момент времени. На рис. б-г приведены изображения кончика острия с застывшей жидкостью при изменении температуры с Tstart до Tend за время 50 мс, 100 мс и 300 мс, соответственно. Видно, что интервал времени 50 мс (рис. б) недостаточен для «обнажения» кончика острия из-за слишком быстрого застывания апиезона; интервал 300 мс (рис. г) слишком велик в виду
излишнего «сползания» изолирующего материала, а близким к оптимальному режиму является интервал в 100 мс (рис. в), при котором неизолированный кончик имеет минимальную площадь.
Таким образом проведенные исследования показали целесообразность применения метода гидродинамики сглаженных частиц для моделирования процесса нанесения изолирующего покрытия апиезона Ь на острие зонда ЭСТМ. Результаты расчетов хорошо согласуются с экспериментальными данными и доказывают реалистичность применяемой модели в задаче об изоляции зондирующих острий ЭСТМ. При этом в дальнейших исследованиях перспективным представляется более точный учет распределения температуры в моделируемой системе, поскольку именно этот параметр является (в отношении заготовки нанозонда) наиболее поддающимся контролю и, одновременно, наиболее сильно влияющим на геометрические характеристики получаемого неизолированного кончика зонда ЭСТМ.
СПИСОК ЛИТЕРАТУРЫ
1. Полежаев В. И., Бунэ А. В., Верезуби А. В., Глушко Г. С., Грязнов В. Л., Дубовик К. Г., Никитин С. А., Простомолотов А. И. Математическое моделирование конвективного тепломассообмена на основе уравнений Навье-Стокса. М.: Наука, 1987. 272 с.
2. Тюриков А. В., Шелковников Е. Ю., Гуляев П. В., Жуйков Б. Л., Липанов С. И. Исследование влияния макроскопической формы перетравливаемой шейки заготовки СТМ-зонда на его нанотопологию // Химическая физика и мезоскопия. 2015. Т. 17, № 3. C. 482-487.
3. Тюриков А. В., Шелковников Е. Ю., Гуляев П. В., Жуйков Б. Л., Липанов С. И. Моделирование электрохимической стадии процесса формирования острий СТМ-зондов // Химическая физика и мезоскопия, 2016. Т. 18, № 2. С. 323-330.
4. Шелковников Е.Ю., Тюриков А.В., Кизнерцев С.Р., Липанов С.И. Методика моделирования процесса комбинированного травления нанозондов // Химическая физика и мезоскопия. 2012. Т. 14, № 4. С. 650-654.
5. Desbrun M., Cascuel M.-P. Smoothed Particles: A new paradigm for animating highly deformable bodies // Computer Animation and Simulation'96: Proceedings of the Eurographics Workshop in Poitiers, France, August 31-September 1, 1996. NewYork: Springer Wien, 1996, pp. 61-76. http://ru.bookzz.org/book/2118951/bf82ee
6. Monaghan J. J. Smoothed Particle Hydrodynamics // Annual Review of Astronomy and Astrophysics, 1992, vol. 30, pp. 543-574. DOI: 10.1146/annurev.aa.30.090192.002551
7. Muller M., Charypar D., Gross M. Particle-Based Fluid Simulation for Interactive Applications // Proceedings of 2003 ACM SIGGRAPH Symposium on Computer Animation, 2003. pp. 154-159. http://dx.doi.org/10.2312/SCA03/154-159
8. Фогельсон Р. Л., Лихачев Е. Р. Температурная зависимость вязкости // Журнал технической физики. 2001. Т. 71, вып. 8. C. 128-131.
9. Шелковников Е. Ю., Тюриков А. В., Гуляев П. В., Жуйков Б. Л., Липанов С. И. Особенности применения метода частиц для моделирования процесса изоляции зондирующих острий электрохимического СТМ // Химическая физика и мезоскопия. 2016. Т.18. №1. С. 154-159.
APPLICABITITY OF THE NONGRID SMOOTHED PARTICLES HYDRODYNAMICS METHOD FOR MODELLING THE PROCESS OF ISOLATION OF THE PROBES OF ELECTROCHEMICAL STM
Shelkovnikov E. Yu., Tyurikov A. V., Zhuikov B. L., Kiznertsev S. R.
Institute of Mechanics, Ural Branch of the Russian Academy of Sciences, Izhevsk, Russian
SUMMARY. In the paper questions of application of the nongrid smoothed particles hydrodynamics method for a process of modeling the insulation of the electrochemical STM probe tip with an apiezone L are considered. The methods of a hydrodynamics based on grid computation were popular within last decades, but in case of simulation of processes in systems with the changing boundaries and the free surface the decision based on a grid isn't optimum. Such the process, in particular, is insulation of probe tip for electrochemical STM which is irreplaceable means of a research of surfaces which are in solutions of electrolytes with resolution close to atomic. Feature of such microscope is the additional processing of probes consisting in their insulation by nonconducting viscous liquid material for preventing Faraday currents from their surface which complicate receiving atomic resolution on the STM-topographies. The small tip of the probe must remain uninsulated for a possibility of registration of a tunnel current. In the paper
features of smoothed particles hydrodynamics method and its difference from classical grid methods are described. Its advantages in case of simulation of liquids with the free boundaries and variable viscosity and surface tension are also shown for researches of processes of formation of probe tip of the electrochemical STM. The specialized program using the modern graphics adapter and DirectCompute technologies is developed for carrying out a computer experiment. The analysis showed that when carrying out the experiments on the bench with the GeForce GTX 960 GPU accelerator, productivity gain reached several of orders in comparison with the version on the central processor.
KEYWORDS: the electrochemical scanning tunnel microscope, the probe, an insulation, nongrid smoothed particles hydrodynamics method.
REFERENCES
1. Polezhaev V. I., Bune A. V., Verezubi A. V., Glushko G. S., Gryaznov V. L., Dubovik K. G., Nikitin S. A., Prostomolotov A. I. Matematicheskoe modelirovanie konvektivnogo teplomassoobmena na osnove uravneniy Nav'e-Stoksa [Mathematical modeling of convective heat and mass transfer on the basis of the Navier-Stokes equations]. Moscow: Nauka Publ., 1987. 272 p.
2. Tyurikov A. V., Shelkovnikov E. Yu., Gulyaev P. V., Zhuykov B. L., Lipanov S. I. Issledovanie vliyaniya makroskopicheskoy formy peretravlivaemoy sheyki zagotovki STM-zonda na ego nanotopologiyu [Research of the influence of the macroscopic shape of the etching neck of the STM probe blank at its nanoscale topology]. Khimicheskaya fizika i mezoskopiya [Chemical Physics & Mesoscopy], 2015, vol. 17, no. 3, pp. 482-487.
3. Tyurikov A. V., Shelkovnikov E. Yu., Gulyaev P. V., Zhuykov B. L., Lipanov S. I. Modelirovanie elektrokhimicheskoy stadii protsessa formirovaniya ostriy STM-zondov [Modelling of the electrochemical stage of process of formation the probes of STM]. Khimicheskaya fizika i mezoskopiya [Chemical Physics & Mesoscopy], 2016, vol. 18, no. 2, pp. 323-330.
4. Shelkovnikov E. Yu., Tyurikov A. V., Kiznertsev S. R., Lipanov S. I. Metodika modelirovaniya protsessa kombinirovannogo travleniya nanozondov [Method simulation of combined etching nanoprobes]. Khimicheskaya fizika i mezoskopiya [Chemical Physics & Mesoscopy], 2012, vol. 14, no. 4, pp. 650-654.
5. Desbrun M., Cascuel M.-P. Smoothed Particles: A new paradigm for animating highly deformable bodies. Computer Animation and Simulation'96: Proceedings of the Eurographics Workshop in Poitiers, France, August 31-September 1, 1996. NewYork: Springer Wien, 1996, pp. 61-76. http://ru.bookzz.org/book/2118951/bf82ee
6. Monaghan J. J. Smoothed Particle Hydrodynamics. Annual Review of Astronomy and Astrophysics, 1992, vol. 30, pp. 543-574. DOI: 10.1146/annurev.aa.30.090192.002551
7. Muller M., Charypar D., Gross M. Particle-Based Fluid Simulation for Interactive Applications. Proceedings of 2003 ACM SIGGRAPH Symposium on Computer Animation, 2003. pp. 154-159. http://dx.doi.org/10.2312/SCA03/154-159
8. Fogel'son R. L., Likhachev E. R. Temperature dependence of viscosity. Technical Physics. The Russian Journal of Applied Physics, 2001, vol. 46, no. 8, pp. 1056-1059.
9. Shelkovnikov E. Yu., Tyurikov A. V., Gulyaev P. V., Zhuykov B. L., Lipanov S. I. Osobennosti primeneniya metoda chastits dlya modelirovaniya protsessa izolyatsii zondiruyushchikh ostriy elektrokhimicheskogo STM [Features of application of the method of particles for modelling the process of isolation probes for electrochemical STM]. Khimicheskaya fizika i mezoskopiya [Chemical Physics & Mesoscopy], 2016, vol. 18, no. 1, pp. 154-159.
Тюриков Александр Валерьевич, кандидат физико-математических наук, старший научный сотрудник ИМ УрО РАН, e-mail: alex. tyurikov@mail. ru
Шелковников Евгений Юрьевич, доктор технических наук, профессор, заведующий лабораторией ИМ УрО РАН, профессор кафедры «Вычислительная техника» ИжГТУ имениМ.Т. Калашникова, e-mail: evshelk@mail.ru
Жуйков Богдан Леонидович, аспирант ИМ УрО РАН, e-mail: fastblood@mail.ru
Кизнерцев Станислав Рафаилович, кандидат технических наук, старший научный сотрудник ИМ УрО РАН, e-mail: iit@udman.ru
