Молекулярно-динамическое моделирование системы заряженных частиц

Основы метода молекулярной динамики. Способы повышения производительности молекулярно-динамического моделирования при помощи различных вычислительных методов и графических ускорителей. Моделирование системы заряженных частиц в виде кластерной наноплазмы.

Рубрика Физика и энергетика
Вид курсовая работа
Язык русский
Дата добавления 09.02.2017
Размер файла 1,8 M

Отправить свою хорошую работу в базу знаний просто. Используйте форму, расположенную ниже

Студенты, аспиранты, молодые ученые, использующие базу знаний в своей учебе и работе, будут вам очень благодарны.

Размещено на http://www.allbest.ru

Москва 2015

ВВЕДЕНИЕ

В современной физике применение вычислительных методов занимает все более важное место и используется при решении широкого класса задач в самых разных областях. В данной работе анализируются данные, полученные при помощи расчётов методом молекулярной динамики (МД), который является одним из наиболее мощных подходов, используемых в моделировании систем многих частиц. Метод МД позволяет на атомных и молекулярных масштабах изучать свойства различных объектов, а также моделировать происходящие в них сложные динамические процессы. Этот метод требует высокой производительности вычислительных систем, поэтому очень важным является изучение и развитие алгоритмов и технических средств, позволяющих производить большие объёмы вычислений.

Основной целью работы является анализ результатов молекулярно-динамического моделирования системы заряженных частиц, представляющей собой плазму, образованную в результате ионизации наноразмерного кластера коротким лазерным импульсом. Данное плазменное образование зачастую называют кластерной наноплазмой. От пространственно однородной плазмы её отличает важная роль размерных эффектов. В частности, в данной работе анализируется явно выраженная зависимость спектра электронных колебаний от размера кластера в интервале от 102 до 105 атомов.

При определенных условиях кластерная плазма может оказаться неидеальной. Изучение неидеальной плазмы занимает особую роль среди фундаментальных задач физики плазмы. Такая плазма отличается от идеальной (газовой) тем, что средняя энергия взаимодействия частиц сопоставима с их кинетической энергией. В то время как для идеальной плазмы, широко распространенной в атмосфере, космосе и в повседневных электронных устройствах, существуют достаточно точные теоретические модели, неидеальная плазма исследована сравнительно мало. Здесь незаменимым инструментом исследования становится компьютерное моделирование методом МД.

Исследование ионизованных наноразмерных кластеров важно для построения теории неидеальной плазмы, а также для практических аспектов, связанных с уникальными свойствами кластерных мишеней, в частности с высоким поглощением энергии лазерного излучения при генерации рентгеновских импульсов или потоков заряженных частиц.

Также в работе оценивается применимость метода Particle-Particle-Particle-Mesh для моделирования указанного типа систем заряженных частиц.

В первой главе излагаются основы метода МД. Во второй главе описаны способы повышения производительности МД моделирования при помощи различных вычислительных методов, а так же использования графических ускорителей. В третьей главе описывается модель моделирования кластерной наноплазмы. В четвёртой главе проводится анализ результатов МД моделирования и доработка существующей теоретической модели ионизированной наноплазмы.

молекулярный динамика заряженный частица кластерный наноплазма

1. МЕТОД МОЛЕКУЛЯРНОЙ ДИНАМИКИ

1.1 Решение уравнений движения

Метод молекулярной динамики - это вычислительный метод, в котором эволюция системы взаимодействующих атомов или частиц во времени определяется численным интегрированием их уравнений движения. В нашем случае используются классические ньютоновские уравнения движения, которые для упрощения построения численной схемы записываются в виде системы дифференциальных уравнений первого порядка:

, (1.1)

, (1.2)

где - координаты частиц, - скорости частиц, - силы, действующие на частицы, - массы частиц.

Уравнения требуют задания начальных условий: координат частиц и их скоростей:

, , i = 1,…,N. (1.3)

Затем уравнения движения интегрируются с заданным шагом по времени. Каждый шаг требуется вычисления сил, действующих на частицы, которые выражаются через заданный потенциал взаимодействия, определяющий свойства системы [1-3]:

. (1.4)

Метод классической молекулярной динамики позволяет рассматривать системы, состоящие из более чем триллионов атомов на временах до нескольких пикосекунд. Однако, максимальный размер системы сильно зависит от сложности потенциала взаимодействия и применения различных методов оптимизации, которые будут рассмотрены в дальнейшем. Системы больших масштабов и на больших временах могут быть рассмотрены в рамках так называемой многомасштабной модели [12].

1.2 Потенциалы взаимодействия

Для каждой задачи требуется индивидуальный выбор потенциала взаимодействия. Например, часто используется потенциал Леннард-Джонса применяется для моделирования взаимодействия инертных газов, других слабовзаимодействующих незаряженных атомов или молекул, а также в качестве модельного потенциала для исследования фундаментальных свойств молекулярно-динамических систем. Потенциал Леннард-Джонса - это короткодействующий потенциал, для которого на сравнительно небольшом расстоянии можно сделать «отсечку», то есть, считать, что на расстоянии больше этого частицы не взаимодействуют.

