Updated example input scripts, data files and README

This commit is contained in:
Trung Nguyen
2021-06-02 12:03:46 -05:00
parent 1ce50e0b1d
commit b9e9dda1ab
9 changed files with 39 additions and 5817 deletions

View File

@ -2,17 +2,11 @@ This folder contains some example data and input scripts for the USER-DIELECTRIC
Nguyen TD, Li H, Bagchi D, Solis FJ, Olvera de la Cruz, Incorporating surface polarization effects into large-scale coarse-grained molecular dynamics simulation, Computer Physics Communications 2019, 241, 80--91.
- data.minimal : two point opposite charges facing a flat interface
- data.plane : two point opposite charges facing a larger flat interface
- data.confined : two point opposite charges confined between two interfaces (epsilon1=2/epsilon2=10/epsilon2=2)
- data.sphere : two point opposite charges outside a spherical interface (epsilon_in=1/epsilon2=10)
- data.cylinder : single point charge facing a cylindrical interface
- data.confined : two point opposite charges confined between two interfaces (epsilon1=2/epsilon2=10/epsilon2=2)
- data.sphere : two point opposite charges outside a spherical interface (epsilon_in=1/epsilon2=10)
- in.minimal : read in data.minimal
- in.nopbc : read in data.* files, using non-periodic boundary conditions, with a large cutoff
using "-v data data.foo" from the LAMMPS command line arguments
- in.pbc : read in data.* files using periodic boundary conditions
with a slab correction for kspace_style pppm and ewald; not required for msm
- in.confined : read in data.confined
- in.nopbc : read in data.* files, using non-periodic boundary conditions, with a large cutoff
For "atom_style dielectric" the Atoms section in the data file contains 15 following columns:
@ -20,28 +14,31 @@ id mol type q x y z normx normy normz area_per_patch ed em epsilon curvature
where
id, mol, type, q, x, y and z are similar to those in atom_style full
* id, mol, type, q, x, y and z are similar to those in atom_style full
normx, normy and normz are the three components of the normal unit vector
of the interface at the patch (vertex). For real charges (ions), these 3 values
are irrelevant and can be anything (e.g. 0,0,1).
* normx, normy and normz are the three components of the normal unit vector
of the interface at the boundary element (also called vertex, or patch).
For real charges (ions), these 3 values are irrelevant,
and can be anything (e.g. 0,0,1). normx, normy, and normz can be
accessed through mux, muy and muz as if they were dipole components.
ed = dielectric difference at the vertex along the normal vector direction.
For example, if (normx,normy,normz) points from medium with epsilon_in to medium with epsilon_out,
then ed = epsilon_out - epsilon_in
* ed = dielectric difference at the vertex along the normal vector direction.
For example, if (normx,normy,normz) points from medium with epsilon_in
to medium with epsilon_out, then ed = epsilon_out - epsilon_in
em = (epsilon_out + epsilon_in)/2: the mean dielectric value
epsilon = the local epsilon value at the vertex or at the ion.
For real charges, epsilon is the medium dielectric constant, q is the actual real charges.
For interface particles, if q is zero (zero surface charges), epsilon is set to be 1.0;
if q is nonzero (charged surfaces), epsilon is set to be em (the mean dielectric value above).
area_per_patch: surface area of the patch. For real charges, this value is irrelevant, can be 1.0.
curvature: surface curvature at the patch.
For example, for spherical interfaces, curvature = 1/spherical radius.
For planar interfaces, curvature = 0.
* em = (epsilon_out + epsilon_in)/2: the mean dielectric value
* epsilon = the local epsilon value at the vertex or at the ion.
For real charges, epsilon is the medium dielectric constant,
and q is the real (unscaled) charges.
For interface particles:
+ if q is zero (zero surface charges), epsilon is set to be 1.0;
+ if q is nonzero (charged surfaces), epsilon is set to be em
(the mean dielectric value above).
* area_per_patch: the surface area of the patch (element).
For real charges, this value is irrelevant, can be 1.0.
* curvature: surface mean curvature at the patch.
For example, for spherical interfaces, curvature = 1/spherical radius.
For planar interfaces, curvature = 0.

View File

