Фазовые переходы МД — различия между версиями

Материал из Department of Theoretical and Applied Mechanics
Перейти к: навигация, поиск
(Задача I)
 
(не показано 20 промежуточных версий 7 участников)
Строка 9: Строка 9:
 
Написать графический интерфейс, позволяющий наблюдать движение частиц. Предусмотреть возможность отключаемого отображения: температуры (цветом), скорости (светом и отрезком), связей между частицами (отрезком). Реализовать возможность выбора частицы мышкой и вывода подробной информации (номер, скорость, сила).
 
Написать графический интерфейс, позволяющий наблюдать движение частиц. Предусмотреть возможность отключаемого отображения: температуры (цветом), скорости (светом и отрезком), связей между частицами (отрезком). Реализовать возможность выбора частицы мышкой и вывода подробной информации (номер, скорость, сила).
 
Список Группы:
 
Список Группы:
* [http://tm.spbstu.ru/Нарядчиков_Александр Нарядчиков Александр]
+
* [[Нарядчиков Александр]]
* Лебедев Станислав
+
* [[Лебедев Станислав]]
* Демченко Артем
+
* [[Демченко Артем]]
* Киселев Лев
+
* [[Киселёв Лев]]
  
 
==Задача II==
 
==Задача II==
Строка 19: Строка 19:
  
 
Список Группы:
 
Список Группы:
* Абрамов Игорь
+
* [[Абрамов Игорь]]
* Ляжков Сергей
+
* [[Ляжков Сергей]]
* Сенников Иван
+
* [[Сенников Иван]]
* Степаняц Степан
+
* [[Степаняц Степан]]
* Лосева Татьяна
+
* [[Лосева Татьяна]]
  
 
==Задача III==
 
==Задача III==
Строка 30: Строка 30:
  
 
Список Группы:
 
Список Группы:
* Давыдова Алена
+
* [[Давыдова Алена]]
* Бальцер Анастасия
+
* [[Бальцер Анастасия]]
* Васильева Анастасия
+
* [[Васильева Анастасия]]
* Иванова Яна
+
* [[Иванова Яна]]
* Лобанов Илья
+
* [[Лобанов Илья]]
  
 
==Задача IV==
 
==Задача IV==
Строка 41: Строка 41:
  
 
Список Группы:
 
Список Группы:
*[http://tm.spbstu.ru/%D0%A0%D1%83%D0%B1%D0%B8%D0%BD%D0%BE%D0%B2%D0%B0_%D0%A0%D0%B0%D0%B8%D1%81%D0%B0 Рубинова Раиса ]
+
*[[Рубинова Раиса ]]
*Андреева Полина
+
*[[Андреева Полина]]
*Белоусова Екатерина
+
*[[Белоусова Екатерина]]
*Тимошенко Валентина
+
*[[Тимошенко Валентина]]
*Уманский Александр
+
*[[Уманский Александр]]
  
 
==Решение задачи==
 
==Решение задачи==
 +
Открывать лучше в Mozile FireFox, либо настраивать аппаратное ускорение самому (если программа не открывается)
 
{{#widget:Iframe |url=http://tm.spbstu.ru/htmlets/Lebedev/Phase_Transition_MD_2017/Main.html |width=1300 |height=750 |border=0 }}
 
{{#widget:Iframe |url=http://tm.spbstu.ru/htmlets/Lebedev/Phase_Transition_MD_2017/Main.html |width=1300 |height=750 |border=0 }}
 +
 +
== Потенциал Бреннера второго поколения ==
 +
 +
Потенциал Бреннера второго поколения позволяет представить энергию связи
 +
в виде
 +
 +
<math>
 +
    E_b = \sum_i \sum_{j (> i)} \left[ V^R (r_{ij}) - b_{ij} V^A (r_{ij}) \right].
 +
</math><br />
 +
 +
Силу, действующую на частицу с номером I можно рассчитать как минус градиент энергии (производная по радиус-вектору частицы i)
 +
 +
<math>
 +
\vec{F}_i = - \sum \limits_{j \neq i} \left[ \frac{\partial V^R(r_{ij})}{\partial \vec{r}_i}-\frac{\partial b_{ij}}{\partial \vec{r}_i}V^A(r_{ij})-b_{ij}\frac{\partial V^A(r_{ij})}{\partial \vec{r}_i} \right],
 +
</math>
 +
 +
Между атомами углерода функции отталкивания и притяжения имеют вид:
 +
 +
<math>
 +
\frac{\partial V^R(r_{ij})}{\partial \vec{r}_i} = \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} \left( 1+\frac{Q}{r} \right) A e^{-ar_{ij}}+f^c(r_{ij}) Q \frac{\vec{r}_{ij}}{r_{ij}^3}Ae^{-ar_{ij}}+f^c(r_{ij}) \left( 1+\frac{Q}{r} \right) A a e^{-ar_{ij}}\frac{\vec{r}_{ij}}{r_{ij}}.
 +
</math>
 +
 +
<math>
 +
    V^A (r) = f^c (r) \sum_{n = 1,3} B_n e^{-\beta_n r},
 +
</math>
 +
 +
где
 +
 +
<math>
 +
    f^c (r) = \left\{
 +
    \begin{array}{l}
 +
        1,  \\
 +
        \left[ 1 + \cos(\pi(r - D_{\min}) / (D_{\max} - D_{\min})) \right] / 2,\\
 +
        0, \\
 +
    \end{array} \right.
 +
    \begin{array}{l}
 +
        r < D_{\min}, \\
 +
        D_{\min} < r < D_{\max}, \\
 +
        r > D_{\max},
 +
    \end{array}
 +
</math>
 +
 +
<math>
 +
\frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} = \left\{
 +
\begin{array}{l}
 +
\vec{0}, \ \ \ \ r_{ij} < D_{min} \\
 +
\sin{\left( \pi \frac{(r_{ij}-D_{min})}{(D_{max} - D_{min})} \right)} \cdot \frac{\pi \vec{r}_{ij}}{(D_{min}- D_{max})r_{ij}}, \ \ \ \ \mbox{при} D_{min} < r_{ij} < D_{max} \\
 +
\vec{0}, \ \ \ \ r_{ij} > D_{max}
 +
\end{array} \right.
 +
</math>
 +
 +
<math>
 +
\frac{\partial V^A(r_{ij})}{\partial \vec{r}_i} = \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} \cdot \sum \limits_{n=1,3} B_n e^{-\beta_nr_{ij}} + f^c(r_{ij}) \cdot \sum \limits_{n=1,3}B_n \beta_n e^{-\beta_nr_{ij}}\frac{\vec{r}_{ij}}{r_{ij}}.
 +
</math>
 +
 +
 +
 +
Параметры имеют вид:
 +
 +
<math>
 +
    \begin{array}{l}
 +
    B_1  = 12 388.791 977 98 \,\mbox{eV},\; \beta_1  = 4.720 452 3127 \,\mbox{Å}^{-1},\; Q = 0.313 460 296
 +
    0833 \,\mbox{Å},\\
 +
    B_2  = 17.567 406 465 09 \,\mbox{eV},\; \beta_2 = 1.433 213 2499 \,\mbox{Å}^{-1},\; A = 10 953.544 162 170
 +
    \,\mbox{eV},\\
 +
    B_3  = 30.714 932 080 65 \,\mbox{eV},\; \beta_3  = 1.382 691 2506 \,Å^{-1},\; \alpha = 4.746 539 060 6595
 +
    \,\mbox{Å}^{-1},\\
 +
    D_{\min}  = 1.7 \,\mbox{Å},\; D_{\max}  = 2.0 \,\mbox{Å}.
 +
  \end{array}
 +
</math>
 +
 +
Множитель <math>b_{ij}</math> равен <math>b_{ij} = (B_{ij} + B_{ji}) / 2</math>
 +
 +
а, соответственно его производная
 +
<math>
 +
\frac{\partial b_{ij}}{\partial \vec{r}_i} = \frac{1}{2} \left( \frac{\partial B_{ij}}{\vec{r}_i}+ \frac{\partial B_{ji}}{\partial \vec{r}_i} \right)
 +
</math>, где
 +
 +
<math>
 +
    B_{ij} = \left[ 1 + \sum_{k (\neq i, j)} f^c (r_{ik}) G(\cos(\theta_{ijk}))
 +
    \right]^{-1/2},
 +
</math>
 +
 +
А производная<math>B_{ij}</math> считается по следующей формуле
 +
 +
<math>
 +
\frac{\partial B_{ij}}{\partial \vec{r}_i} = \sum \limits_{k \neq i,j} \left[ \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} \cdot G(\cos{(\theta_{ijk})}) + f^c(r_{ik}) \cdot \frac{\partial G(\cos{(\theta_{ijk})})}{\partial \vec{r}_i} \right] .
 +
</math>
 +
 +
<math>
 +
\frac{\partial G(\cos{\theta_{ijk}})}{\partial \vec{r}_i} = \frac{\partial G(\cos{\theta_{ijk}})}{\partial \cos{\theta_{ijk}}} \cdot \frac{\partial \cos{\theta_{ijk}}}{\partial \vec{r}_i}.
 +
</math>
 +
 +
где <math>\theta_{ijk}</math> – угол между связями, соединяющими атомы
 +
<math>i,j</math> и <math>i,k</math>. Функция <math>G(\cos\theta)</math> строится как полином через значения функции и ее
 +
производных в точках, соответствующих равновесным конфигурациям алмаза (<math>\theta =
 +
\arccos(-1/3)</math>) и графена (<math>\theta = 2 \pi / 3</math>):
 +
 +
{| class="wikitable"
 +
|  <math>\theta(rad)</math>
 +
|  <math>G(\cos \theta)</math>
 +
|  <math>dG(\cos \theta) / d\cos \theta</math>
 +
|  <math>d^2 G(\cos \theta) / d\cos \theta^2</math>
 +
|-
 +
|  <math>0.6082\pi</math>
 +
|  <math>0.097 33</math>
 +
|  <math>0.400 00</math>
 +
|  <math>1.980 00</math>
 +
|-
 +
|  <math>2\pi / 3</math>
 +
|  <math>0.052 80</math>
 +
|  <math>0.170 00</math>
 +
|  <math>0.370 00</math>
 +
|}
 +
 +
Производные от косинуса по радиус-векторам i-ой и j-ой частицы  высчитываются так (где i – вершина угла):
 +
 +
<math>
 +
\frac{\partial \cos{\theta_{ijk}}}{\partial \vec{r}_i} = \frac{\vec{r}_{ij} \times \left( \vec{r}_{ij} \times \vec{r}_{ik} \right)}{r_{ij}^3 r_{ik}} + \frac{\vec{r}_{ik} \times \left( \vec{r}_{ik} \times \vec{r}_{ij} \right)}{r_{ik}^3 r_{ij}},\ \ \ \ \mbox{i —- вершина}
 +
</math>
 +
 +
<math>
 +
\frac{\partial \cos{\theta_{ijk}}}{\partial \vec{r}_j} = \frac{\vec{r}_{ij} \times \left( \vec{r}_{ik} \times \vec{r}_{ij} \right)}{r_{ij}^3 r_{ik}}.
 +
</math>
 +
 +
==Потенциал погруженного атома==
 +
Модель погружённого атома (англ. embedded atom model, EAM) используется для приближенного описания энергии взаимодействия между двумя атомами. Полная энергия системы состоит из двух слагаемых – энергии парного взаимодействия атомов и энергии взаимодействия каждого атома с электронной плотностью, создаваемой другими атомами.
 +
 +
Для расчета энергии парного взаимодействия используется следующая формула:
 +
<math>
 +
E(\vec{r}_1,...,\vec{r}_N)= \frac{1}{2} \sum \limits_{i \neq j} \left[ φ(r_{ij}) \right]
 +
</math>
 +
 +
где <math> φ(r_{ij}) </math> −потенциал взаимодействия i−го и j−го атомов, находящихся на расстоянии <math> r_{ij} </math>.
 +
 +
Расчет энергии взаимодействия каждого атома с электронной плотностью, создаваемой другими атомами, идет по формуле
 +
<math>
 +
E = F \sum \limits_{i \neq j} \left[ ρ_{α}(r_{ij}) \right]   
 +
</math>
 +
 +
где <math> r_{ij} </math> — расстояние между i−м и j−м атомами, <math> ρ_{α} </math>  — вклад в плотность заряда электронов от j−го атома в месте расположения i−го атома и F — это функция «погружения», которая представляет энергию,необходимую для помещения i−го атома в электронное облако.
 +
 +
Таким образом, энергия i-го атома равна
 +
<math>
 +
E_{i} = F_{α} (\sum \limits_{i \neq j} \left[ ρ_{α}(r_{ij}) \right]  + \frac{1}{2} \sum \limits_{i \neq j} \left[ φ(r_{ij}) \right]
 +
</math>
 +
 +
Для расчета силы от функции погружения используется дифференцирование по плотности:
 +
<math>
 +
F_{pogr} = − \frac{\partial F_{α}}{\partial ρ} (\sum \limits_{i} \left[ \frac{\partial ρ_{i}}{\partial r_{ij}} \frac{\partial r_{ij}}{\partial r_{i}}  \right] 
 +
</math>
 +
 +
Для расчета силы от парного потенциала используется дифференцирование по расстоянию:
 +
<math>
 +
F_{parn} = \frac{\partial Π}{\partial ρ} \frac{\partial ρ}{\partial R} \frac{\vec{R}}{\partial R} 
 +
</math>
 +
 +
Общая сила равна сумме сил, действующих со стороны обоих потенциалов:
 +
<math>
 +
F = F_{pogr} + F_{parn}
 +
</math>

Текущая версия на 00:19, 25 мая 2018

Виртуальная лаборатория > Фазовые переходы МД

Задача[править]

Переход от кристаллической структуры к газу. В направлении абсцисс используются периодические ГУ, в направлении оси ординат один ряд частиц фиксирован, с другой стороны несколько рядов частиц (3-5) нагреваются посредством термостата Берендсена (регулируемые параметры). Частицы взаимодействуют посредством потенциала. Уравнения движения интегрируются методом Leapfrog. Система забывает об улетевших частицах.

Задача I[править]

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

Задача II[править]

Основные элементы расчетной части: Запуск расчета, создание образца с треугольной решеткой, задание начальных условий, определение связей, интегрирование уравнений движения методом Leapfrog, расчет сил парным потенциалом. Удаление улетевших частиц.

Список Группы:

Задача III[править]

Расчет сил потенциалом Бреннера второго поколения, создание решетки графена, расчет связей, термостат Берендсена.

Список Группы:

Задача IV[править]

Расчет сил потенциалом погруженного атома для Железа. Задание периодических граничных условий.

Список Группы:

Решение задачи[править]

Открывать лучше в Mozile FireFox, либо настраивать аппаратное ускорение самому (если программа не открывается)

Потенциал Бреннера второго поколения[править]

Потенциал Бреннера второго поколения позволяет представить энергию связи в виде

[math] E_b = \sum_i \sum_{j (\gt i)} \left[ V^R (r_{ij}) - b_{ij} V^A (r_{ij}) \right]. [/math]

Силу, действующую на частицу с номером I можно рассчитать как минус градиент энергии (производная по радиус-вектору частицы i)

[math] \vec{F}_i = - \sum \limits_{j \neq i} \left[ \frac{\partial V^R(r_{ij})}{\partial \vec{r}_i}-\frac{\partial b_{ij}}{\partial \vec{r}_i}V^A(r_{ij})-b_{ij}\frac{\partial V^A(r_{ij})}{\partial \vec{r}_i} \right], [/math]

Между атомами углерода функции отталкивания и притяжения имеют вид:

[math] \frac{\partial V^R(r_{ij})}{\partial \vec{r}_i} = \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} \left( 1+\frac{Q}{r} \right) A e^{-ar_{ij}}+f^c(r_{ij}) Q \frac{\vec{r}_{ij}}{r_{ij}^3}Ae^{-ar_{ij}}+f^c(r_{ij}) \left( 1+\frac{Q}{r} \right) A a e^{-ar_{ij}}\frac{\vec{r}_{ij}}{r_{ij}}. [/math]

[math] V^A (r) = f^c (r) \sum_{n = 1,3} B_n e^{-\beta_n r}, [/math]

где

[math] f^c (r) = \left\{ \begin{array}{l} 1, \\ \left[ 1 + \cos(\pi(r - D_{\min}) / (D_{\max} - D_{\min})) \right] / 2,\\ 0, \\ \end{array} \right. \begin{array}{l} r \lt D_{\min}, \\ D_{\min} \lt r \lt D_{\max}, \\ r \gt D_{\max}, \end{array} [/math]

[math] \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} = \left\{ \begin{array}{l} \vec{0}, \ \ \ \ r_{ij} \lt D_{min} \\ \sin{\left( \pi \frac{(r_{ij}-D_{min})}{(D_{max} - D_{min})} \right)} \cdot \frac{\pi \vec{r}_{ij}}{(D_{min}- D_{max})r_{ij}}, \ \ \ \ \mbox{при} D_{min} \lt r_{ij} \lt D_{max} \\ \vec{0}, \ \ \ \ r_{ij} \gt D_{max} \end{array} \right. [/math]

[math] \frac{\partial V^A(r_{ij})}{\partial \vec{r}_i} = \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} \cdot \sum \limits_{n=1,3} B_n e^{-\beta_nr_{ij}} + f^c(r_{ij}) \cdot \sum \limits_{n=1,3}B_n \beta_n e^{-\beta_nr_{ij}}\frac{\vec{r}_{ij}}{r_{ij}}. [/math]


Параметры имеют вид:

[math] \begin{array}{l} B_1 = 12 388.791 977 98 \,\mbox{eV},\; \beta_1 = 4.720 452 3127 \,\mbox{Å}^{-1},\; Q = 0.313 460 296 0833 \,\mbox{Å},\\ B_2 = 17.567 406 465 09 \,\mbox{eV},\; \beta_2 = 1.433 213 2499 \,\mbox{Å}^{-1},\; A = 10 953.544 162 170 \,\mbox{eV},\\ B_3 = 30.714 932 080 65 \,\mbox{eV},\; \beta_3 = 1.382 691 2506 \,Å^{-1},\; \alpha = 4.746 539 060 6595 \,\mbox{Å}^{-1},\\ D_{\min} = 1.7 \,\mbox{Å},\; D_{\max} = 2.0 \,\mbox{Å}. \end{array} [/math]

Множитель [math]b_{ij}[/math] равен [math]b_{ij} = (B_{ij} + B_{ji}) / 2[/math]

а, соответственно его производная [math] \frac{\partial b_{ij}}{\partial \vec{r}_i} = \frac{1}{2} \left( \frac{\partial B_{ij}}{\vec{r}_i}+ \frac{\partial B_{ji}}{\partial \vec{r}_i} \right) [/math], где

[math] B_{ij} = \left[ 1 + \sum_{k (\neq i, j)} f^c (r_{ik}) G(\cos(\theta_{ijk})) \right]^{-1/2}, [/math]

А производная[math]B_{ij}[/math] считается по следующей формуле

[math] \frac{\partial B_{ij}}{\partial \vec{r}_i} = \sum \limits_{k \neq i,j} \left[ \frac{\partial f^c(r_{ij})}{\partial \vec{r}_i} \cdot G(\cos{(\theta_{ijk})}) + f^c(r_{ik}) \cdot \frac{\partial G(\cos{(\theta_{ijk})})}{\partial \vec{r}_i} \right] . [/math]

[math] \frac{\partial G(\cos{\theta_{ijk}})}{\partial \vec{r}_i} = \frac{\partial G(\cos{\theta_{ijk}})}{\partial \cos{\theta_{ijk}}} \cdot \frac{\partial \cos{\theta_{ijk}}}{\partial \vec{r}_i}. [/math]

где [math]\theta_{ijk}[/math] – угол между связями, соединяющими атомы [math]i,j[/math] и [math]i,k[/math]. Функция [math]G(\cos\theta)[/math] строится как полином через значения функции и ее производных в точках, соответствующих равновесным конфигурациям алмаза ([math]\theta = \arccos(-1/3)[/math]) и графена ([math]\theta = 2 \pi / 3[/math]):

[math]\theta(rad)[/math] [math]G(\cos \theta)[/math] [math]dG(\cos \theta) / d\cos \theta[/math] [math]d^2 G(\cos \theta) / d\cos \theta^2[/math]
[math]0.6082\pi[/math] [math]0.097 33[/math] [math]0.400 00[/math] [math]1.980 00[/math]
[math]2\pi / 3[/math] [math]0.052 80[/math] [math]0.170 00[/math] [math]0.370 00[/math]

Производные от косинуса по радиус-векторам i-ой и j-ой частицы высчитываются так (где i – вершина угла):

[math] \frac{\partial \cos{\theta_{ijk}}}{\partial \vec{r}_i} = \frac{\vec{r}_{ij} \times \left( \vec{r}_{ij} \times \vec{r}_{ik} \right)}{r_{ij}^3 r_{ik}} + \frac{\vec{r}_{ik} \times \left( \vec{r}_{ik} \times \vec{r}_{ij} \right)}{r_{ik}^3 r_{ij}},\ \ \ \ \mbox{i —- вершина} [/math]

[math] \frac{\partial \cos{\theta_{ijk}}}{\partial \vec{r}_j} = \frac{\vec{r}_{ij} \times \left( \vec{r}_{ik} \times \vec{r}_{ij} \right)}{r_{ij}^3 r_{ik}}. [/math]

Потенциал погруженного атома[править]

Модель погружённого атома (англ. embedded atom model, EAM) используется для приближенного описания энергии взаимодействия между двумя атомами. Полная энергия системы состоит из двух слагаемых – энергии парного взаимодействия атомов и энергии взаимодействия каждого атома с электронной плотностью, создаваемой другими атомами.

Для расчета энергии парного взаимодействия используется следующая формула: [math] E(\vec{r}_1,...,\vec{r}_N)= \frac{1}{2} \sum \limits_{i \neq j} \left[ φ(r_{ij}) \right] [/math]

где [math] φ(r_{ij}) [/math] −потенциал взаимодействия i−го и j−го атомов, находящихся на расстоянии [math] r_{ij} [/math].

Расчет энергии взаимодействия каждого атома с электронной плотностью, создаваемой другими атомами, идет по формуле [math] E = F \sum \limits_{i \neq j} \left[ ρ_{α}(r_{ij}) \right] [/math]

где [math] r_{ij} [/math] — расстояние между i−м и j−м атомами, [math] ρ_{α} [/math] — вклад в плотность заряда электронов от j−го атома в месте расположения i−го атома и F — это функция «погружения», которая представляет энергию,необходимую для помещения i−го атома в электронное облако.

Таким образом, энергия i-го атома равна [math] E_{i} = F_{α} (\sum \limits_{i \neq j} \left[ ρ_{α}(r_{ij}) \right] + \frac{1}{2} \sum \limits_{i \neq j} \left[ φ(r_{ij}) \right] [/math]

Для расчета силы от функции погружения используется дифференцирование по плотности: [math] F_{pogr} = − \frac{\partial F_{α}}{\partial ρ} (\sum \limits_{i} \left[ \frac{\partial ρ_{i}}{\partial r_{ij}} \frac{\partial r_{ij}}{\partial r_{i}} \right] [/math]

Для расчета силы от парного потенциала используется дифференцирование по расстоянию: [math] F_{parn} = \frac{\partial Π}{\partial ρ} \frac{\partial ρ}{\partial R} \frac{\vec{R}}{\partial R} [/math]

Общая сила равна сумме сил, действующих со стороны обоих потенциалов: [math] F = F_{pogr} + F_{parn} [/math]