В общем случае с ростом числа частиц сложность вычислений растёт как O(N2). Однако, для короткодействующих потенциалов существует много методов оптимизации вычислений, например, вычисление списков ближайших соседей. Суть метода в том, что для каждой частицы создаётся список, в который включаются частицы, находящиеся в определённом расстоянии вокруг неё. Считается, что с другими частицами она не взаимодействует. На каждом шаге проводятся расчёты взаимодействий только внутри списка, а сами списки пересчитываются через определённые интервалы по времени. Применение такого метода (часто в совокупности с методом связанных списков) позволяет добиться роста сложности как О(N).

В случае короткодействующего потенциала достаточно просто использовать параллельные вычисления, так как наличие отсечки позволяет разделить область вычислений на несколько более мелких областей, и расчёт на каждом из участков можно проводить отдельным процессором. Грамотно подобранное разбиение на области и составление списков ближайших соседей позволяет сократить время передачи информации между процессорами и оптимизировать вычисления.

Тем не менее, большое количество задач требует применения дальнодействующих потенциалов. Например, расчёт электростатических или гравитационных взаимодействий. В этом случае отсечку нужно устанавливать на достаточно большом расстоянии, и разбиение на области становится неэффективным в силу того, что тратится много времени на коммуникацию между процессорами. Сложность вычисления сильно возрастает с ростом числа частиц, и для таких задач необходимы иные методы оптимизации [1], которые, однако, не являются столь же эффективными, как в случае короткодействующих потенциалов, поэтому размеры моделируемых системы обычно ограничиваются тысячами или миллионами частиц.

2. ОПТИМИЗАЦИЯ МОДЕЛИРОВАНИЯ СИСТЕМ ЗАРЯЖЕННЫХ ЧАСТИЦ

2.1 Приближенные методы моделирования систем заряженных частиц

Для ускорения моделирования заряженных частиц, взаимодействующих посредством дальнодействующего куловноского потенциала, в больших системах существуют различные вычислительные методы. Например, метод сумм Эвальда [1, 13] очень популярен в современной молекулярной динамике, однако он применим только к периодическим системам. Этот метод является частным случаем формулы суммирования Пуассона, в которой суммирование энергий взаимодействия в реальном пространстве заменяется на эквивалентное суммирование в пространстве Фурье.

Применение сеток или решёток для ускорения расчётов для N тел долгое время являлось стандартной практикой за пределами стандартной молекулярной динамики, например, в астрофизике и физике плазмы, где необходимо отслеживать макроскопические тенденции на длинных отрезках времени. Было обнаружено, что в случаях, когда необходимо узнать общие параметры системы и не требуется высокая точность, данные методы работают очень хорошо. В связи с этим, методы вычислений на решётках во многих ситуациях заменяют прямые вычисления.

Метод PPPM (Particle-Particle, Particle-Mesh) содержит в себе плюсы точных молекулярно-динамических расчётов и вычислений на сетках [1]. Сила взаимодействия разделяется на два слагаемых: Fij = FijPP + FijPM , где РР (Particle-Particle) - непосредственное межчастичное взаимодействие, вычисляемая до определённого радиуса отсечки, а РМ (Particle-Mesh) считается достаточно гладкой функцией для вычисления на сетке. Вычисления на сетке производятся путем численного решения уравнения Пуассона в Фурье-пространстве.

Для проверки применимости метода PPPM к задаче моделирования кластерной наноплазмы был использован пакет программ LAMMPS (англ. Large-scale Atomic/Molecular Massively Parallel Simulator) [14]. Это свободно распространяемый пакет программ для классической молекулярной динамики, разработанный в Сандийской национальной лаборатории США. LAMMPS может применяться как для небольших расчётов, так и для изучения взаимодействия миллионов частиц. В этом пакете реализовано большое количество потенциалов и методов, применяемых для оптимизации вычислений, например, пересчёт списков ближайших соседей. Так же LAMMPS поддерживает параллельные вычисления. В частности, в пакете LAMMPS реализован метод PPPM. В данной работе был проведён тестовый расчет на LAMMPS с целью выяснения, подходят ли использованные в нем приближённые вычисления для задачи моделирования наноплазмы.

Также следует отметить метод TreeMD [1], использованный, в частности, и для задач моделирования кластерной наноплазмы15]. Тем не менее, результаты, полученные этим методом имеют явные отличия, от результатов, полученных точным МД моделированием [8,11]

2.2 Расчет взаимодействия частиц методом PPPM

При помощи LAMMPS было смоделировано движение двух отрицательно заряженных частиц с кулоновским потенциалом взаимодействия. За образец взят расчёт без применения PPPM. По траектории движения частиц построена зависимость потенциальной энергии от расстояния (рис. 2.1).

Рис. 2.1 Кулоновский потенциал для двух отрицательно заряженных частиц

Для вычислений был использован потенциал coul/cut с отсечкой на расстоянии 75 атомных единиц. Во всех расчётах сделано 2000 шагов по 0.007фс.