@ -1,4 +1,4 @@
LAMMPS data file: using 2001 atoms for single point charge (no pbc), or 2002 atoms for two opposite charges (to have charge neutrality for pbc)
LAMMPS data file: two oppositely charged ions confined between two walls
4002 atoms
3 atom types

File diff suppressed because it is too large Load Diff

View File

@ -1,50 +0,0 @@
LAMMPS data file via write_data, version 22 Nov 2016, timestep = 0
32 atoms
3 atom types
0.0000000000000000e+00 5.0000007843669128e+00 xlo xhi
0.0000000000000000e+00 5.1961532151380023e+00 ylo yhi
0.0000000000000000e+00 1.0000000000000000e+01 zlo zhi
Masses
1 1
2 1
3 1
Atoms # dielectric: id mol type q x y z normx normy normz area_per_patch ed em epsilon curvature
1 1 1 0 0 0 5 0 0 1 1 8 6 1 0
2 1 1 0 0.5 0.866026 5 0 0 1 1 8 6 1 0
3 1 1 0 1 0 5 0 0 1 1 8 6 1 0
4 1 1 0 1.5 0.866026 5 0 0 1 1 8 6 1 0
5 1 1 0 2 0 5 0 0 1 1 8 6 1 0
6 1 1 0 2.5 0.866026 5 0 0 1 1 8 6 1 0
7 1 1 0 3 0 5 0 0 1 1 8 6 1 0
8 1 1 0 3.5 0.866026 5 0 0 1 1 8 6 1 0
9 1 1 0 4 0 5 0 0 1 1 8 6 1 0
10 1 1 0 4.5 0.866026 5 0 0 1 1 8 6 1 0
11 1 1 0 0 1.73205 5 0 0 1 1 8 6 1 0
12 1 1 0 0.5 2.59808 5 0 0 1 1 8 6 1 0
13 1 1 0 1 1.73205 5 0 0 1 1 8 6 1 0
14 1 1 0 1.5 2.59808 5 0 0 1 1 8 6 1 0
15 1 1 0 2 1.73205 5 0 0 1 1 8 6 1 0
16 1 1 0 2.5 2.59808 5 0 0 1 1 8 6 1 0
17 1 1 0 3 1.73205 5 0 0 1 1 8 6 1 0
18 1 1 0 3.5 2.59808 5 0 0 1 1 8 6 1 0
19 1 1 0 4 1.73205 5 0 0 1 1 8 6 1 0
20 1 1 0 4.5 2.59808 5 0 0 1 1 8 6 1 0
21 1 1 0 0 3.4641 5 0 0 1 1 8 6 1 0
22 1 1 0 0.5 4.33013 5 0 0 1 1 8 6 1 0
23 1 1 0 1 3.4641 5 0 0 1 1 8 6 1 0
24 1 1 0 1.5 4.33013 5 0 0 1 1 8 6 1 0
25 1 1 0 2 3.4641 5 0 0 1 1 8 6 1 0
26 1 1 0 2.5 4.33013 5 0 0 1 1 8 6 1 0
27 1 1 0 3 3.4641 5 0 0 1 1 8 6 1 0
28 1 1 0 3.5 4.33013 5 0 0 1 1 8 6 1 0
29 1 1 0 4 3.4641 5 0 0 1 1 8 6 1 0
30 1 1 0 4.5 4.33013 5 0 0 1 1 8 6 1 0
31 2 2 1 1 1 7 0 0 1 1 8 6 10 0
32 2 3 -1 3 3 7 0 0 1 1 8 6 10 0

File diff suppressed because it is too large Load Diff

View File

@ -14,11 +14,11 @@ units lj
atom_style dielectric
atom_modify map array
dimension 3
boundary f f f
boundary p p f
variable method index gmres # gmres = BEM/GMRES
# icc = BEM/ICC*
# dof = Variational
# dof = Direct optimization of the functional
# none
variable ed equal "v_epsilon2 - v_epsilon1"
@ -41,10 +41,13 @@ set group cations charge ${qscale}
variable qscale equal "-1.0 / v_epsilon2"
set group anions charge ${qscale}
pair_style lj/cut/coul/cut/dielectric 1.122 20.0
pair_style lj/cut/coul/long/dielectric 1.122 10.0
pair_coeff * * 1.0 1.0
pair_coeff 1 1 0.0 1.0
kspace_style pppm/dielectric 0.0001
kspace_modify slab 3.0
neigh_modify every 1 delay 0 check yes one 5000
#compute ef all efield/atom

