1-
21Prepare and relax
32=================
43
54
65The first step is to relax the structure with ReaxFF, which which will be achieved using
76molecular dynamics. To ensure the system equilibrates properly, we will monitor certain
87parameters over time, such as the system volume. To set up this
9- tutorial, select * Start Tutorial 5 * from the
10- * Tutorials * menu of LAMMPS--GUI and follow the instructions.
11- The editor should display the following content corresponding to *relax.lmp *:
8+ tutorial, select `` Start Tutorial 5 `` from the
9+ `` Tutorials `` menu of LAMMPS--GUI and follow the instructions.
10+ The editor should display the following content corresponding to ** relax.lmp * *:
1211
1312.. code-block :: lammps
1413
@@ -17,15 +16,14 @@ The editor should display the following content corresponding to *relax.lmp*:
1716
1817 read_data silica.data
1918
20-
2119 So far, the input is very similar to what was seen in the previous tutorials.
22- Some basic parameters are defined (* units * and * atom_style * ),
23- and a *.data * file is imported by the * read \ data * command.
20+ Some basic parameters are defined (`` units `` and `` atom_style `` ),
21+ and a ** .data ** file is imported by the `` read_data `` command.
2422
2523The initial topology given by \h ref{\f ilepath tutorial5/silica.data}{\d wlcmd{silica.data}}
2624is a small amorphous silica structure. This structure was created using a force field called
27- Vashishta :cite: `vashishta1990interaction `. If you open the *silica.data *
28- file, you will find in the * Atoms * section that all silicon atoms have a
25+ Vashishta :cite: `vashishta1990interaction `. If you open the ** silica.data * *
26+ file, you will find in the `` Atoms `` section that all silicon atoms have a
2927charge of :math: `q = 1.1 \,\text {e}`, and all oxygen atoms have a charge of :math: `q = -0.55 \,\text {e}`.
3028
3129.. admonition :: Note
@@ -36,28 +34,28 @@ charge of :math:`q = 1.1\,\text{e}`, and all oxygen atoms have a charge of :math
3634 changes once ReaxFF is used: the charge of each atom will adjust to its local
3735 environment.
3836
39- Next, copy the following three crucial lines into the *relax.lmp * file:
37+ Next, copy the following three crucial lines into the ** relax.lmp * * file:
4038
4139.. code-block :: lammps
4240
4341 pair_style reaxff NULL safezone 3.0 mincap 150
4442 pair_coeff * * reaxCHOFe.inc Si O
4543 fix myqeq all qeq/reaxff 1 0.0 10.0 1.0e-6 reaxff maxiter 400
4644
47- In this case, the * pair \ style reaxff* is used without a control file. The
48- * safezone * and * mincap * keywords are added to prevent
45+ In this case, the `` pair_style reaxff`` is used without a control file. The
46+ `` safezone `` and `` mincap `` keywords are added to prevent
4947allocation issues, which sometimes can trigger segmentation faults and
50- * bondchk * errors. The * pair \ coeff * command uses the
48+ `` bondchk `` errors. The `` pair_coeff `` command uses the
5149\h ref{\f ilepath tutorial5/reaxCHOFe.inc}{\d wlcmd{reaxCHOFe.inc}}
5250file, which should have been downloaded during the tutorial set up. Finally, the
53- * fix qeq/reaxff * is used to perform charge equilibration :cite: `rappe1991charge `,
51+ `` fix qeq/reaxff `` is used to perform charge equilibration :cite: `rappe1991charge `,
5452which occurs at every step. The values 0.0 and 10.0 represent the
5553low and the high cutoffs, respectively, and :math: `1.0 \text {e} -6 ` is the tolerance.
56- The * maxiter * sets an upper limit to the number of attempts to
54+ The `` maxiter `` sets an upper limit to the number of attempts to
5755equilibrate the charge.
5856
5957
60- Next, add the following commands to the *relax.lmp * file to track the
58+ Next, add the following commands to the ** relax.lmp * * file to track the
6159evolution of the charges during the simulation:
6260
6361.. code-block :: lammps
@@ -68,20 +66,18 @@ evolution of the charges during the simulation:
6866 variable qO equal charge(grpO)/count(grpO)
6967 variable vq atom q
7068
71- To print the averaged charges * qSi * and * qO * using the
72- * thermo_style * command, and create images of the system. Add the
73- following lines to *relax.lmp *:
69+ To print the averaged charges `` qSi `` and `` qO `` using the
70+ `` thermo_style `` command, and create images of the system. Add the
71+ following lines to ** relax.lmp * *:
7472
7573.. code-block :: lammps
7674
7775 thermo 100
7876 thermo_style custom step temp etotal press vol v_qSi v_qO
79- dump viz all image 100 myimage-*.ppm q &
80- type shiny 0.1 box no 0.01 view 180 90 zoom 2.3 size 1200 500
81- dump_modify viz adiam Si 2.6 adiam O 2.3 backcolor white &
82- amap -1 2 ca 0.0 3 min royalblue 0 green max orangered
77+ dump viz all image 100 myimage-*.ppm q type shiny 0.1 box no 0.01 view 180 90 zoom 2.3 size 1200 500
78+ dump_modify viz adiam Si 2.6 adiam O 2.3 backcolor white amap -1 2 ca 0.0 3 min royalblue 0 green max orangered
8379
84- Here, the atoms are colored by their charges * q * , ranging from royal blue
80+ Here, the atoms are colored by their charges `` q `` , ranging from royal blue
8581(when :math: `q=-1 \,\text {e}`) to orange-red (when :math: `q=2 \,\text {e}`).
8682
8783ADD FIGURE SIO-deformed Amorphous silicon oxide after deformation.
@@ -91,24 +87,22 @@ about :math:`1.8~\text{e}` appear in red/orange, and bulk O atoms with a charge
9187about :math: `-0.9 ~ \text {e}` appear in blue.
9288
9389We can generate histograms of the charges for each atom type using
94- * fix ave/histo * commands:
90+ `` fix ave/histo `` commands:
9591
9692.. code-block :: lammps
9793
98- fix myhis1 grpSi ave/histo 10 500 5000 -1.5 2.5 1000 v_vq &
99- file relax-Si.histo mode vector
100- fix myhis2 grpO ave/histo 10 500 5000 -1.5 2.5 1000 v_vq &
101- file relax-O.histo mode vector
94+ fix myhis1 grpSi ave/histo 10 500 5000 -1.5 2.5 1000 v_vq file relax-Si.histo mode vector
95+ fix myhis2 grpO ave/histo 10 500 5000 -1.5 2.5 1000 v_vq file relax-O.histo mode vector
10296
103- We can also use the * fix reaxff/species * to evaluate what species are
97+ We can also use the `` fix reaxff/species `` to evaluate what species are
10498present within the simulation. It will be useful later when the system is deformed,
10599and bonds are broken:
106100
107101.. code-block :: lammps
108102
109103 fix myspec all reaxff/species 5 1 5 relax.species element Si O
110104
111- Here, the information will be printed every 5 steps in a file called *relax.species *.
105+ Here, the information will be printed every 5 steps in a file called ** relax.species * *.
112106Let us perform a very short run using the anisotropic NPT command and relax the
113107density of the system:
114108
@@ -122,8 +116,8 @@ density of the system:
122116
123117 write_data relax.data
124118
125- Run the *relax.lmp * file using LAMMPS. As seen from *relax.species *,
126- only one species is detected, called * O384Si192 * , representing the entire system.
119+ Run the ** relax.lmp ** file using LAMMPS. As seen from ** relax.species * *,
120+ only one species is detected, called `` O384Si192 `` , representing the entire system.
127121
128122As the simulation progresses, the charge of every atom fluctuates
129123because it is adjusting to the local environment of the atom (Fig.~\r ef{fig:SIO-charge}\, a).
@@ -133,7 +127,7 @@ to a rapid change in the system volume, which causes interatomic distances to
133127shift quickly (Fig.~\r ef{fig:SIO-charge}\, b). The atoms with the
134128most extreme charges are located at structural defects,
135129such as dangling oxygen groups (Fig.~\r ef{fig:SIO-slice}).
136- Finally, the generated *.histo * files can be used to
130+ Finally, the generated ** .histo * * files can be used to
137131plot the probability distributions, :math: `P(q)` (see Fig.~\r ef{fig:SIO-distribution}\, a).
138132
139133FIGURE SIO-charge a) Average charge per atom of the silicon, :math: `q_\text {Si}`, atoms as
@@ -155,7 +149,7 @@ Deform the structure
155149
156150
157151Let us apply a deformation to the structure to force some :math: `\text {Si}-\text {O}`
158- bonds to break (and eventually re-assemble). Open the *deform.lmp *
152+ bonds to break (and eventually re-assemble). Open the ** deform.lmp * *
159153file, which must contain the following lines:
160154
161155.. code-block :: lammps
@@ -177,21 +171,17 @@ file, which must contain the following lines:
177171
178172 thermo 200
179173 thermo_style custom step temp etotal press vol v_qSi v_qO
180- dump viz all image 100 myimage-*.ppm q &
181- type shiny 0.1 box no 0.01 view 180 90 zoom 2.3 size 1200 500
182- dump_modify viz adiam Si 2.6 adiam O 2.3 backcolor white &
183- amap -1 2 ca 0.0 3 min royalblue 0 green max orangered
184-
185- fix myhis1 grpSi ave/histo 10 500 5000 -1.5 2.5 1000 v_vq &
186- file deform-Si.histo mode vector
187- fix myhis2 grpO ave/histo 10 500 5000 -1.5 2.5 1000 v_vq &
188- file deform-O.histo mode vector
174+ dump viz all image 100 myimage-*.ppm q type shiny 0.1 box no 0.01 view 180 90 zoom 2.3 size 1200 500
175+ dump_modify viz adiam Si 2.6 adiam O 2.3 backcolor white amap -1 2 ca 0.0 3 min royalblue 0 green max orangered
176+
177+ fix myhis1 grpSi ave/histo 10 500 5000 -1.5 2.5 1000 v_vq file deform-Si.histo mode vector
178+ fix myhis2 grpO ave/histo 10 500 5000 -1.5 2.5 1000 v_vq file deform-O.histo mode vector
189179 fix myspec all reaxff/species 5 1 5 deform.species element Si O
190180
191- The only difference with the previous *relax.lmp * file is the path to
192- the *relax.data * file.
181+ The only difference with the previous ** relax.lmp * * file is the path to
182+ the ** relax.data * * file.
193183
194- Next, let us use * fix nvt * instead of * fix npt * to apply a
184+ Next, let us use `` fix nvt `` instead of `` fix npt `` to apply a
195185Nosé-Hoover thermostat without a barostat:
196186
197187.. code-block :: lammps
@@ -200,7 +190,7 @@ Nosé-Hoover thermostat without a barostat:
200190 timestep 0.5
201191
202192 Here, no barostat is used because the change in the box volume will be imposed
203- by the * fix deform * .
193+ by the `` fix deform `` .
204194
205195ADD FIGURE SIO-deformed-charge a) Average charge per atom of the silicon, :math: `q_\text {Si}`, atoms as
206196a function of time, :math: `t`, during deformation of the :math: `\text {SiO}_2 ` system.
@@ -209,9 +199,9 @@ silica structure occurs near :math:`t = 11`\,ps. b) Temperature, :math:`T`, of
209199system as a function of :math: `t`.
210200
211201
212- Let us run for 5000 steps without deformation, then apply the * fix deform *
202+ Let us run for 5000 steps without deformation, then apply the `` fix deform ``
213203to progressively elongate the box along the :math: `x`-axis during 25000 steps. Add
214- the following line to *deform.lmp *:
204+ the following line to ** deform.lmp * *:
215205
216206.. code-block :: lammps
217207
@@ -223,7 +213,7 @@ the following line to *deform.lmp*:
223213
224214 write_data deform.data
225215
226- Run the *deform.lmp * file using LAMMPS. During the deformation, the charge
216+ Run the ** deform.lmp * * file using LAMMPS. During the deformation, the charge
227217values progressively evolve until the structure eventually breaks down. After the
228218structure breaks down, the charges equilibrate near new average values that differ
229219from the initial averages (Fig.~\r ef{fig:SIO-deformed-charge}\, a). The difference
@@ -236,7 +226,7 @@ the material (Fig.~\ref{fig:SIO-deformed-charge}\,b).
236226You can examine the charge distribution after deformation, as well as during
237227deformation (Fig.~\r ef{fig:SIO-distribution}\, b). As expected, the final
238228charge distribution slightly differs from the previously calculated one. If
239- no new species were formed during the simulation, the *deform.species * file
229+ no new species were formed during the simulation, the ** deform.species * * file
240230should look like this:
241231
242232.. code-block :: lammps
@@ -248,7 +238,7 @@ should look like this:
248238 30000 1 1 1
249239
250240 Sometimes, :math: `\text {O}_2 ` molecules are formed during the deformation. If this occurs,
251- a new column * O2 * appears in the *deform.species * file.
241+ a new column `` O2 `` appears in the ** deform.species * * file.
252242
253243Decorate the surface
254244--------------------
@@ -257,9 +247,9 @@ Under ambient conditions, some of the surface :math:`\text{SiO}_2` atoms become
257247passivated by forming covalent bonds with hydrogen (H) atoms :cite: `sulpizi2012silica `.
258248We will add hydrogen atoms randomly to the cracked silica and observe how the
259249system evolves. To do so, we first need to modify the previously generated data
260- file *deform.data * and make space for a third atom type.
261- Copy *deform.data *, name the copy *deform-mod.data *, and modify the
262- first lines of *deform-mod.data * as follows:
250+ file ** deform.data * * and make space for a third atom type.
251+ Copy ** deform.data ** , name the copy ** deform-mod.data * *, and modify the
252+ first lines of ** deform-mod.data * * as follows:
263253
264254.. code-block :: lammps
265255
@@ -282,7 +272,7 @@ first lines of *deform-mod.data* as follows:
282272
283273 (...)
284274
285- Open the *decorate.lmp * file, which must contain the following lines:
275+ Open the ** decorate.lmp * * file, which must contain the following lines:
286276
287277.. code-block :: lammps
288278
@@ -296,12 +286,12 @@ Open the *decorate.lmp* file, which must contain the following lines:
296286 pair_coeff * * reaxCHOFe.inc Si O H
297287 fix myqeq all qeq/reaxff 1 0.0 10.0 1.0e-6 reaxff maxiter 400
298288
299- The * displace_atoms * command is used to move the center of the
289+ The `` displace_atoms `` command is used to move the center of the
300290crack near the center of the box. This step is optional but makes for a nicer
301291visualization. A different value for the shift may be needed in
302292your case, depending on the location of the crack. A difference with the previous
303- input is that three atom types are specified in the * pair_coeff * command, i.e.
304- * Si O H * .
293+ input is that three atom types are specified in the `` pair_coeff `` command, i.e.
294+ `` Si O H `` .
305295
306296ADD FIGURE SIO-decorated Cracked silicon oxide after the addition of hydrogen atoms.
307297The atoms are colored by their charges, with the newly added hydrogen atoms appearing as small
@@ -322,16 +312,12 @@ types of atoms, and output the charge values into log files:
322312 thermo 5
323313 thermo_style custom step temp etotal press v_qSi v_qO v_qH
324314
325- dump viz all image 100 myimage-*.ppm q &
326- type shiny 0.1 box no 0.01 view 180 90 zoom 2.3 size 1200 500
327- dump_modify viz adiam Si 2.6 adiam O 2.3 adiam H 1.0 &
328- backcolor white amap -1 2 ca 0.0 3 min royalblue &
329- 0 green max orangered
315+ dump viz all image 100 myimage-*.ppm q type shiny 0.1 box no 0.01 view 180 90 zoom 2.3 size 1200 500
316+ dump_modify viz adiam Si 2.6 adiam O 2.3 adiam H 1.0 backcolor white amap -1 2 ca 0.0 3 min royalblue 0 green max orangered
330317
331- fix myspec all reaxff/species 5 1 5 decorate.species &
332- element Si O H
318+ fix myspec all reaxff/species 5 1 5 decorate.species element Si O H
333319
334- Here, the :math: `+1 \mathrm {e}{-10 }` was added to the denominator of the * variable qH *
320+ Here, the :math: `+1 \mathrm {e}{-10 }` was added to the denominator of the `` variable qH ``
335321to avoid dividing by 0 at the beginning of the simulation. Finally, let us
336322create a loop with 10 steps, and create two hydrogen atoms at random locations at
337323every step:
@@ -353,7 +339,7 @@ every step:
353339 jump SELF loop
354340
355341 Run the simulation with LAMMPS. When the simulation is over,
356- it can be seen from the *decorate.species * file that
342+ it can be seen from the ** decorate.species * * file that
357343all the created hydrogen atoms reacted with the :math: `\text {SiO}_{2 }` structure to
358344form surface groups (such as hydroxyl (-OH) groups).
359345
0 commit comments