При использовании PPPM на результат влияют несколько параметров. Одним из них является размер ячейки, в котором проводится расчёт. При уменьшении ячейки возрастает погрешность относительно исходного потенциала, так как меньшее количество частиц попадают в зону непосредственного расчёта взаимодействия. Для проверки взяты ячейки размером от 100 до 1000, измерялась средняя ошибка. Результаты представлены на рис. 2.2.

Рис. 2.2 Зависимость средней погрешности от размера ячейки

Ещё одним параметром, влияющим на результат, является точность Фурье-преобразования, указываемая в параметре kspace. Зависимость погрешности определения потенциала взаимодействия от этого параметра приведена на рис. 2.3. На рисунке видно, что в точке 0.0001 погрешность увеличивается. Это связано с автоматическим выбором размера сетки. При 0.0001 программа автоматически выбирает размер сетки равным 4 а.е., а при других параметрах больше. Это связано с тем, что параметр размера сетки должен быть кратным степеням чисел 2, 3, 5. То есть, если на основе остальных параметров программой автоматически был выбран неподходящий размер, резко возрастает погрешность.

Рис. 2.3 Зависимость средней погрешности от точности kspace.

Расчёты на PPPM занимали достаточно много времени: взаимодействие двух частиц с точностью 0. 001 занимало вплоть до 34.3807с. Это связано с тем, что изучаемая система не была периодической, частицы могут вылетать достаточно далеко и не должны возвращаться с другой стороны расчётной области, поэтому приходится брать размер области очень большим. При этом укрупнение сетки приводило к большой погрешности.

Вывод: применение методов приближённого вычисления в случае задачи, требующей крупного размера вычислительной области нецелесообразно. Следовательно, для расчётов, необходимых для данной работы, следует применять точные методы вычисления.

2.3 Применение графических ускорителей

Метод МД, в особенности применительно к системам с дальнодействующими потенциалами, требует высокой производительности вычислительных систем, и это приводит к тому, что решение задач МД моделирования требует большой вычислительной мощности компьютеров. В связи с этим очень важной становится эффективность распараллеливания. На скорость работы МД программы вообще и на эффективность распараллеливания в частности влияет число частиц в системе, а так же характер их взаимодействий [4].

В последнее время все чаще используются гибридные системы, в которые помимо обычных процессоров входят графические ускорители (GPU). Изначально созданные как часть видеокарты, графические ускорители продемонстрировали эффективность в научных расчётах. Высокая вычислительная мощность GPU достигается за счёт особой многоядерной архитектуры.

Из разницы в архитектуре следует разница в принципах работы. Если архитектура CPU предполагает последовательную обработку информации, то GPU рассчитаны на массивно параллельные вычисления. Особенностью графических ускорителей является большое количество простых вычислительных ядер. В отличие от ядер центрального процессора, графические ядра способны совершать только простые математические операции и не способны к поддержке операционной системы, таким образом, GPU фактически является математическим сопроцессором к центральному процессору.

Гибридные кластеры на основе GPU зачастую показывают быстродействие, во много раз превышающее быстродействие систем на универсальных центральных процессорах (см. пример расчетов для Леннард-Джоносоновской жидкости на рис. 3.1).

Рис 3.1. Отношение производительности GPU Nvidia X2050 к одному ядру CPU Intel Xeon X5670 при моделировании леннард-джонсовской жидкости с параметрами: ? = 0.19, T=1.0, rcut= 3, rskin= 0.8. Число 9, тра: 106.

Несмотря на преимущества вычислений на GPU, нужно помнить, что создание программ, эффективно использующих их вычислительные возможности - это сложный процесс, требующий учета особенностей архитектуры конкретного устройства применительно к конкретной задаче, без которого производительность может оказаться недостаточно высокой. На рисунке 3.1 можно видеть, что наибольшая эффективность в сравнении с вычислениями на СPU при расчётах МД достигается на системах с большим количеством частиц, когда загружено большое количество ядер GPU, а так же в случае применения оптимизированного алгоритма [5].

Одним из условий эффективного использования GPU является хорошая распараллеливаемость программы. Метод молекулярной динамики подходит под это условие. В настоящее время существует много пакетов программ, позволяющих проводить МД моделирование на графических ускорителях, однако набор потенциалов взаимодействия частиц по-прежнему очень невелик и ограничивается, как правило, простейшими типа потенциала Леннард-Джонса. Для моделирования большинства реальных веществ требуется большее количество потенциалов, чем доступно в существующих пакетах, поэтому адаптация необходимого потенциала под вычисления на графических ускорителях становится отдельной задачей на начальных этапах исследований.

Максимальное ускорение для алгоритма моделирования кластерной наноплазмы, достигнутое на GPU Nvidia Tesla M2050 по сравнению с одним ядром центрального процессора (CPU) Intel Xeon E5520, составило 100 раз (рис. 3.2). Как видно из рисунка, эффективность применения GPU существенным образом снижается при переходе к малому числу частиц N < 104. Это объясняется малой загрузкой GPU при использовании одного потока команд (thread) на частицу.