View File

@ -1,39 +0,0 @@
# Interface
newton off
units lj
atom_style dielectric
atom_modify map array
dimension 3
boundary p p f
read_data data.minimal
group interface type 1
group ions type 2 3
pair_style coul/long/dielectric 5.0
#pair_style coul/cut/dielectric 20.0
pair_coeff * *
thermo 1000
thermo_style custom step evdwl ecoul elong epair
thermo_modify flush yes
#compute ef all efield/atom
dump 1 all custom 100 dump.txt id type q x y z #fx fy fz c_ef[1] c_ef[2] c_ef[3]
dump_modify 1 sort id
fix 1 ions nve
fix 3 interface polarize/bem/icc 1 1.0e-4
kspace_style pppm/dielectric 0.0001
kspace_modify slab 3.0
velocity all create 1.0 2947942 mom yes dist gaussian
run 0

View File

@ -8,25 +8,25 @@ boundary f f f
variable method index gmres # gmres = BEM/GMRES
# icc = BEM/ICC*
# dof = Variational
# dof = Direct optimization of the functional
# none
variable data index data.sphere #data.plane # data.minimal
variable data index data.sphere
read_data ${data}
group interface type 1
group ions type 2 3
pair_style lj/cut/coul/cut/dielectric 50.0
pair_style lj/cut/coul/cut/dielectric 1.122 20.0
pair_coeff * * 1.0 1.0
pair_coeff 1 1 0.0 1.0
neigh_modify one 5000
#compute ef all efield/atom
dump 1 all custom 100 all.dump id mol type q x y z #fx fy fz c_ef[1] c_ef[2] c_ef[3]
dump 2 interface custom 100 interface.dump id mol type q x y z #fx fy fz c_ef[1] c_ef[2] c_ef[3]
dump 1 all custom 100 all.dump id mol type q x y z #c_ef[1] c_ef[2] c_ef[3]
dump 2 interface custom 100 interface.dump id mol type q x y z #c_ef[1] c_ef[2] c_ef[3]
dump_modify 1 sort id
@ -42,7 +42,7 @@ else &
"print 'Unsupported method for polarization' "
thermo 1000
thermo_style custom step evdwl ecoul elong epair #f_3
thermo_style custom step evdwl ecoul elong epair
thermo_modify flush yes
run 0

View File

@ -1,54 +0,0 @@
# Interface
newton off
units lj
atom_style dielectric
atom_modify map array
dimension 3
boundary p p f
variable method index gmres # gmres = BEM/GMRES
# icc = BEM/ICC*
# dof = Variational
# none
variable data index data.sphere #data.plane # data.minimal
read_data ${data}
group interface type 1
group ions type 2 3
pair_style lj/cut/coul/long/dielectric 1.122 10.0
pair_coeff * * 1.0 1.0
pair_coeff 1 1 0.0 1.0
kspace_style pppm/dielectric 0.0001
kspace_modify slab 3.0
neigh_modify every 1 delay 0 check yes one 5000
#compute ef all efield/atom
dump 1 all custom 100 all.dump id mol type q x y z #fx fy fz c_ef[1] c_ef[2] c_ef[3]
dump 2 interface custom 100 interface.dump id mol type q x y z #fx fy fz c_ef[1] c_ef[2] c_ef[3]
dump_modify 1 sort id
fix 1 ions nve
if "${method} == gmres" then &
"fix 3 interface polarize/bem/gmres 1 1.0e-4" &
elif "${method} == icc"&
"fix 3 interface polarize/bem/icc 1 1.0e-4 itr_max 50" &
elif "${method} == dof" &
"fix 3 interface polarize/functional 1 0.0001" &
else &
"print 'Unsupported method for polarization' "
thermo 100
thermo_style custom step evdwl ecoul elong epair cpu
thermo_modify flush yes
run 0