Skip to content

Commit 552d4f3

Browse files
committed
prepared umbrella sampling
1 parent ac66598 commit 552d4f3

3 files changed

Lines changed: 103 additions & 20 deletions

File tree

Lines changed: 39 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,39 @@
1+
2+
variable sigma equal 3.405
3+
variable epsilon equal 0.238
4+
variable U0 equal 10*${epsilon}
5+
variable dlt equal 1.0
6+
variable x0 equal 10.0
7+
8+
units real
9+
atom_style atomic
10+
pair_style lj/cut 3.822
11+
pair_modify shift yes
12+
boundary p p p
13+
14+
region myreg block -50 50 -15 15 -50 50
15+
create_box 2 myreg
16+
create_atoms 2 single 0 0 0
17+
create_atoms 1 random 199 34134 myreg overlap 3 maxtry 50
18+
19+
mass * 39.948
20+
pair_coeff * * ${epsilon} ${sigma}
21+
group topull type 2
22+
23+
variable U atom ${U0}*atan((x+${x0})/${dlt})-${U0}*atan((x-${x0})/${dlt})
24+
variable F atom ${U0}/((x-${x0})^2/${dlt}^2+1)/${dlt}-${U0}/((x+${x0})^2/${dlt}^2+1)/${dlt}
25+
fix myadf all addforce v_F 0.0 0.0 energy v_U
26+
27+
fix mynve all nve
28+
fix mylgv all langevin 119.8 119.8 500 30917
29+
30+
thermo 2000
31+
32+
dump viz all image 50000 myimage-*.ppm type type &
33+
shiny 0.1 box yes 0.01 view 180 90 zoom 6 &
34+
size 1600 500 fsaa yes
35+
dump_modify viz backcolor white acolor 1 cyan &
36+
acolor 2 red adiam 1 3 adiam 2 3 boxcolor black
37+
38+
timestep 2.0
39+
run 50000
Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,13 @@
1+
2+
variable sigma equal 3.405
3+
variable epsilon equal 0.238
4+
variable U0 equal 10*${epsilon}
5+
variable dlt equal 1.0
6+
variable x0 equal 10.0
7+
8+
units real
9+
atom_style atomic
10+
pair_style lj/cut 3.822
11+
pair_modify shift yes
12+
boundary p p p
13+

lammps-tutorials.tex

Lines changed: 51 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -3504,7 +3504,7 @@ \subsubsection{Adding water}
35043504
\centering
35053505
\includegraphics[width=\linewidth]{GCMC-number}
35063506
\caption{Number of water molecules $N_\text{H2O}$ as a function of the time $t$
3507-
as extracgted from \hyperref[gcmc-silica-label]{Tutorial 6}.}
3507+
as extracted from \hyperref[gcmc-silica-label]{Tutorial 6}.}
35083508
\label{fig:GCMC-number}
35093509
\end{figure}
35103510

@@ -3633,7 +3633,7 @@ \subsubsection{Method 1: Free sampling}
36333633
\centering
36343634
\includegraphics[width=\linewidth]{US-potential}
36353635
\caption{Potential $U$ given in Eq.\,\eqref{eq:U} (a) and force $F$ given in
3636-
Eq.\,\eqref{eq:F} (b) as a function of the coordinate $x$. Here,
3636+
Eq.\,\eqref{eq:F} (b) as functions of the coordinate $x$. Here,
36373637
$U_0 = 0.36~\text{kcal/mol}$, $\delta = 1.0~\mathrm{\AA{}}$, and $x_0 = 10~\mathrm{\AA{}}$.}
36383638
\label{fig:potential}
36393639
\end{figure}
@@ -3774,43 +3774,74 @@ \subsubsection{Method 2: Umbrella sampling}
37743774
unbiased free energy profile.
37753775