Таким образом, применение GPU оказывается эффективным в случае моделирования систем заряженных частиц. Это избавляет от необходимости применять приближенные методы типа PPPM, которые, как показано в предыдущем параграфе, плохо применимы для расчета кластерной наноплазмы.

Рис 3.2 Отношение быстродействия GPU Nvidia Tesla M2050 к одному ядру CPU Intel Xeon E5520 при моделировании неидеальной плазмы в зависимости от количества частиц в системе. Крестики - алгоритм с одним потоком на частицу, треугольники - алгоритм с несколькими потоками на частицу.

Малое суммарное число потоков не позволяет задействовать внутренний механизм GPU, скрывающий задержки доступа к памяти за счет быстрого аппаратного переключения между ждущими и готовыми к выполнению потоками. Для оптимизации полученного алгоритма применяет метод расчета, увеличивающий число потоков на частицу для малых N (верхний ряд точек на рис. 3.2).

3. МОДЕЛИРОВАНИЕ НЕИДЕАЛЬНОЙ ПЛАЗМЫ, ОБРАЗОВАННОЙ В РЕЗУЛЬТАТЕ ИОНИЗАЦИИ НАНОРАЗМЕРНЫХ КЛАСТЕРОВ

3.1 Постановка задачи

Исследование ионизованных наноразмерных кластеров важно для построения теории неидеальной плазмы, а так же для практического применения, связанного с уникальными свойствами кластерных мишеней, в частности с высоким поглощением энергии лазерного излучения при генерации рентгеновских импульсов или потоков заряженных частиц.

При помощи фемтосекундных лазерных установок можно изучать колебания электронов в наноплазме, образующейся в результате быстрой ионизации кластеров короткими лазерными импульсами [6].

Электронные колебания в наноплазме исследовались экспериментально. Кластеры облучались двумя последовательными лазерными импульсами с контролируемой задержкой. В результате первого импульса происходит ионизация атомов кластера, а также его зарядка за счет вылетевших электронов. Иллюстрация разлёта ионов показана на рисунке 4.1. Вследствие суммарного отрицательного заряда кластера происходит разлёт мишени (кулоновский взрыв). Следующий импульс позволяет определить коэффициент резонансного поглощения лазерного излучения в образовавшейся наноплазме на различных этапах кулоновского взрыва. При этом в экспериментах было показано, что резонансная частота колебаний электронов может отличаться от известной частоты Ми.

Рис. 4.1 Ионизация и разлет кластера. Лазерный импульс: I = 5·1012 Вт/см2, t = 50 фс. Чёрным обозначены траектории ионов, красным - электронов.

Для построения теоретической модели указанных выше процессов была написана программа молекулярно-динамического моделирования неидеальной электрон-ионной плазмы [7,11]. В ней были смоделированы колебания электронов в неидеальной наноплазме при неподвижных ионах.

3.2 Численная модель моделирования кластерной наноплазмы

Классический кулоновский потенциал электрон-ионного взаимодействия на малых расстояниях быстро убывает, из-за этого на короткой дистанции при расчётах получаются очень большие по модулю отрицательные значения энергии. Это противоречит квантово-механическим представлениям о наличии минималного уровня энергии (основного состояния) в системе электрона и иона. Для приближенного учета этого явления куловновский потенциал взаимодействия принято заменять псевдопотенциалом (сравнение на рис. 4.2), который при нулевом расстоянии между атомами имеет определённое значение. Использованный в [9,10] потенциал имеет вид:

, (3.1)

где , . Значение параметра выбрано так, что при стремлении к нулю значение энергии получается равным энергии ионизации изолированного атома натрия .

Для интегрирования уравнений движения применялась схема с перешагиванием (Leap-frog) с шагом по времени три аттосекунды. В начале расчёта ионы располагались в узлах икосаэдрической кристаллической решетки, обрезанной по сфере радиуса R.

Разлёт ионов не моделировался, считалось, что ионы неподвижны в определённых конфигурациях. Это оправдано, так как время расширения ионной сферы намного больше, чем время установления равновесия в электронной подсистеме. Для изучения разлёта электронов рассматривались только ранние стадии разлета кластера.

Рис. 4.2. Потенциалы электрон-ионного взаимодействия на коротких расстояниях: 1 - кулоновский, 2 -

Рассматривалась изначально ионизованная кластерная наноплазма при определенной температуре, задаваемой на начальном этапе термостатом. После достижения стационарных значений температуры электронов и заряда кластера термостат выключался, и основной расчет проводился для консервативной системы с сохранением полной энергии. Полученные в результате моделирования результаты усреднялись по ансамблю начальных конфигураций ионов.

3.3 Спектры электронных колебаний в ионизованном нанокластере

Анализ колебаний электронов в электронной подсистеме проводился путем вычисления спектра автокорреляционной функции тока [8]:

, (3.2)

, (3.3)

где - полный электронный ток в кластере, - нормировочный объем, - диэлектрическая проницаемость.

Рис. 4.3. Спектры автокорреляционной функции тока для различного числа ионов в кластере Ni (размера кластера R). Пунктирные линии показывают теоретические значения частоты Ми и плазменной частоты .

