forked from google-deepmind/deepmind-research
-
Notifications
You must be signed in to change notification settings - Fork 0
/
lammps.dat
129 lines (104 loc) · 5.3 KB
/
lammps.dat
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
# Simulate a WCA solute particle immersed into a Lennard-Jones solvent.
#
# We first generate an initial configuration of particles located at random
# positions. We then minimise the system and equilibrate it to the target
# temperature. During the production run, we sample configurations at regular
# intervals.
################################################################################
# Specify input parameters.
################################################################################
# We use the same box dimensions and temperature as in Jarzynski (2002),
# Phys. Rev. E 65, 046122, but convert these quantities to reduced units.
# T=300 K, eps=0.1854 kcal/mol, k_boltzmann=0.001987204 kcal/mol/K.
variable solute_radius equal 9.2/3.542
variable box_length equal 22.28/3.542
variable temperature_reference equal 300*0.001987204/0.1854
variable temperature_damp equal 0.5
variable seed equal 1234
variable num_solvent_particles equal 125
variable dump_frequency equal 5000
variable timestep equal 0.002
variable runtime_equi equal 50000
variable runtime_prod equal 500000
variable steps_equilibration equal floor(${runtime_equi}/${timestep})
variable steps_production equal floor(${runtime_prod}/${timestep})
################################################################################
# Create solvent particles at random positions and place solute at the origin.
################################################################################
units lj
atom_style atomic
dimension 3
boundary p p p
variable box_length_half equal ${box_length}/2.0
region domain block -${box_length_half} ${box_length_half} &
-${box_length_half} ${box_length_half} &
-${box_length_half} ${box_length_half} &
units box
create_box 2 domain
create_atoms 2 random 1 ${seed} domain
create_atoms 1 random ${num_solvent_particles} ${seed} domain
set type 2 x 0.0 y 0.0 z 0.0
group solute type 2
group solvent type 1
# Initialise particle masses.
mass 1 1
mass 2 10000000
################################################################################
# Define pair-style.
################################################################################
variable sigma_solute equal ${solute_radius}
variable radius_cutoff_solute equal 1.12246*${sigma_solute}
variable radius_cutoff equal ${box_length}/2.0
pair_style lj/cut ${radius_cutoff}
# Lennard-Jones interaction between solvent particles.
pair_coeff 1 1 1.0 1.0
# WCA interaction between solute and solvent particles.
pair_coeff 1 2 1.0 ${sigma_solute} ${radius_cutoff_solute}
# No interaction between solute particles.
pair_coeff 2 2 0.0 0.0
# Shift energies to be zero at the cutoff and update neighbor list settings.
pair_modify shift yes
neighbor 0.3 bin
neigh_modify delay 5
################################################################################
# Perform minimisation.
################################################################################
fix freeze solute setforce 0.0 0.0 0.0
minimize 1.0e-4 1.0e-6 100 1000
################################################################################
# Perform equilibration run.
################################################################################
reset_timestep 0
compute temperature_solvent solvent temp
compute kinetic_energy all ke
compute potential_energy all pe
variable total_energy equal c_kinetic_energy+c_potential_energy
# Compute centre-of-mass velocity and solute position for monitoring.
variable vcmx equal "vcm(all,x)"
variable vcmy equal "vcm(all,y)"
variable vcmz equal "vcm(all,z)"
variable vcm2 equal v_vcmx*v_vcmx+v_vcmy*v_vcmy+v_vcmz*v_vcmz
compute position_solute solute com
# Specify terminal output.
thermo_style custom step temp c_temperature_solvent c_potential_energy &
c_kinetic_energy v_total_energy press v_vcm2 &
c_position_solute[1] c_position_solute[2] c_position_solute[3]
thermo_modify norm no
thermo 1000
# Specify integration timestep.
timestep ${timestep}
# Apply Langevin thermostat to solvent particles.
fix flangevin solvent langevin ${temperature_reference} &
${temperature_reference} ${temperature_damp} ${seed} zero yes
fix fnve all nve
run ${steps_equilibration}
################################################################################
# Production run.
################################################################################
reset_timestep 0
dump fDumpTrajectory all custom ${dump_frequency} &
trajectory.dat id type x y z
dump_modify fDumpTrajectory sort id format float %.8g
fix fDumpEnergy all ave/time ${dump_frequency} 1 ${dump_frequency} &
c_potential_energy file energy.dat format %.8g
run ${steps_production}