37763776
\paragraph{LAMMPS input script}
3777-
Create a new folder called \flrcmd{BiasedSampling/}, and create a new input file
3778-
called \flecmd{input.lmp} in it, and copy the following lines:
3777+
3778+
Open the file named \lmpcmd{umbrella-sampling.lmp}, which should
3779+
contain the following lines:
37793780
\begin{lstlisting}
37803781
variable sigma equal 3.405
37813782
variable epsilon equal 0.238
37823783
variable U0 equal 10*${epsilon}
3783-
variable dlt equal 0.5
3784-
variable x0 equal 5.0
3784+
variable dlt equal 1.0
3785+
variable x0 equal 10
37853786
variable k equal 1.5
37863787

37873788
units real
37883789
atom_style atomic
37893790
pair_style lj/cut 3.822
37903791
pair_modify shift yes
37913792
boundary p p p
3793+
\end{lstlisting}
3794+
The only difference with the previous case is the larger value
3795+
for the repulsive potential $U_0$, which will make the central area
3796+
of the system very unliky to be visited by free particles.
37923797

3793-
region myreg block -25 25 -5 5 -25 25
3798+
Let us create a box with 2 atom types, with a single particle of type 2,
3799+
by adding the following lines to \lmpcmd{umbrella-sampling.lmp}:
3800+
\begin{lstlisting}
3801+
region myreg block -50 50 -15 15 -50 50
37943802
create_box 2 myreg
3795-
labelmap atom 1 A1 2 A2
3796-
create_atoms A2 single 0 0 0
3797-
create_atoms A1 random 5 34134 myreg &
3798-
overlap 1.5 maxtry 50
3799-
3803+
create_atoms 2 single 0 0 0
3804+
create_atoms 1 random 199 34134 myreg overlap 3 maxtry 50
3805+
\end{lstlisting}
3806+
Then, let us give the same mass and LJ parameters to both atoms of types
3807+
1 and 2, and let us place the atoms of type 2 into a group named \lmpcmd{topull}:
3808+
\begin{lstlisting}
38003809
mass * 39.948
38013810
pair_coeff * * ${epsilon} ${sigma}
3802-
neigh_modify every 1 delay 4 check yes
38033811
group topull type 2
3804-
3805-
minimize 1e-4 1e-6 100 1000
3806-
reset_timestep 0
3807-
3812+
\end{lstlisting}
3813+
Then, the same potential $U$ and force $F$ is applied to all the atoms,
3814+
toguether with the same \lmpcmd{fix nve} and \lmpcmd{fix langevin}:
3815+
\begin{lstlisting}
38083816
variable U atom ${U0}*atan((x+${x0})/${dlt})&
38093817
-${U0}*atan((x-${x0})/${dlt})
3810-
variable F atom &
3811-
${U0}/((x-${x0})^2/${dlt}^2+1)/${dlt}&
3818+
variable F atom ${U0}/((x-${x0})^2/${dlt}^2+1)/${dlt}&
38123819
-${U0}/((x+${x0})^2/${dlt}^2+1)/${dlt}
3813-
fix pot all addforce v_F 0.0 0.0 energy v_U
3820+
fix myadf all addforce v_F 0.0 0.0 energy v_U
3821+
3822+
fix mynve all nve
3823+
fix mylgv all langevin 119.8 119.8 500 30917
3824+
\end{lstlisting}
3825+
And then let us perform a short equilibration in the nvt ensemble
3826+
to prepare for the umbrella sampling run:
3827+
\begin{lstlisting}
3828+
thermo 2000
3829+
3830+
dump viz all image 50000 myimage-*.ppm type type &
3831+
shiny 0.1 box yes 0.01 view 180 90 zoom 6 &
3832+
size 1600 500 fsaa yes
3833+
dump_modify viz backcolor white acolor 1 cyan &
3834+
acolor 2 red adiam 1 3 adiam 2 3 boxcolor black
3835+
3836+
timestep 2.0
3837+
run 50000
3838+
\end{lstlisting}
3839+
3840+
3841+
\clearpage
3842+
3843+
3844+
\begin{lstlisting}
38143845

38153846
fix mynve all nve
38163847
fix mylgv all langevin 119.8 119.8 50 1530917

0 commit comments

Comments
 (0)