Теоретическое значение частоты Ми имеет вид . Это значение получено на основе модели, в которой электроны представляют собой сферическое облако, которое колеблется около примерно такого же по размеру облака ионов. Однако в результате удара лазером некоторое количество электронов вылетело, следовательно, кластер имеет положительный заряд, поэтому радиус электронного облака меньше. Это объясняет сдвиг резонанса Ми для маленьких кластеров [5].

Вставка на рисунке 4.3 показывает, что амплитуда колебаний электронов, соответствующих второму максимуму, в разных частях кластера имеет разный знак. Далее будем называть этот пик плазменным резонансом. Теоретическое значение его частоты имеет вид . С увеличением размера кластера амплитуда этого резонанса резко возрастает, и при он становится доминирующим, в то время как резонанс Ми затухает. Как видно из рисунка, частоты первого и второго резонансов стремятся к частоте Ми при увеличении размера кластера, а для частоты третьего резонанса предельным значением является плазменная частота.

Положения основных резонансов в зависимости от размера кластера показаны на рис. 4.4. Определение плазменной частоты становится неточным для кластеров с .

Рис. 4.4. Зависимость частот основных мод колебаний электронов от размера кластера (количества ионов Ni). Квадраты и треугольники соответствуют двум максимумам на кривой вблизи частоты Ми; круги - положение максимума вблизи плазменной частоты

Полученные отклонения собственных частот электронных колебаний от теоретических значений потребовали дополнительного анализа, который и является основным предметом исследования в данной работе. В следующей главе подробно рассматриваются причины сдвига частоты Ми.

4. АНАЛИЗ КОЛЕБАНИЙ ЭЛЕКТРОНОВ В КЛАСТЕРНОЙ НАНОПЛАЗМЕ

4.1 Теоретическая оценка частоты электронных колебаний

Результаты молекулярно-динамического моделирования неидеальной плазмы потребовали дополнительного теоретического обоснования. Был проведён анализ и проверка вывода теоретических формул, а так же анализ результатов молекулярно-динамического моделирования и дополнение теоретической базы по проведённому эксперименту.

Выведем формулу коллективных колебаний электронов в кластере. Пусть облако электронов имеет фиксированную форму. Допустим, что плотность электронов n(r) и потенциальная энергия электрона в кластере V(r)=-(r) сферически симметричны и их центры находятся в равновесии. Предположим сдвиг z относительно минимума потенциала. Необходимая частота в гармоническом приближении:

,

где U - энергия всей системы, а m=Neme - это масса Ne электронов в облаке. В то же время ДU может быть выражено через V(r) и n(r).

Если потенциал непрерывный, можно продифференцировать его перед интегрированием

Полагая z=0, с=r, можно получить:

(4.1)

Полученная в результате формула (4.1) - это формула электронных колебаний для произвольного потенциала. В данной задаче используется потенциал (3.3).

В работе [8] сделано предположение о том, что распределения электронов и ионов представляют собой однородные сферы разного радиуса (за счет вылета части электронов). На основе этого предположения выведена формула:

(4.2)

При помощи программы Wolfram была проверена корректность вывода этой формулы. График частоты, вычисленной по этой формуле, изображён на рис. 3.4 сплошной линией. На нём видно, что на эту кривую не попадают точки, обозначающие результаты молекулярной динамики, следовательно, требуется ввести в неё дополнительные параметры, улучшающие точность.

Для того чтобы аналитические расчёты давали результат, лучше согласованный с результатами молекулярной динамики, было предложено учитывать реальный профиль электронной плотности. Это означает, что электронное облако считается шаром с плотностью, зависящей от радиуса, в то время как ионы остаются замороженными, и их распределение по-прежнему считается однородной сферой. Необходимый для этого профиль электронной плотности брался из результатов МД моделирования путем усреднения средней плотности электронов в различных точках на протяжении достаточно длинной МД траектории, пример его сравнения с постоянной плотностью на рис. 4.1.

Рис. 4.1. Пример сравнения реального профиля электронной плотности с постоянной электронной плотностью. Размер кластера - 3150 частиц.

Расчёты проводились в программе Wolfram Mathematica. Для подстановки в формулу колебаний (5.1) был взят коллективный потенциал кластера: , где б = e, i - тип частиц. После подстановки erf-like потенциала 3.1 получаем (здесь и далее код из программы Mathematica):

(4.3)

Для того чтобы программа могла корректно аналитически проинтегрировать потенциал, он был разбит на две составляющие: первая для случая r > r1, вторая для r1 > r:

Затем аналитически было получено выражение fV(r) = и взят его интеграл.

В результате этого получилась следующая функция:

(4.4)

Следовательно, выражение для частоты приняло следующий вид:

, (4.5)

где n(r) - электронная плотность.

Для определения новых значений плазменной частоты электронная плотность считывалась как таблица из результатов расчётов молекулярной динамики для различных размеров кластеров. По таблице строилась интерполяционная функция (использовалась сплайн-интерполяция), а затем численно брался интеграл (4.5). Рисунок 4.2 показывает, что новые результаты лучше приближают данные, полученные из молекулярно-динамических расчётов.

Рис. 4.2 Зависимость резонансных частот от числа электронов в кластере. Ромбы - результаты численного интегрирования (4.5) с учётом реальной электронной плотности, полученной из МД моделирования.

Данный результат указывает на то, что смещение частоты Ми в зависимости от размера кластера действительно может быть объяснено только лишь деформацией профиля электронной плотности.

4.2 Параметризация профиля электронной плотности

Предположение о важности роли реального профиля электронной плотности в вычислении плазменного резонанса подтвердилось, следовательно, потребовалось дальнейшее его изучение. В результате молекулярно-динамических расчётов были получены массивы данных, соответствующие профилям плотности для кластеров различных размеров при температуре 1эВ. Для всех профилей были построены соответствующие графики.

При изучении графиков можно увидеть закономерность, с которой меняется форма профиля при увеличении числа ионов в кластере. Было выдвинуто предположение, что возможно вывести аналитическую формулу, описывающую зависимость электронной плотности в кластере от радиуса и температуры. Вывод осуществлялся при помощи программы Wolfram Mathematica, которая позволяет наглядно отслеживать зависимость искомой функции от параметров и минимизировать отклонения от референтных данных МД моделирования.

Сначала была создана тестовая функция, зависящая от нескольких параметров, определяющих положение последнего пика, его амплитуду, плавности возрастания и убывания, а так же положение точки перехода к затухающей асимптоте на больших расстояниях. Тестовая функция, приближающая профиль плотности, имела вид:

, (4.6)

где fFermi(x,x0,s) - это функции Ферми вида . Здесь и далее расстояния измеряются в нм, концентрация электронов в нм-3, частота в фс-1.

Рис. 4.3 Сравнение выбранной функции с исходными данными. Подбор параметров при помощи опции Manipulate в Wolfram Mathematica. Зеленая линяя - функция, синяя - данные из массива результатов моделирования.

N = 55.

Требовалось сократить количество параметров для того, чтобы можно было определять профиль плотности, зная только число частиц. Для этого были подобраны оптимальные наборы параметров такие, что отклонение функции от исходного набора данных было минимальным. Подбор проводился для результатов расчётов в кластерах размером 55 ионов (рис. 4.4), 147 (рис. 4.5), 309 (рис. 4.6), 1000 (рис. 4.7), 3150 (рис. 4.8), 7100 (рис. 4.9).

Рис. 4.4 Сравнение выбранной функции с исходными данными.

N = 147.

Рис. 4.5 Сравнение выбранной функции с исходными данными. N = 309.

Рис. 4.6 Сравнение выбранной функции с исходными данными.

N = 1000.

Рис. 4.7 Сравнение выбранной функции с исходными данными.

N = 3150.

Рис. 4.8 Сравнение выбранной функции с исходными данными.

N = 7100.

Таб. 4.1 Параметры функции.

N

N1/3

c1,

вершина

c2,

высота пика

st, степень сглаживания

параметр наклона функции Ферми

частота

per, точка

перехода

55

3,803

0,359

0,900

0,990

0,08

8

0,405

147

5,278

0,545

0,412

2,860

0,08

8

0,670

309

6,761

0,800

0,474

3,420

0,08

8

0,900

1000

10,000

1,210

0,260

5,010

0,08

8

1,375

3150

14,659

1,718

0,130

6,090

0,08

8

2,015

7100

19,220

2,220

0,125

7,560

0,08

8

2,600

14100

24,159

3,150

0,120

10,430

0,08

8

3,490

После подбора оптимальных параметров была установлена их зависимость друг от друга или от размера кластера и температуры. Параметр c1 отвечал за сдвиг вершины последнего пика и являлся, по сути, главным для повторения формы реального профиля, зависимость получилась квадратичной (рис. 4.9). Точка переключения между функциями Ферми отвечала за переход к асимптоте на больших расстояниях, рис. 4.10. Гладкость переключения и частота возмущения для всех размеров кластеров совпали со значениями 0.08 и 8 соответственно. Параметр амплитуды пика отвечал за интенсивность уплотнения на краях кластера (рис. 4.11), а степень сглаживания определяла плавность перехода от относительно равномерного распределения плотности в центре кластера к неравномерной плотности по краям (рис. 4.12).

Рис. 4.9 Зависимость положения вершины от кубического корня из числа частиц.

Рис. 4.10 Зависимость положения т. переключения между функциями Ферми от N1/3

.

Рис. 4.11 Зависимость амплитуды пика от N1/3.

Рис. 4.12 Зависимость степени сглаживания от N1/3.

Затем была построена итоговая функция

,

где параметры выражены через следующим образом:

Эта функция была подставлена в уравнение (4.5), в результате чего получилось окончательное выражение для численного интегрирования:

Рис. 4.13 Зависимость частоты резонанса от числа ионов.

В результате интегрирования получилась кривая значений частоты, нанесённая на график вместе с предыдущими результатами на рис. 4.14. Эта кривая гораздо лучше описывает результаты МД моделирования.

Рис. 4.14 Зависимость резонансных частот от числа электронов в кластере. Оранжевая линия - результат.

ВЫВОДЫ

Проверена гипотеза о том, что смещение частоты резонанса Ми в сторону меньших значений связано с неоднородностью профиля электронной плотности и суммарным зарядом кластера. Предложена интерполяционная формула для профиля электронной плотности в кластере произвольного размера. Результаты применения этой формулы показывают хорошее совпадение с данными МД моделирования по частоте колебаний Ми.

Проведено сравнение двух методик ускорения МД моделирования систем заряженных частиц: приближение PPPM и полноатомная модель с расчетом сил и решением уравнений движения на графических ускорителях. Показаны явные преимущества использования полноатомной модели.

СПИСОК ЦИТИРОВАННОЙ ЛИТЕРАТУРЫ

[1] Gibbon P., Sutmann G. Long-Range Interactions in Many-Particle Simulation. In: Quantum Simulations of Complex Many-Body Systems: From Theory to Algorithms eds. J. Grotendorst, et al, Julich: NIC. 2002. V. 10. P. 467-506.

[2] Allen M.P., Tildesley D.J. Computer Simulation of Liquids. Oxford : Clarendon Press, 1989.

[3] Frenkel D., Smit B. Understanding Molecular Simulation: From Algorithms to Applications. San Diego: Academic Press, 2002.

[4] I.V. Morozov, A.M. Kazennov, R.G. Bystryi, G.E. Norman, V.V. Pisarev, V.V. Stegailov. Comput. Phys. Commun. 2011. V. 182. P. 1974.

[5] T. Raitza, H. Reinholz, G. Rцpke, I. Morozov, E. Suraud. Laser excited expanding small clusters: Single time distribution functions // Contributions to Plasma Physics. 2009. V. 49. P. 496-506.

[6] T. Doeppner, T. Fennel, P. Radcliffe, J. Tiggesbaumker, and K.-H. Meiwes-Broer, Phys. Rev. A 73, 031202 (2006).

[7] H. Reinholz, T. Raitza, G. Roepke, and I. V. Morozov, Int. J. Mod. Phys. B 22, 4627 (2008).

[8] T. Raitza, G. Roepke, G. Reinholz, and I. V. Morozov, Phys. Rev. E 84, 036406 (2011).

[9] G. Zwicknagel and T. Pschiwul, Contrib. Plasma Phys. 43, 393 (2003).

[10] A. Filinov, M. Bonitz, and W. Ebeling, J. Phys. A 36, 5957 (2003).

[11] R.G. Bystryi, I.V. Morozov. Electronic oscillations in ionized sodium nanoclusters // J. Phys. B: At. Mol. Opt. Phys. 2015. V. 48. P. 015401.

[12] Тимофеев А.В., Морозов И.В., Стегайлов В.В., Норман Г.Э., Куксин А.Ю., Ланкин А.В., Орехов Н.Д., Писарев В.В., Смирнов Г.С., Стариков С.В. Зачем и какие нужны суперкомпьютеры эксафлопсного класса? Предсказательное моделирование свойств и многомасштабных процессов в материаловедении // Программные системы: теория и приложения. 2014. Т. 5. № 1. С. 191-244.»

[13] Ewald, P (1921). "Die Berechnung optischer und elektrostatischer Gitterpotentiale". Ann. Phys. 369 (3): 253-287.

[14] S. Plimpton, Fast Parallel Algorithms for Short-Range Molecular Dynamics, J Comp Phys, 117, 1-19 (1995).

[15] Winkel M., Gibbon P. Spatially Resolved Electronic Correlations in Nanoclusters // Contrib. Plasma Phys. 2013. V. 53. P. 254-262.

Размещено на Allbest.ru

...

Подобные документы

  • Ускорители заряженных частиц — устройства для получения заряженных частиц больших энергий, один из основных инструментов современной физики. Проектирование и испытание предшественников адронного коллайдера, поиск возможности увеличения мощности систем.

    реферат [685,8 K], добавлен 01.12.2010

  • Ускорители заряженных частиц как устройства, в которых под действием электрических и магнитных полей создаются и управляются пучки высокоэнергетичных заряженных частиц. Общая характеристика высоковольтного генератора Ван-де-Граафа, знакомство с функциями.

    презентация [4,2 M], добавлен 14.03.2016

  • Взаимодействие заряженных частиц и со средой. Детектирование. Определение граничной энергии бета-спектра методом поглощения. Взаимодействие заряженных частиц со средой. Пробег заряженных частиц в веществе. Ядерное взаимодействие. Тормозное излучение.

    курсовая работа [1,1 M], добавлен 06.02.2008

  • Понятие и принцип работы ускорителей, их внутреннее устройство и основные элементы. Ускорение пучков частиц с высокой энергией в электрическом поле как способ их получения. Типы ускорителей и их функциональные особенности. Генератор Ван де Граафа.

    контрольная работа [276,8 K], добавлен 18.09.2015

  • Динамика частиц, захваченных геомагнитным полем, ее роль в механизме динамики космического изучения в околоземном пространстве. Геометрия радиационных поясов Земли. Ускорение частиц космического излучения. Происхождение галактических космических лучей.

    дипломная работа [1,2 M], добавлен 24.06.2015

  • Сущность молекулярно-динамического моделирования. Обзор методов моделирования. Анализ дисперсионного взаимодействия между твердой стенкой и жидкостью. Использование результатов исследования для анализа адсорбции, микроскопических свойств течения жидкости.

    контрольная работа [276,7 K], добавлен 20.12.2015

  • Изучение процессов рассеяния заряженных и незаряженных частиц как один из основных экспериментальных методов исследования строения атомов, атомных ядер и элементарных частиц. Борновское приближение и формула Резерфорда. Фазовая теория рассеяния.

    курсовая работа [555,8 K], добавлен 03.05.2011

  • Коэффициенты диффузии, ступенчатые поверхности. Алгоритм Метраполиса, метод Монте-Карло, парциальное и среднее покрытие, термодинамический фактор. Диффузия системы взаимодействующих частиц. Зависимость среднего покрытия от химического потенциала.

    дипломная работа [1,2 M], добавлен 10.12.2013

  • Создание большого адронного коллайдера, ускорителя заряженных частиц на встречных пучках. Предназначение его для разгона протонов и ионов, изучение продуктов их соударений. Изучение космических лучей, моделируемых с помощью несталкивающихся частиц.

    презентация [1,1 M], добавлен 16.04.2015

  • Анализ естественных и искусственных радиоактивных веществ. Методы анализа, основанные на взаимодействии излучения с веществами. Радиоиндикаторные методы анализа. Метод анализа, основанный на упругом рассеянии заряженных частиц, на поглощении P-частиц.

    реферат [23,4 K], добавлен 10.03.2011

  • Энергетическое разрешение полупроводникового детектора. Механизмы взаимодействия альфа-частиц с веществом. Моделирование прохождения элементарных частиц через вещество с использованием методов Монте–Карло. Потери энергии на фотоядерные взаимодействия.

    курсовая работа [502,5 K], добавлен 07.12.2015

  • Создание физической модели деформации материала. Система кластеров структурированных частиц. Описание механики процесса пластической деформации металла при обработке давлением и разрушения материала при гидрорезке на основе кавитации, резонансных явлений.

    статья [794,6 K], добавлен 07.02.2014

  • Понятие электрического тока как упорядоченного движения заряженных частиц. Виды электрических батарей и способы преобразования энергии. Устройство гальванического элемента, особенности работы аккумуляторов. Классификация источников тока и их применение.

    презентация [2,2 M], добавлен 18.01.2012

  • Понятие и основные положения молекулярно-кинетической теории. Диффузия как самопроизвольное перемешивание соприкасающихся веществ. Броуновское движение – беспорядочное движение частиц. Молекула - система из небольшого числа связанных друг с другом атомов.

    презентация [123,0 K], добавлен 06.06.2012

  • Математическая модель и решение задачи очистки технических жидкостей от твердых частиц в роторной круговой центрифуге. Система дифференциальных уравнений, описывающих моделирование процесса движения твердой частицы. Физические характеристики жидкости.

    презентация [139,6 K], добавлен 18.10.2015

  • Физические свойства висмута и его полиморфных модификаций. Исследование влияния мощных пучков заряженных частиц на микроструктуры и свойства мишеней. Преимущества применения методов рентгеноструктурного фазового анализа для расчета дифракционных картин.

    курсовая работа [5,2 M], добавлен 13.08.2013

  • История открытия рентгеновского излучения. Источники рентгеновских лучей, их основные свойства и способы регистрации. Рентгеновская трубка, ускорители заряженных частиц. Естественная и искусственная радиоактивность. Применение рентгеновского излучения.

    презентация [427,3 K], добавлен 28.11.2013

  • Определения молекулярной физики и термодинамики. Понятие давления, основное уравнение молекулярно-кинетической теории. Температура и средняя кинетическая энергия теплового движения молекул. Уравнение состояния идеального газа (Менделеева - Клапейрона).

    презентация [972,4 K], добавлен 06.12.2013

  • Основные понятия и определения молекулярной физики и термодинамики. Основное уравнение молекулярно-кинетической теории. Температура и средняя кинетическая энергия теплового движения молекул. Состояние идеального газа (уравнение Менделеева-Клапейрона).

    презентация [1,1 M], добавлен 13.02.2016

  • Движение несвободной частицы. Силы реакции и динамика частиц. Движение центра масс, закон сохранения импульса системы. Закон сохранения кинетического момента системы. Закон сохранения и превращения механической энергии системы частиц. Теорема Кёнига.

    доклад [32,7 K], добавлен 30.04.2009

Работы в архивах красиво оформлены согласно требованиям ВУЗов и содержат рисунки, диаграммы, формулы и т.д.
PPT, PPTX и PDF-файлы представлены только в архивах.
Рекомендуем скачать работу.