Merge pull request #3787 from oywg11/ilp-water-graphene

Registry-Dependent Potential for Interfaces of Water with Graphene
This commit is contained in:
Axel Kohlmeyer
2023-06-07 11:48:43 -04:00
committed by GitHub
33 changed files with 3867 additions and 71 deletions

View File

@ -37,6 +37,7 @@ OPT.
*
* :doc:`adp (ko) <pair_adp>`
* :doc:`agni (o) <pair_agni>`
* :doc:`aip/water/2dm (t) <pair_aip_water_2dm>`
* :doc:`airebo (io) <pair_airebo>`
* :doc:`airebo/morse (io) <pair_airebo>`
* :doc:`amoeba (g) <pair_amoeba>`

View File

@ -0,0 +1,167 @@
.. index:: pair_style aip/water/2dm
.. index:: pair_style aip/water/2dm/opt
pair_style aip/water/2dm command
===================================
Accelerator Variant: *aip/water/2dm/opt*
Syntax
""""""
.. code-block:: LAMMPS
pair_style [hybrid/overlay ...] aip/water/2dm cutoff tap_flag
* cutoff = global cutoff (distance units)
* tap_flag = 0/1 to turn off/on the taper function
Examples
""""""""
.. code-block:: LAMMPS
pair_style hybrid/overlay aip/water/2dm 16.0 1
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw
pair_style hybrid/overlay aip/water/2dm 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw
Description
"""""""""""
.. versionadded:: TBD
The *aip/water/2dm* style computes the anisotropic interfacial potential
(AIP) potential for interfaces of water with two-dimensional (2D)
materials as described in :ref:`(Feng) <Feng>`.
.. math::
E = & \frac{1}{2} \sum_i \sum_{j \neq i} V_{ij} \\
V_{ij} = & {\rm Tap}(r_{ij})\left \{ e^{-\alpha (r_{ij}/\beta -1)}
\left [ \epsilon + f(\rho_{ij}) + f(\rho_{ji})\right ] -
\frac{1}{1+e^{-d\left [ \left ( r_{ij}/\left (s_R \cdot r^{eff} \right ) \right )-1 \right ]}}
\cdot \frac{C_6}{r^6_{ij}} \right \}\\
\rho_{ij}^2 = & r_{ij}^2 - ({\bf r}_{ij} \cdot {\bf n}_i)^2 \\
\rho_{ji}^2 = & r_{ij}^2 - ({\bf r}_{ij} \cdot {\bf n}_j)^2 \\
f(\rho) = & C e^{ -( \rho / \delta )^2 } \\
{\rm Tap}(r_{ij}) = & 20\left ( \frac{r_{ij}}{R_{cut}} \right )^7 -
70\left ( \frac{r_{ij}}{R_{cut}} \right )^6 +
84\left ( \frac{r_{ij}}{R_{cut}} \right )^5 -
35\left ( \frac{r_{ij}}{R_{cut}} \right )^4 + 1
Where :math:`\mathrm{Tap}(r_{ij})` is the taper function which provides
a continuous cutoff (up to third derivative) for interatomic separations
larger than :math:`r_c` :doc:`pair_style ilp_graphene_hbn
<pair_ilp_graphene_hbn>`.
.. note::
This pair style uses the atomic normal vector definition from
:ref:`(Feng) <Feng>`), where the atomic normal vectors of the
hydrogen atoms are assumed to lie along the corresponding
oxygen-hydrogen bonds and the normal vector of the central oxygen
atom is defined as their average.
The provided parameter file, ``COH.aip.water.2dm``, is intended for use
with *metal* :doc:`units <units>`, with energies in meV. Two additional
parameters, *S*, and *rcut* are included in the parameter file. *S* is
designed to facilitate scaling of energies; *rcut* is the cutoff for an
internal, short distance neighbor list that is generated for speeding up
the calculation of the normals for all atom pairs.
.. note::
The parameters presented in the provided parameter file,
``COH.aip.water.2dm``, are fitted with the taper function enabled by
setting the cutoff equal to 16.0 Angstrom. Using a different cutoff
or taper function setting should be carefully checked as they can
lead to significant errors. These parameters provide a good
description in both short- and long-range interaction regimes. This
is essential for simulations in high pressure regime (i.e., the
interlayer distance is smaller than the equilibrium distance).
This potential must be used in combination with hybrid/overlay. Other
interactions can be set to zero using :doc:`pair_coeff settings
<pair_coeff>` with the pair style set to *none*\ .
This pair style tallies a breakdown of the total interlayer potential
energy into sub-categories, which can be accessed via the :doc:`compute
pair <compute_pair>` command as a vector of values of length 2. The 2
values correspond to the following sub-categories:
1. *E_vdW* = vdW (attractive) energy
2. *E_Rep* = Repulsive energy
To print these quantities to the log file (with descriptive column
headings) the following commands could be included in an input script:
.. code-block:: LAMMPS
compute 0 all pair aip/water/2dm
variable Evdw equal c_0[1]
variable Erep equal c_0[2]
thermo_style custom step temp epair v_Erep v_Evdw
----------
.. include:: accel_styles.rst
----------
Mixing, shift, table, tail correction, restart, rRESPA info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
This pair style does not support the pair_modify mix, shift, table, and
tail options.
This pair style does not write their information to binary restart
files, since it is stored in potential files. Thus, you need to
re-specify the pair_style and pair_coeff commands in an input script
that reads a restart file.
Restrictions
""""""""""""
This pair style is part of the INTERLAYER package. It is only enabled
if LAMMPS was built with that package. See the :doc:`Build package
<Build_package>` page for more info.
This pair style requires the newton setting to be *on* for pair
interactions.
The ``COH.aip.water.2dm`` potential file provided with LAMMPS is
parameterized for *metal* units. You can use this pair style with any
LAMMPS units, but you would need to create your own potential file with
parameters in the appropriate units, if your simulation does not use
*metal* units.
Related commands
""""""""""""""""
:doc:`pair_coeff <pair_coeff>`,
:doc:`pair_none <pair_none>`,
:doc:`pair_style hybrid/overlay <pair_hybrid>`,
:doc:`pair_style drip <pair_drip>`,
:doc:`pair_style ilp_tmd <pair_ilp_tmd>`,
:doc:`pair_style saip_metal <pair_saip_metal>`,
:doc:`pair_style ilp_graphene_hbn <pair_ilp_graphene_hbn>`,
:doc:`pair_style pair_kolmogorov_crespi_z <pair_kolmogorov_crespi_z>`,
:doc:`pair_style pair_kolmogorov_crespi_full <pair_kolmogorov_crespi_full>`,
:doc:`pair_style pair_lebedeva_z <pair_lebedeva_z>`,
:doc:`pair_style pair_coul_shield <pair_coul_shield>`.
Default
"""""""
tap_flag = 1
----------
.. _Feng:
**(Feng)** Z. Feng and W. Ouyang et al., J. Phys. Chem. C. 127, 8704-8713 (2023).

View File

@ -155,8 +155,8 @@ interactions.
The BNCH.ILP potential file provided with LAMMPS (see the potentials
directory) are parameterized for *metal* units. You can use this
potential with any LAMMPS units, but you would need to create your
BNCH.ILP potential file with coefficients listed in the appropriate
potential with any LAMMPS units, but you would need to create your own
custom BNCH.ILP potential file with coefficients listed in the appropriate
units, if your simulation does not use *metal* units.
Related commands
@ -181,6 +181,14 @@ tap_flag = 1
----------
.. _Ouyang1:
**(Ouyang1)** W. Ouyang, D. Mandelli, M. Urbakh and O. Hod, Nano Lett. 18, 6009-6016 (2018).
.. _Ouyang2:
**(Ouyang2)** W. Ouyang et al., J. Chem. Theory Comput. 16(1), 666-676 (2020).
.. _Leven1:
**(Leven1)** I. Leven, I. Azuri, L. Kronik and O. Hod, J. Chem. Phys. 140, 104106 (2014).
@ -196,11 +204,3 @@ tap_flag = 1
.. _Kolmogorov2:
**(Kolmogorov)** A. N. Kolmogorov, V. H. Crespi, Phys. Rev. B 71, 235415 (2005).
.. _Ouyang1:
**(Ouyang1)** W. Ouyang, D. Mandelli, M. Urbakh and O. Hod, Nano Lett. 18, 6009-6016 (2018).
.. _Ouyang2:
**(Ouyang2)** W. Ouyang et al., J. Chem. Theory Comput. 16(1), 666-676 (2020).

View File

@ -135,8 +135,8 @@ interactions.
The TMD.ILP potential file provided with LAMMPS (see the potentials
directory) are parameterized for *metal* units. You can use this
potential with any LAMMPS units, but you would need to create your
BNCH.ILP potential file with coefficients listed in the appropriate
potential with any LAMMPS units, but you would need to create your own
custom TMD.ILP potential file with coefficients listed in the appropriate
units, if your simulation does not use *metal* units.
Related commands

View File

@ -129,7 +129,7 @@ interactions.
The CH.KC potential file provided with LAMMPS (see the potentials
folder) is parameterized for metal units. You can use this pair style
with any LAMMPS units, but you would need to create your own custom
CC.KC potential file with all coefficients converted to the appropriate
CH.KC potential file with all coefficients converted to the appropriate
units.
Related commands

View File

@ -134,8 +134,8 @@ interactions.
The CHAu.ILP potential file provided with LAMMPS (see the potentials
directory) are parameterized for *metal* units. You can use this
potential with any LAMMPS units, but you would need to create your
BNCH.ILP potential file with coefficients listed in the appropriate
potential with any LAMMPS units, but you would need to create your own
custom CHAu.ILP potential file with coefficients listed in the appropriate
units, if your simulation does not use *metal* units.
Related commands

View File

@ -114,6 +114,7 @@ accelerated styles exist.
* :doc:`adp <pair_adp>` - angular dependent potential (ADP) of Mishin
* :doc:`agni <pair_agni>` - AGNI machine-learning potential
* :doc:`aip/water/2dm <pair_aip_water_2dm>` - anisotropic interfacial potential for water in 2d geometries
* :doc:`airebo <pair_airebo>` - AIREBO potential of Stuart
* :doc:`airebo/morse <pair_airebo>` - AIREBO with Morse instead of LJ
* :doc:`amoeba <pair_amoeba>` -

View File

@ -55,6 +55,7 @@ Ai
Aidan
aij
aimd
aip
airebo
Aj
ajaramil
@ -528,6 +529,7 @@ collisional
Columic
colvars
Colvars
COH
COLVARS
comID
Commun
@ -1090,6 +1092,7 @@ Fellinger
femtosecond
femtoseconds
fene
Feng
Fennell
fep
FEP

View File

@ -0,0 +1 @@
../../../../potentials/COH.aip.water.2dm

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,51 @@
# Initialization
units metal
boundary p p p
atom_style full
processors * * 1 # domain decomposition over x and y
read_data ./gra_water.data
mass 1 12.0107 # carbon mass (g/mole)
mass 2 15.9994 # oxygen mass (g/mole)
mass 3 1.008 # hydrogen mass (g/mole)
# Separate atom groups
group gr molecule 1
group water molecule 2
######################## Potential defition ##############################
# Interlayer potential
pair_style hybrid/overlay aip/water/2dm 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
####################################################################
pair_coeff 1 1 none
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw # C-H2O
# bond and angle
bond_style harmonic
bond_coeff 1 0.0 0.9572
angle_style harmonic
angle_coeff 1 0.0 104.52
# define kspace calculation
kspace_style pppm/tip4p 1E-5
# Neighbor update settings
neighbor 2.0 bin
neigh_modify every 1 delay 5 check yes page 1000000 one 100000
####################################################################
# Calculate pair energy
compute 1 all pair lj/cut/tip4p/long
compute 2 all pair aip/water/2dm
compute wt water temp
variable TIP4P equal c_1
variable EILP equal c_2 # total interlayer energy
variable temp_wt equal c_wt
############# Output ##############
thermo_style custom step etotal pe ke v_TIP4P v_EILP v_temp_wt
thermo 100
thermo_modify lost error
fix subf gr setforce 0.0 0.0 0.0
fix 1 water shake 0.0001 20 100 b 1 a 1
timestep 1e-3
velocity water create 300.0 12345 dist gaussian mom yes rot yes
fix 2 water nve
run 1000

View File

@ -0,0 +1,51 @@
# Initialization
units metal
boundary p p p
atom_style full
processors * * 1 # domain decomposition over x and y
read_data ./gra_water.data
mass 1 12.0107 # carbon mass (g/mole)
mass 2 15.9994 # oxygen mass (g/mole)
mass 3 1.008 # hydrogen mass (g/mole)
# Separate atom groups
group gr molecule 1
group water molecule 2
######################## Potential defition ##############################
# Interlayer potential
pair_style hybrid/overlay aip/water/2dm/opt 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
####################################################################
pair_coeff 1 1 none
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm/opt COH.aip.water.2dm C Ow Hw # C-H2O
# bond and angle
bond_style harmonic
bond_coeff 1 0.0 0.9572
angle_style harmonic
angle_coeff 1 0.0 104.52
# define kspace calculation
kspace_style pppm/tip4p 1E-5
# Neighbor update settings
neighbor 2.0 bin
neigh_modify every 1 delay 5 check yes page 1000000 one 100000
####################################################################
# Calculate pair energy
compute 1 all pair lj/cut/tip4p/long
compute 2 all pair aip/water/2dm/opt
compute wt water temp
variable TIP4P equal c_1
variable EILP equal c_2 # total interlayer energy
variable temp_wt equal c_wt
############# Output ##############
thermo_style custom step etotal pe ke v_TIP4P v_EILP v_temp_wt
thermo 100
thermo_modify lost error
fix subf gr setforce 0.0 0.0 0.0
fix 1 water shake 0.0001 20 100 b 1 a 1
timestep 1e-3
velocity water create 300.0 12345 dist gaussian mom yes rot yes
fix 2 water nve
run 1000

View File

@ -0,0 +1,239 @@
LAMMPS (23 Jun 2022 - Update 4)
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:98)
using 1 OpenMP thread(s) per MPI task
# Initialization
units metal
boundary p p p
atom_style full
processors * * 1 # domain decomposition over x and y
read_data ./gra_water.data
Reading data file ...
orthogonal box = (0 0 0) to (46.92336 44.331078 200)
1 by 1 by 1 MPI processor grid
reading atoms ...
936 atoms
reading velocities ...
936 velocities
scanning bonds ...
2 = max bonds/atom
scanning angles ...
1 = max angles/atom
reading bonds ...
96 bonds
reading angles ...
48 angles
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 0 0
special bond factors coul: 0 0 0
2 = max # of 1-2 neighbors
1 = max # of 1-3 neighbors
1 = max # of 1-4 neighbors
2 = max # of special neighbors
special bonds CPU = 0.000 seconds
read_data CPU = 0.012 seconds
mass 1 12.0107 # carbon mass (g/mole)
mass 2 15.9994 # oxygen mass (g/mole)
mass 3 1.008 # hydrogen mass (g/mole)
# Separate atom groups
group gr molecule 1
792 atoms in group gr
group water molecule 2
144 atoms in group water
######################## Potential defition ##############################
# Interlayer potential
pair_style hybrid/overlay aip/water/2dm 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
####################################################################
pair_coeff 1 1 none
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw # C-H2O
Reading aip/water/2dm potential file COH.aip.water.2dm with DATE: 2022-12-02
# bond and angle
bond_style harmonic
bond_coeff 1 0.0 0.9572
angle_style harmonic
angle_coeff 1 0.0 104.52
# define kspace calculation
kspace_style pppm/tip4p 1E-5
# Neighbor update settings
neighbor 2.0 bin
neigh_modify every 1 delay 5 check yes page 1000000 one 100000
####################################################################
# Calculate pair energy
compute 1 all pair lj/cut/tip4p/long
compute 2 all pair aip/water/2dm
compute wt water temp
variable TIP4P equal c_1
variable EILP equal c_2 # total interlayer energy
variable temp_wt equal c_wt
############# Output ##############
thermo_style custom step etotal pe ke v_TIP4P v_EILP v_temp_wt
thermo 100
thermo_modify lost error
fix subf gr setforce 0.0 0.0 0.0
fix 1 water shake 0.0001 20 100 b 1 a 1
0 = # of size 2 clusters
0 = # of size 3 clusters
0 = # of size 4 clusters
48 = # of frozen angles
find clusters CPU = 0.000 seconds
timestep 1e-3
velocity water create 300.0 12345 dist gaussian mom yes rot yes
fix 2 water nve
run 1000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- ilp/graphene/hbn potential doi:10.1021/acs.nanolett.8b02848
@Article{Ouyang2018
author = {W. Ouyang and D. Mandelli and M. Urbakh and O. Hod},
title = {Nanoserpents: Graphene Nanoribbon Motion on Two-Dimensional Hexagonal Materials},
journal = {Nano Letters},
volume = 18,
pages = 6009,
year = 2018,
}
- ilp/tmd potential doi:10.1021/acs.jctc.1c00782
@Article{Ouyang2021
author = {W. Ouyang and R. Sofer and X. Gao and J. Hermann and
A. Tkatchenko and L. Kronik and M. Urbakh and O. Hod},
title = {Anisotropic Interlayer Force Field for Transition
Metal Dichalcogenides: The Case of Molybdenum Disulfide},
journal = {J.~Chem.\ Theory Comput.},
volume = 17,
pages = {7237--7245}
year = 2021,
}
- ilp/water/2dm potential doi/10.1021/acs.jpcc.2c08464
@Article{Feng2023
author = {Z. Feng, Y. Yao, J. Liu, B. Wu, Z. Liu, and W. Ouyang},
title = {Registry-Dependent Potential for Interfaces of Water with Graphene},
journal = {J. Phys. Chem. C},
volume = 127,
pages = {8704-8713}
year = 2023,
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
PPPM initialization ...
extracting TIP4P info from pair style
using 12-bit tables for long-range coulomb (../kspace.cpp:342)
G vector (1/distance) = 0.28684806
grid = 25 24 80
stencil order = 5
estimated absolute RMS force accuracy = 0.0001640931
estimated relative force accuracy = 1.1395635e-05
using single precision MKL FFT
3d grid and FFT values/proc = 84320 48000
WARNING: Using a manybody potential with bonds/angles/dihedrals and special_bond exclusions (../pair.cpp:239)
WARNING: Communication cutoff 0 is shorter than a bond length based estimate of 3.4358. This may lead to errors. (../comm.cpp:727)
WARNING: Increasing communication cutoff to 11.6118 for TIP4P pair style (../pair_lj_cut_tip4p_long.cpp:484)
Neighbor list info ...
update every 1 steps, delay 5 steps, check yes
max neighbors/atom: 100000, page size: 1000000
master list distance cutoff = 18
ghost atom cutoff = 18
binsize = 9, bins = 6 5 23
3 neighbor lists, perpetual/occasional/extra = 3 0 0
(1) pair aip/water/2dm, perpetual
attributes: full, newton on, ghost
pair build: full/bin/ghost
stencil: full/ghost/bin/3d
bin: standard
(2) pair lj/cut/tip4p/long, perpetual, skip from (3)
attributes: half, newton on
pair build: skip
stencil: none
bin: none
(3) neighbor class addition, perpetual
attributes: half, newton on
pair build: half/bin/newton
stencil: half/bin/3d
bin: standard
WARNING: Communication cutoff adjusted to 18 (../comm.cpp:736)
SHAKE stats (type/ave/delta/count) on step 0
Bond: 1 0.957201 2.19705e-06 96
Angle: 1 104.52 0.000203056 48
Per MPI rank memory allocation (min/avg/max) = 32.9 | 32.9 | 32.9 Mbytes
Step TotEng PotEng KinEng v_TIP4P v_EILP v_temp_wt
0 -16.131263 -19.815178 3.6839141 189.37246 -1.903543 300
SHAKE stats (type/ave/delta/count) on step 100
Bond: 1 0.9572 9.54949e-07 96
Angle: 1 104.52 6.01522e-05 48
100 -17.494868 -20.796993 3.3021253 188.4955 -1.8981262 268.90898
SHAKE stats (type/ave/delta/count) on step 200
Bond: 1 0.9572 9.63922e-07 96
Angle: 1 104.52 7.7021e-05 48
200 -17.486271 -21.194892 3.7086213 188.14561 -1.9871708 302.01203
SHAKE stats (type/ave/delta/count) on step 300
Bond: 1 0.9572 1.4264e-06 96
Angle: 1 104.52 6.48393e-05 48
300 -17.502844 -20.993704 3.49086 188.23268 -1.8457229 284.27861
SHAKE stats (type/ave/delta/count) on step 400
Bond: 1 0.9572 1.33728e-06 96
Angle: 1 104.52 7.6239e-05 48
400 -17.495287 -20.828353 3.3330658 188.48002 -1.8429075 271.42862
SHAKE stats (type/ave/delta/count) on step 500
Bond: 1 0.9572 1.14685e-06 96
Angle: 1 104.52 8.58621e-05 48
500 -17.491435 -20.443044 2.9516084 188.7589 -1.8566335 240.36459
SHAKE stats (type/ave/delta/count) on step 600
Bond: 1 0.9572 9.17601e-07 96
Angle: 1 104.52 8.24516e-05 48
600 -17.505684 -20.608457 3.1027731 188.72078 -1.9560796 252.67471
SHAKE stats (type/ave/delta/count) on step 700
Bond: 1 0.9572 9.50422e-07 96
Angle: 1 104.52 5.62423e-05 48
700 -17.496703 -21.072663 3.5759596 188.2777 -1.9833956 291.20871
SHAKE stats (type/ave/delta/count) on step 800
Bond: 1 0.9572 1.15262e-06 96
Angle: 1 104.52 7.02157e-05 48
800 -17.478623 -20.819504 3.3408809 188.37868 -1.9112996 272.06505
SHAKE stats (type/ave/delta/count) on step 900
Bond: 1 0.9572 9.14138e-07 96
Angle: 1 104.52 6.98742e-05 48
900 -17.48086 -20.728495 3.2476349 188.59022 -1.8922102 264.47155
SHAKE stats (type/ave/delta/count) on step 1000
Bond: 1 0.9572 1.00586e-06 96
Angle: 1 104.52 0.000111712 48
1000 -17.498465 -20.331545 2.8330804 188.87473 -1.812177 230.71225
Loop time of 20.3334 on 1 procs for 1000 steps with 936 atoms
Performance: 4.249 ns/day, 5.648 hours/ns, 49.180 timesteps/s
99.5% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 16.179 | 16.179 | 16.179 | 0.0 | 79.57
Bond | 0.00021103 | 0.00021103 | 0.00021103 | 0.0 | 0.00
Kspace | 3.3118 | 3.3118 | 3.3118 | 0.0 | 16.29
Neigh | 0.79017 | 0.79017 | 0.79017 | 0.0 | 3.89
Comm | 0.026379 | 0.026379 | 0.026379 | 0.0 | 0.13
Output | 0.00046496 | 0.00046496 | 0.00046496 | 0.0 | 0.00
Modify | 0.017013 | 0.017013 | 0.017013 | 0.0 | 0.08
Other | | 0.008835 | | | 0.04
Nlocal: 936 ave 936 max 936 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 5242 ave 5242 max 5242 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 0 ave 0 max 0 min
Histogram: 1 0 0 0 0 0 0 0 0 0
FullNghs: 431382 ave 431382 max 431382 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 431382
Ave neighs/atom = 460.87821
Ave special neighs/atom = 0.30769231
Neighbor list builds = 28
Dangerous builds = 0
Total wall time: 0:00:20

View File

@ -0,0 +1,239 @@
LAMMPS (23 Jun 2022 - Update 4)
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:98)
using 1 OpenMP thread(s) per MPI task
# Initialization
units metal
boundary p p p
atom_style full
processors * * 1 # domain decomposition over x and y
read_data ./gra_water.data
Reading data file ...
orthogonal box = (0 0 0) to (46.92336 44.331078 200)
2 by 2 by 1 MPI processor grid
reading atoms ...
936 atoms
reading velocities ...
936 velocities
scanning bonds ...
2 = max bonds/atom
scanning angles ...
1 = max angles/atom
reading bonds ...
96 bonds
reading angles ...
48 angles
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 0 0
special bond factors coul: 0 0 0
2 = max # of 1-2 neighbors
1 = max # of 1-3 neighbors
1 = max # of 1-4 neighbors
2 = max # of special neighbors
special bonds CPU = 0.001 seconds
read_data CPU = 0.017 seconds
mass 1 12.0107 # carbon mass (g/mole)
mass 2 15.9994 # oxygen mass (g/mole)
mass 3 1.008 # hydrogen mass (g/mole)
# Separate atom groups
group gr molecule 1
792 atoms in group gr
group water molecule 2
144 atoms in group water
######################## Potential defition ##############################
# Interlayer potential
pair_style hybrid/overlay aip/water/2dm 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
####################################################################
pair_coeff 1 1 none
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm COH.aip.water.2dm C Ow Hw # C-H2O
Reading aip/water/2dm potential file COH.aip.water.2dm with DATE: 2022-12-02
# bond and angle
bond_style harmonic
bond_coeff 1 0.0 0.9572
angle_style harmonic
angle_coeff 1 0.0 104.52
# define kspace calculation
kspace_style pppm/tip4p 1E-5
# Neighbor update settings
neighbor 2.0 bin
neigh_modify every 1 delay 5 check yes page 1000000 one 100000
####################################################################
# Calculate pair energy
compute 1 all pair lj/cut/tip4p/long
compute 2 all pair aip/water/2dm
compute wt water temp
variable TIP4P equal c_1
variable EILP equal c_2 # total interlayer energy
variable temp_wt equal c_wt
############# Output ##############
thermo_style custom step etotal pe ke v_TIP4P v_EILP v_temp_wt
thermo 100
thermo_modify lost error
fix subf gr setforce 0.0 0.0 0.0
fix 1 water shake 0.0001 20 100 b 1 a 1
0 = # of size 2 clusters
0 = # of size 3 clusters
0 = # of size 4 clusters
48 = # of frozen angles
find clusters CPU = 0.000 seconds
timestep 1e-3
velocity water create 300.0 12345 dist gaussian mom yes rot yes
fix 2 water nve
run 1000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- ilp/graphene/hbn potential doi:10.1021/acs.nanolett.8b02848
@Article{Ouyang2018
author = {W. Ouyang and D. Mandelli and M. Urbakh and O. Hod},
title = {Nanoserpents: Graphene Nanoribbon Motion on Two-Dimensional Hexagonal Materials},
journal = {Nano Letters},
volume = 18,
pages = 6009,
year = 2018,
}
- ilp/tmd potential doi:10.1021/acs.jctc.1c00782
@Article{Ouyang2021
author = {W. Ouyang and R. Sofer and X. Gao and J. Hermann and
A. Tkatchenko and L. Kronik and M. Urbakh and O. Hod},
title = {Anisotropic Interlayer Force Field for Transition
Metal Dichalcogenides: The Case of Molybdenum Disulfide},
journal = {J.~Chem.\ Theory Comput.},
volume = 17,
pages = {7237--7245}
year = 2021,
}
- ilp/water/2dm potential doi/10.1021/acs.jpcc.2c08464
@Article{Feng2023
author = {Z. Feng, Y. Yao, J. Liu, B. Wu, Z. Liu, and W. Ouyang},
title = {Registry-Dependent Potential for Interfaces of Water with Graphene},
journal = {J. Phys. Chem. C},
volume = 127,
pages = {8704-8713}
year = 2023,
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
PPPM initialization ...
extracting TIP4P info from pair style
using 12-bit tables for long-range coulomb (../kspace.cpp:342)
G vector (1/distance) = 0.28684806
grid = 25 24 80
stencil order = 5
estimated absolute RMS force accuracy = 0.0001640931
estimated relative force accuracy = 1.1395635e-05
using single precision MKL FFT
3d grid and FFT values/proc = 30685 12480
WARNING: Using a manybody potential with bonds/angles/dihedrals and special_bond exclusions (../pair.cpp:239)
WARNING: Communication cutoff 0 is shorter than a bond length based estimate of 3.4358. This may lead to errors. (../comm.cpp:727)
WARNING: Increasing communication cutoff to 11.6118 for TIP4P pair style (../pair_lj_cut_tip4p_long.cpp:484)
Neighbor list info ...
update every 1 steps, delay 5 steps, check yes
max neighbors/atom: 100000, page size: 1000000
master list distance cutoff = 18
ghost atom cutoff = 18
binsize = 9, bins = 6 5 23
3 neighbor lists, perpetual/occasional/extra = 3 0 0
(1) pair aip/water/2dm, perpetual
attributes: full, newton on, ghost
pair build: full/bin/ghost
stencil: full/ghost/bin/3d
bin: standard
(2) pair lj/cut/tip4p/long, perpetual, skip from (3)
attributes: half, newton on
pair build: skip
stencil: none
bin: none
(3) neighbor class addition, perpetual
attributes: half, newton on
pair build: half/bin/newton
stencil: half/bin/3d
bin: standard
WARNING: Communication cutoff adjusted to 18 (../comm.cpp:736)
SHAKE stats (type/ave/delta/count) on step 0
Bond: 1 0.957201 2.19705e-06 96
Angle: 1 104.52 0.000203056 48
Per MPI rank memory allocation (min/avg/max) = 25.22 | 25.25 | 25.29 Mbytes
Step TotEng PotEng KinEng v_TIP4P v_EILP v_temp_wt
0 -16.131263 -19.815178 3.6839141 189.37246 -1.903543 300
SHAKE stats (type/ave/delta/count) on step 100
Bond: 1 0.9572 9.54949e-07 96
Angle: 1 104.52 6.01522e-05 48
100 -17.494869 -20.796995 3.3021253 188.4955 -1.8981262 268.90898
SHAKE stats (type/ave/delta/count) on step 200
Bond: 1 0.9572 9.63922e-07 96
Angle: 1 104.52 7.7021e-05 48
200 -17.48627 -21.194892 3.7086213 188.14561 -1.9871708 302.01203
SHAKE stats (type/ave/delta/count) on step 300
Bond: 1 0.9572 1.4264e-06 96
Angle: 1 104.52 6.48393e-05 48
300 -17.502843 -20.993703 3.4908599 188.23268 -1.8457229 284.27861
SHAKE stats (type/ave/delta/count) on step 400
Bond: 1 0.9572 1.33728e-06 96
Angle: 1 104.52 7.6239e-05 48
400 -17.495285 -20.82835 3.333065 188.48003 -1.8429074 271.42856
SHAKE stats (type/ave/delta/count) on step 500
Bond: 1 0.9572 1.14685e-06 96
Angle: 1 104.52 8.58621e-05 48
500 -17.491436 -20.443043 2.9516075 188.7589 -1.8566335 240.36452
SHAKE stats (type/ave/delta/count) on step 600
Bond: 1 0.9572 9.17601e-07 96
Angle: 1 104.52 8.24517e-05 48
600 -17.505683 -20.608456 3.1027734 188.72078 -1.9560795 252.67474
SHAKE stats (type/ave/delta/count) on step 700
Bond: 1 0.9572 9.50425e-07 96
Angle: 1 104.52 5.62422e-05 48
700 -17.496706 -21.072664 3.575958 188.2777 -1.9833951 291.20858
SHAKE stats (type/ave/delta/count) on step 800
Bond: 1 0.9572 1.15256e-06 96
Angle: 1 104.52 7.02177e-05 48
800 -17.478628 -20.819507 3.340879 188.37868 -1.9113009 272.06489
SHAKE stats (type/ave/delta/count) on step 900
Bond: 1 0.9572 9.14163e-07 96
Angle: 1 104.52 6.98849e-05 48
900 -17.480865 -20.728504 3.2476386 188.5902 -1.8922108 264.47185
SHAKE stats (type/ave/delta/count) on step 1000
Bond: 1 0.9572 1.00568e-06 96
Angle: 1 104.52 0.000111707 48
1000 -17.498474 -20.331607 2.833133 188.87466 -1.8121689 230.71654
Loop time of 9.05361 on 4 procs for 1000 steps with 936 atoms
Performance: 9.543 ns/day, 2.515 hours/ns, 110.453 timesteps/s
99.6% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 2.462 | 4.2266 | 7.1075 | 89.3 | 46.68
Bond | 0.00018424 | 0.00019878 | 0.00022125 | 0.0 | 0.00
Kspace | 1.4698 | 4.338 | 6.1001 | 87.8 | 47.91
Neigh | 0.39462 | 0.39489 | 0.39518 | 0.0 | 4.36
Comm | 0.043753 | 0.055826 | 0.062746 | 3.1 | 0.62
Output | 0.00048755 | 0.00053971 | 0.00062785 | 0.0 | 0.01
Modify | 0.027255 | 0.028664 | 0.030278 | 0.7 | 0.32
Other | | 0.008904 | | | 0.10
Nlocal: 234 ave 302 max 198 min
Histogram: 2 0 0 1 0 0 0 0 0 1
Nghost: 2876.5 ave 3122 max 2632 min
Histogram: 1 0 1 0 0 0 0 1 0 1
Neighs: 0 ave 0 max 0 min
Histogram: 4 0 0 0 0 0 0 0 0 0
FullNghs: 107846 ave 150684 max 82181 min
Histogram: 2 0 0 0 1 0 0 0 0 1
Total # of neighbors = 431382
Ave neighs/atom = 460.87821
Ave special neighs/atom = 0.30769231
Neighbor list builds = 28
Dangerous builds = 0
Total wall time: 0:00:09

View File

@ -0,0 +1,256 @@
LAMMPS (23 Jun 2022 - Update 4)
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:98)
using 1 OpenMP thread(s) per MPI task
# Initialization
units metal
boundary p p p
atom_style full
processors * * 1 # domain decomposition over x and y
read_data ./gra_water.data
Reading data file ...
orthogonal box = (0 0 0) to (46.92336 44.331078 200)
1 by 1 by 1 MPI processor grid
reading atoms ...
936 atoms
reading velocities ...
936 velocities
scanning bonds ...
2 = max bonds/atom
scanning angles ...
1 = max angles/atom
reading bonds ...
96 bonds
reading angles ...
48 angles
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 0 0
special bond factors coul: 0 0 0
2 = max # of 1-2 neighbors
1 = max # of 1-3 neighbors
1 = max # of 1-4 neighbors
2 = max # of special neighbors
special bonds CPU = 0.000 seconds
read_data CPU = 0.012 seconds
mass 1 12.0107 # carbon mass (g/mole)
mass 2 15.9994 # oxygen mass (g/mole)
mass 3 1.008 # hydrogen mass (g/mole)
# Separate atom groups
group gr molecule 1
792 atoms in group gr
group water molecule 2
144 atoms in group water
######################## Potential defition ##############################
# Interlayer potential
pair_style hybrid/overlay aip/water/2dm/opt 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
####################################################################
pair_coeff 1 1 none
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm/opt COH.aip.water.2dm C Ow Hw # C-H2O
Reading aip/water/2dm potential file COH.aip.water.2dm with DATE: 2022-12-02
# bond and angle
bond_style harmonic
bond_coeff 1 0.0 0.9572
angle_style harmonic
angle_coeff 1 0.0 104.52
# define kspace calculation
kspace_style pppm/tip4p 1E-5
# Neighbor update settings
neighbor 2.0 bin
neigh_modify every 1 delay 5 check yes page 1000000 one 100000
####################################################################
# Calculate pair energy
compute 1 all pair lj/cut/tip4p/long
compute 2 all pair aip/water/2dm/opt
compute wt water temp
variable TIP4P equal c_1
variable EILP equal c_2 # total interlayer energy
variable temp_wt equal c_wt
############# Output ##############
thermo_style custom step etotal pe ke v_TIP4P v_EILP v_temp_wt
thermo 100
thermo_modify lost error
fix subf gr setforce 0.0 0.0 0.0
fix 1 water shake 0.0001 20 100 b 1 a 1
0 = # of size 2 clusters
0 = # of size 3 clusters
0 = # of size 4 clusters
48 = # of frozen angles
find clusters CPU = 0.000 seconds
timestep 1e-3
velocity water create 300.0 12345 dist gaussian mom yes rot yes
fix 2 water nve
run 1000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- ilp/graphene/hbn potential doi:10.1021/acs.nanolett.8b02848
@Article{Ouyang2018
author = {W. Ouyang and D. Mandelli and M. Urbakh and O. Hod},
title = {Nanoserpents: Graphene Nanoribbon Motion on Two-Dimensional Hexagonal Materials},
journal = {Nano Letters},
volume = 18,
pages = 6009,
year = 2018,
}
- ilp/tmd potential doi:10.1021/acs.jctc.1c00782
@Article{Ouyang2021
author = {W. Ouyang and R. Sofer and X. Gao and J. Hermann and
A. Tkatchenko and L. Kronik and M. Urbakh and O. Hod},
title = {Anisotropic Interlayer Force Field for Transition
Metal Dichalcogenides: The Case of Molybdenum Disulfide},
journal = {J.~Chem.\ Theory Comput.},
volume = 17,
pages = {7237--7245}
year = 2021,
}
- ilp/water/2dm potential doi/10.1021/acs.jpcc.2c08464
@Article{Feng2023
author = {Z. Feng, Y. Yao, J. Liu, B. Wu, Z. Liu, and W. Ouyang},
title = {Registry-Dependent Potential for Interfaces of Water with Graphene},
journal = {J. Phys. Chem. C},
volume = 127,
pages = {8704-8713}
year = 2023,
}
- ilp/graphene/hbn/opt potential doi:10.1145/3458817.3476137
@inproceedings{gao2021lmff
author = {Gao, Ping and Duan, Xiaohui and Others},
title = {LMFF: Efficient and Scalable Layered Materials Force Field on Heterogeneous Many-Core Processors},
year = {2021},
isbn = {9781450384421},
publisher = {Association for Computing Machinery},
address = {New York, NY, USA},
url = {https://doi.org/10.1145/3458817.3476137},
doi = {10.1145/3458817.3476137},
booktitle = {Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis},
articleno = {42},
numpages = {14},
location = {St. Louis, Missouri},
series = {SC'21},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
PPPM initialization ...
extracting TIP4P info from pair style
using 12-bit tables for long-range coulomb (../kspace.cpp:342)
G vector (1/distance) = 0.28684806
grid = 25 24 80
stencil order = 5
estimated absolute RMS force accuracy = 0.0001640931
estimated relative force accuracy = 1.1395635e-05
using single precision MKL FFT
3d grid and FFT values/proc = 84320 48000
WARNING: Using a manybody potential with bonds/angles/dihedrals and special_bond exclusions (../pair.cpp:239)
WARNING: Communication cutoff 0 is shorter than a bond length based estimate of 3.4358. This may lead to errors. (../comm.cpp:727)
WARNING: Increasing communication cutoff to 11.6118 for TIP4P pair style (../pair_lj_cut_tip4p_long.cpp:484)
Neighbor list info ...
update every 1 steps, delay 5 steps, check yes
max neighbors/atom: 100000, page size: 1000000
master list distance cutoff = 18
ghost atom cutoff = 18
binsize = 9, bins = 6 5 23
3 neighbor lists, perpetual/occasional/extra = 3 0 0
(1) pair aip/water/2dm/opt, perpetual
attributes: full, newton on
pair build: full/bin
stencil: full/bin/3d
bin: standard
(2) pair lj/cut/tip4p/long, perpetual, skip from (3)
attributes: half, newton on
pair build: skip
stencil: none
bin: none
(3) neighbor class addition, perpetual, half/full from (1)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
WARNING: Communication cutoff adjusted to 18 (../comm.cpp:736)
SHAKE stats (type/ave/delta/count) on step 0
Bond: 1 0.957201 2.19705e-06 96
Angle: 1 104.52 0.000203056 48
Per MPI rank memory allocation (min/avg/max) = 25.27 | 25.27 | 25.27 Mbytes
Step TotEng PotEng KinEng v_TIP4P v_EILP v_temp_wt
0 -16.131263 -19.815178 3.6839141 189.37246 -1.903543 300
SHAKE stats (type/ave/delta/count) on step 100
Bond: 1 0.9572 9.54949e-07 96
Angle: 1 104.52 6.01522e-05 48
100 -17.494868 -20.796993 3.3021253 188.4955 -1.8981262 268.90898
SHAKE stats (type/ave/delta/count) on step 200
Bond: 1 0.9572 9.63922e-07 96
Angle: 1 104.52 7.7021e-05 48
200 -17.486271 -21.194892 3.7086213 188.14561 -1.9871708 302.01203
SHAKE stats (type/ave/delta/count) on step 300
Bond: 1 0.9572 1.4264e-06 96
Angle: 1 104.52 6.48393e-05 48
300 -17.502844 -20.993704 3.49086 188.23268 -1.8457229 284.27861
SHAKE stats (type/ave/delta/count) on step 400
Bond: 1 0.9572 1.33728e-06 96
Angle: 1 104.52 7.6239e-05 48
400 -17.495287 -20.828353 3.3330658 188.48002 -1.8429075 271.42862
SHAKE stats (type/ave/delta/count) on step 500
Bond: 1 0.9572 1.14685e-06 96
Angle: 1 104.52 8.58621e-05 48
500 -17.491436 -20.443044 2.9516084 188.7589 -1.8566335 240.36459
SHAKE stats (type/ave/delta/count) on step 600
Bond: 1 0.9572 9.17601e-07 96
Angle: 1 104.52 8.24516e-05 48
600 -17.505684 -20.608457 3.1027731 188.72078 -1.9560796 252.67471
SHAKE stats (type/ave/delta/count) on step 700
Bond: 1 0.9572 9.50422e-07 96
Angle: 1 104.52 5.62423e-05 48
700 -17.496701 -21.07266 3.5759595 188.2777 -1.9833956 291.20871
SHAKE stats (type/ave/delta/count) on step 800
Bond: 1 0.9572 1.15262e-06 96
Angle: 1 104.52 7.02158e-05 48
800 -17.478623 -20.819504 3.340881 188.37868 -1.9112996 272.06506
SHAKE stats (type/ave/delta/count) on step 900
Bond: 1 0.9572 9.14138e-07 96
Angle: 1 104.52 6.98742e-05 48
900 -17.480864 -20.728498 3.2476343 188.59022 -1.8922102 264.4715
SHAKE stats (type/ave/delta/count) on step 1000
Bond: 1 0.9572 1.00586e-06 96
Angle: 1 104.52 0.000111711 48
1000 -17.498466 -20.331547 2.8330808 188.87473 -1.8121768 230.71228
Loop time of 8.11929 on 1 procs for 1000 steps with 936 atoms
Performance: 10.641 ns/day, 2.255 hours/ns, 123.163 timesteps/s
99.7% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 4.6849 | 4.6849 | 4.6849 | 0.0 | 57.70
Bond | 0.00017678 | 0.00017678 | 0.00017678 | 0.0 | 0.00
Kspace | 3.2098 | 3.2098 | 3.2098 | 0.0 | 39.53
Neigh | 0.17546 | 0.17546 | 0.17546 | 0.0 | 2.16
Comm | 0.024693 | 0.024693 | 0.024693 | 0.0 | 0.30
Output | 0.00037798 | 0.00037798 | 0.00037798 | 0.0 | 0.00
Modify | 0.015853 | 0.015853 | 0.015853 | 0.0 | 0.20
Other | | 0.007983 | | | 0.10
Nlocal: 936 ave 936 max 936 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 5242 ave 5242 max 5242 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 0 ave 0 max 0 min
Histogram: 1 0 0 0 0 0 0 0 0 0
FullNghs: 431382 ave 431382 max 431382 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 431382
Ave neighs/atom = 460.87821
Ave special neighs/atom = 0.30769231
Neighbor list builds = 28
Dangerous builds = 0
Total wall time: 0:00:08

View File

@ -0,0 +1,256 @@
LAMMPS (23 Jun 2022 - Update 4)
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:98)
using 1 OpenMP thread(s) per MPI task
# Initialization
units metal
boundary p p p
atom_style full
processors * * 1 # domain decomposition over x and y
read_data ./gra_water.data
Reading data file ...
orthogonal box = (0 0 0) to (46.92336 44.331078 200)
2 by 2 by 1 MPI processor grid
reading atoms ...
936 atoms
reading velocities ...
936 velocities
scanning bonds ...
2 = max bonds/atom
scanning angles ...
1 = max angles/atom
reading bonds ...
96 bonds
reading angles ...
48 angles
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 0 0
special bond factors coul: 0 0 0
2 = max # of 1-2 neighbors
1 = max # of 1-3 neighbors
1 = max # of 1-4 neighbors
2 = max # of special neighbors
special bonds CPU = 0.001 seconds
read_data CPU = 0.014 seconds
mass 1 12.0107 # carbon mass (g/mole)
mass 2 15.9994 # oxygen mass (g/mole)
mass 3 1.008 # hydrogen mass (g/mole)
# Separate atom groups
group gr molecule 1
792 atoms in group gr
group water molecule 2
144 atoms in group water
######################## Potential defition ##############################
# Interlayer potential
pair_style hybrid/overlay aip/water/2dm/opt 16.0 lj/cut/tip4p/long 2 3 1 1 0.1546 10 8.5
####################################################################
pair_coeff 1 1 none
pair_coeff 2 2 lj/cut/tip4p/long 8.0313e-3 3.1589 # O-O
pair_coeff 2 3 lj/cut/tip4p/long 0.0 0.0 # O-H
pair_coeff 3 3 lj/cut/tip4p/long 0.0 0.0 # H-H
pair_coeff * * aip/water/2dm/opt COH.aip.water.2dm C Ow Hw # C-H2O
Reading aip/water/2dm potential file COH.aip.water.2dm with DATE: 2022-12-02
# bond and angle
bond_style harmonic
bond_coeff 1 0.0 0.9572
angle_style harmonic
angle_coeff 1 0.0 104.52
# define kspace calculation
kspace_style pppm/tip4p 1E-5
# Neighbor update settings
neighbor 2.0 bin
neigh_modify every 1 delay 5 check yes page 1000000 one 100000
####################################################################
# Calculate pair energy
compute 1 all pair lj/cut/tip4p/long
compute 2 all pair aip/water/2dm/opt
compute wt water temp
variable TIP4P equal c_1
variable EILP equal c_2 # total interlayer energy
variable temp_wt equal c_wt
############# Output ##############
thermo_style custom step etotal pe ke v_TIP4P v_EILP v_temp_wt
thermo 100
thermo_modify lost error
fix subf gr setforce 0.0 0.0 0.0
fix 1 water shake 0.0001 20 100 b 1 a 1
0 = # of size 2 clusters
0 = # of size 3 clusters
0 = # of size 4 clusters
48 = # of frozen angles
find clusters CPU = 0.000 seconds
timestep 1e-3
velocity water create 300.0 12345 dist gaussian mom yes rot yes
fix 2 water nve
run 1000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- ilp/graphene/hbn potential doi:10.1021/acs.nanolett.8b02848
@Article{Ouyang2018
author = {W. Ouyang and D. Mandelli and M. Urbakh and O. Hod},
title = {Nanoserpents: Graphene Nanoribbon Motion on Two-Dimensional Hexagonal Materials},
journal = {Nano Letters},
volume = 18,
pages = 6009,
year = 2018,
}
- ilp/tmd potential doi:10.1021/acs.jctc.1c00782
@Article{Ouyang2021
author = {W. Ouyang and R. Sofer and X. Gao and J. Hermann and
A. Tkatchenko and L. Kronik and M. Urbakh and O. Hod},
title = {Anisotropic Interlayer Force Field for Transition
Metal Dichalcogenides: The Case of Molybdenum Disulfide},
journal = {J.~Chem.\ Theory Comput.},
volume = 17,
pages = {7237--7245}
year = 2021,
}
- ilp/water/2dm potential doi/10.1021/acs.jpcc.2c08464
@Article{Feng2023
author = {Z. Feng, Y. Yao, J. Liu, B. Wu, Z. Liu, and W. Ouyang},
title = {Registry-Dependent Potential for Interfaces of Water with Graphene},
journal = {J. Phys. Chem. C},
volume = 127,
pages = {8704-8713}
year = 2023,
}
- ilp/graphene/hbn/opt potential doi:10.1145/3458817.3476137
@inproceedings{gao2021lmff
author = {Gao, Ping and Duan, Xiaohui and Others},
title = {LMFF: Efficient and Scalable Layered Materials Force Field on Heterogeneous Many-Core Processors},
year = {2021},
isbn = {9781450384421},
publisher = {Association for Computing Machinery},
address = {New York, NY, USA},
url = {https://doi.org/10.1145/3458817.3476137},
doi = {10.1145/3458817.3476137},
booktitle = {Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis},
articleno = {42},
numpages = {14},
location = {St. Louis, Missouri},
series = {SC'21},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
PPPM initialization ...
extracting TIP4P info from pair style
using 12-bit tables for long-range coulomb (../kspace.cpp:342)
G vector (1/distance) = 0.28684806
grid = 25 24 80
stencil order = 5
estimated absolute RMS force accuracy = 0.0001640931
estimated relative force accuracy = 1.1395635e-05
using single precision MKL FFT
3d grid and FFT values/proc = 30685 12480
WARNING: Using a manybody potential with bonds/angles/dihedrals and special_bond exclusions (../pair.cpp:239)
WARNING: Communication cutoff 0 is shorter than a bond length based estimate of 3.4358. This may lead to errors. (../comm.cpp:727)
WARNING: Increasing communication cutoff to 11.6118 for TIP4P pair style (../pair_lj_cut_tip4p_long.cpp:484)
Neighbor list info ...
update every 1 steps, delay 5 steps, check yes
max neighbors/atom: 100000, page size: 1000000
master list distance cutoff = 18
ghost atom cutoff = 18
binsize = 9, bins = 6 5 23
3 neighbor lists, perpetual/occasional/extra = 3 0 0
(1) pair aip/water/2dm/opt, perpetual
attributes: full, newton on
pair build: full/bin
stencil: full/bin/3d
bin: standard
(2) pair lj/cut/tip4p/long, perpetual, skip from (3)
attributes: half, newton on
pair build: skip
stencil: none
bin: none
(3) neighbor class addition, perpetual, half/full from (1)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
WARNING: Communication cutoff adjusted to 18 (../comm.cpp:736)
SHAKE stats (type/ave/delta/count) on step 0
Bond: 1 0.957201 2.19705e-06 96
Angle: 1 104.52 0.000203056 48
Per MPI rank memory allocation (min/avg/max) = 21.4 | 21.44 | 21.48 Mbytes
Step TotEng PotEng KinEng v_TIP4P v_EILP v_temp_wt
0 -16.131263 -19.815178 3.6839141 189.37246 -1.903543 300
SHAKE stats (type/ave/delta/count) on step 100
Bond: 1 0.9572 9.54949e-07 96
Angle: 1 104.52 6.01522e-05 48
100 -17.494869 -20.796995 3.3021253 188.4955 -1.8981262 268.90898
SHAKE stats (type/ave/delta/count) on step 200
Bond: 1 0.9572 9.63922e-07 96
Angle: 1 104.52 7.7021e-05 48
200 -17.48627 -21.194892 3.7086213 188.14561 -1.9871708 302.01203
SHAKE stats (type/ave/delta/count) on step 300
Bond: 1 0.9572 1.4264e-06 96
Angle: 1 104.52 6.48393e-05 48
300 -17.502843 -20.993703 3.4908599 188.23268 -1.8457229 284.27861
SHAKE stats (type/ave/delta/count) on step 400
Bond: 1 0.9572 1.33728e-06 96
Angle: 1 104.52 7.6239e-05 48
400 -17.495285 -20.82835 3.333065 188.48003 -1.8429074 271.42856
SHAKE stats (type/ave/delta/count) on step 500
Bond: 1 0.9572 1.14685e-06 96
Angle: 1 104.52 8.58621e-05 48
500 -17.491436 -20.443043 2.9516075 188.7589 -1.8566335 240.36452
SHAKE stats (type/ave/delta/count) on step 600
Bond: 1 0.9572 9.17601e-07 96
Angle: 1 104.52 8.24517e-05 48
600 -17.505682 -20.608456 3.1027734 188.72078 -1.9560795 252.67474
SHAKE stats (type/ave/delta/count) on step 700
Bond: 1 0.9572 9.50425e-07 96
Angle: 1 104.52 5.62423e-05 48
700 -17.496706 -21.072664 3.575958 188.2777 -1.9833951 291.20858
SHAKE stats (type/ave/delta/count) on step 800
Bond: 1 0.9572 1.15256e-06 96
Angle: 1 104.52 7.02177e-05 48
800 -17.478628 -20.819507 3.340879 188.37868 -1.9113009 272.0649
SHAKE stats (type/ave/delta/count) on step 900
Bond: 1 0.9572 9.14163e-07 96
Angle: 1 104.52 6.98849e-05 48
900 -17.480868 -20.728506 3.2476383 188.5902 -1.8922108 264.47182
SHAKE stats (type/ave/delta/count) on step 1000
Bond: 1 0.9572 1.00568e-06 96
Angle: 1 104.52 0.000111707 48
1000 -17.498472 -20.331605 2.8331335 188.87466 -1.8121689 230.71657
Loop time of 4.24862 on 4 procs for 1000 steps with 936 atoms
Performance: 20.336 ns/day, 1.180 hours/ns, 235.370 timesteps/s
99.6% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.25749 | 1.1806 | 2.6872 | 89.0 | 27.79
Bond | 0.00018656 | 0.00020786 | 0.00025377 | 0.0 | 0.00
Kspace | 1.4259 | 2.9204 | 3.8414 | 56.3 | 68.74
Neigh | 0.057504 | 0.057852 | 0.05818 | 0.1 | 1.36
Comm | 0.041952 | 0.053593 | 0.05876 | 3.0 | 1.26
Output | 0.0004296 | 0.00046809 | 0.00055317 | 0.0 | 0.01
Modify | 0.026204 | 0.027251 | 0.028382 | 0.6 | 0.64
Other | | 0.008209 | | | 0.19
Nlocal: 234 ave 302 max 198 min
Histogram: 2 0 0 1 0 0 0 0 0 1
Nghost: 2876.5 ave 3122 max 2632 min
Histogram: 1 0 1 0 0 0 0 1 0 1
Neighs: 0 ave 0 max 0 min
Histogram: 4 0 0 0 0 0 0 0 0 0
FullNghs: 107846 ave 150684 max 82181 min
Histogram: 2 0 0 0 1 0 0 0 0 1
Total # of neighbors = 431382
Ave neighs/atom = 460.87821
Ave special neighs/atom = 0.30769231
Neighbor list builds = 28
Dangerous builds = 0
Total wall time: 0:00:04

View File

@ -0,0 +1,28 @@
# DATE: 2022-12-02 UNITS: metal CONTRIBUTOR: Wengen Ouyang w.g.ouyang@gmail.com CITATION: Z. Feng, ..., and W. Ouyang, J. Phys. Chem. C 127, 8704 (2023).
# Anisotropic Interfacial Potential (AIP) parameters for water/graphene heterojunctions
# The parameters below are fitted against the DMC reference data that rescaled from PBE+MBD-NL.
#
# ----------------- Repulsion Potential ------------------++++++++++++++ Vdw Potential ++++++++++++++++************
# beta(A) alpha delta(A) epsilon(meV) C(meV) d sR reff(A) C6(meV*A^6) S rcut
# For graphene and hydrocarbons
C C 3.205843 7.511126 1.235334 1.528338E-5 37.530428 15.499947 0.7954443 3.681440 25.714535E3 1.0 2.0
H H 3.974540 6.53799 1.080633 0.6700556 0.8333833 15.022371 0.7490632 2.767223 1.6159581E3 1.0 1.2
C H 2.642950 12.91410 1.020257 0.9750012 25.340996 15.222927 0.8115998 3.887324 5.6874617E3 1.0 1.5
H C 2.642950 12.91410 1.020257 0.9750012 25.340996 15.222927 0.8115998 3.887324 5.6874617E3 1.0 1.5
# For water-graphene
C Ow 4.03686677 15.09817069 1.00000323 1.03838111 0.16372013 9.00010553 1.12057182 2.96484087 21887.14173222 1.0 2.0
C Hw 3.08246994 6.412090180 1.00000265 2.89390420 -1.85748759 9.00009101 1.04574423 2.21642099 4652.78021666 1.0 2.0
Ow C 4.03686677 15.09817069 1.00000323 1.03838111 0.16372013 9.00010553 1.12057182 2.96484087 21887.14173222 1.0 1.2
Hw C 3.08246994 6.412090180 1.00000265 2.89390420 -1.85748759 9.00009101 1.04574423 2.21642099 4652.78021666 1.0 1.2
# # The ILPs for other systems are set to zero
H Ow 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
H Hw 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Ow H 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Hw H 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Ow Ow 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Hw Hw 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Ow Hw 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Hw Ow 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2

View File

@ -0,0 +1,28 @@
# DATE: 2022-12-02 UNITS: metal CONTRIBUTOR: Wengen Ouyang w.g.ouyang@gmail.com CITATION: Z. Feng, ..., and W. Ouyang, J. Phys. Chem. C 127, 8704 (2023).
# Anisotropic Interfacial Potential (AIP) parameters for water/graphene heterojunctions
# The parameters below are fitted against the PBE + MBD-NL DFT reference data from 2.5 A to 15 A.
#
# ----------------- Repulsion Potential ------------------++++++++++++++ Vdw Potential ++++++++++++++++************
# beta(A) alpha delta(A) epsilon(meV) C(meV) d sR reff(A) C6(meV*A^6) S rcut
# For graphene and hydrocarbons
C C 3.205843 7.511126 1.235334 1.528338E-5 37.530428 15.499947 0.7954443 3.681440 25.714535E3 1.0 2.0
H H 3.974540 6.53799 1.080633 0.6700556 0.8333833 15.022371 0.7490632 2.767223 1.6159581E3 1.0 1.2
C H 2.642950 12.91410 1.020257 0.9750012 25.340996 15.222927 0.8115998 3.887324 5.6874617E3 1.0 1.5
H C 2.642950 12.91410 1.020257 0.9750012 25.340996 15.222927 0.8115998 3.887324 5.6874617E3 1.0 1.5
# For water-graphene
C Ow 5.45369612 6.18172364 1.25025450 3.34909245 0.68780636 9.05706482 1.23249498 2.77577173 100226.55503127 1.0 2.0
C Hw 2.55380862 9.68664390 1.96489198 41.77617053 -16.30012807 9.01568534 0.74415463 2.41545571 7409.12856378 1.0 2.0
Ow C 5.45369612 6.18172364 1.25025450 3.34909245 0.68780636 9.05706482 1.23249498 2.77577173 100226.55503127 1.0 1.2
Hw C 2.55380862 9.68664390 1.96489198 41.77617053 -16.30012807 9.01568534 0.74415463 2.41545571 7409.12856378 1.0 1.2
# # The ILPs for other systems are set to zero
H Ow 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
H Hw 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Ow H 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Hw H 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Ow Ow 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Hw Hw 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Ow Hw 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2
Hw Ow 5.45369612 6.18172364 1.25025450 0.00000000 0.00000000 9.05706482 1.23249498 2.77577173 0.00000000 1.0 1.2

4
src/.gitignore vendored
View File

@ -1058,6 +1058,10 @@
/pair_adp.h
/pair_agni.cpp
/pair_agni.h
/pair_aip_water_2dm.cpp
/pair_aip_water_2dm.h
/pair_aip_water_2dm_opt.cpp
/pair_aip_water_2dm_opt.h
/pair_airebo.cpp
/pair_airebo.h
/pair_airebo_morse.cpp

View File

@ -0,0 +1,72 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://www.lammps.org/, Sandia National Laboratories
LAMMPS development team: developers@lammps.org
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Wengen Ouyang (Wuhan University)
e-mail: w.g.ouyang at gmail dot com
This is a full version of the potential described in
[Feng and Ouyang et al, J. Phys. Chem. C 127, 8704-8713 (2023).]
------------------------------------------------------------------------- */
#include "pair_aip_water_2dm.h"
#include "citeme.h"
#include "error.h"
#include "force.h"
#include <cmath>
#include <cstring>
using namespace LAMMPS_NS;
#define MAXLINE 1024
#define DELTA 4
#define PGDELTA 1
static const char cite_aip_water[] =
"aip/water/2dm potential doi/10.1021/acs.jpcc.2c08464\n"
"@Article{Feng2023\n"
" author = {Z. Feng, Y. Yao, J. Liu, B. Wu, Z. Liu, and W. Ouyang},\n"
" title = {Registry-Dependent Potential for Interfaces of Water with Graphene},\n"
" journal = {J. Phys. Chem. C},\n"
" volume = 127,\n"
" pages = {8704-8713}\n"
" year = 2023,\n"
"}\n\n";
/* ---------------------------------------------------------------------- */
PairAIPWATER2DM::PairAIPWATER2DM(LAMMPS *lmp) : PairILPGrapheneHBN(lmp), PairILPTMD(lmp)
{
variant = AIP_WATER_2DM;
single_enable = 0;
// for TMD, each atom have six neighbors
Nnei = 6;
if (lmp->citeme) lmp->citeme->add(cite_aip_water);
}
/* ----------------------------------------------------------------------
global settings
------------------------------------------------------------------------- */
void PairAIPWATER2DM::settings(int narg, char **arg)
{
if (narg < 1 || narg > 2) error->all(FLERR, "Illegal pair_style command");
if (!utils::strmatch(force->pair_style, "^hybrid/overlay"))
error->all(FLERR, "Pair style aip/water/2dm must be used as sub-style with hybrid/overlay");
cut_global = utils::numeric(FLERR, arg[0], false, lmp);
if (narg == 2) tap_flag = utils::numeric(FLERR, arg[1], false, lmp);
}

View File

@ -0,0 +1,39 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://www.lammps.org/, Sandia National Laboratories
LAMMPS development team: developers@lammps.org
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifdef PAIR_CLASS
// clang-format off
PairStyle(aip/water/2dm,PairAIPWATER2DM);
// clang-format on
#else
#ifndef LMP_PAIR_AIP_WATER_2DM_H
#define LMP_PAIR_AIP_WATER_2DM_H
#include "pair_ilp_tmd.h"
namespace LAMMPS_NS {
class PairAIPWATER2DM : virtual public PairILPTMD {
public:
PairAIPWATER2DM(class LAMMPS *);
protected:
void settings(int, char **) override;
};
} // namespace LAMMPS_NS
#endif
#endif

View File

@ -11,7 +11,7 @@
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Wengen Ouyang (Tel Aviv University)
Contributing author: Wengen Ouyang (Wuhan University)
e-mail: w.g.ouyang at gmail dot com
This is a Coulomb potential described in

View File

@ -11,13 +11,11 @@
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Wengen Ouyang (Tel Aviv University)
Contributing author: Wengen Ouyang (Wuhan University)
e-mail: w.g.ouyang at gmail dot com
This is a full version of the potential described in
[Maaravi et al, J. Phys. Chem. C 121, 22826-22835 (2017)]
The definition of normals are the same as that in
[Kolmogorov & Crespi, Phys. Rev. B 71, 235415 (2005)]
[Ouyang et al., J. Chem. Theory Comput. 16(1), 666-676 (2020)]
------------------------------------------------------------------------- */
#include "pair_ilp_graphene_hbn.h"
@ -59,6 +57,7 @@ static const char cite_ilp[] =
static std::map<int, std::string> variant_map = {
{PairILPGrapheneHBN::ILP_GrhBN, "ilp/graphene/hbn"},
{PairILPGrapheneHBN::ILP_TMD, "ilp/tmd"},
{PairILPGrapheneHBN::AIP_WATER_2DM, "aip/water/2dm"},
{PairILPGrapheneHBN::SAIP_METAL, "saip/metal"}};
/* ---------------------------------------------------------------------- */
@ -632,7 +631,7 @@ void PairILPGrapheneHBN::calc_FRep(int eflag, int /* vflag */)
void PairILPGrapheneHBN::ILP_neigh()
{
int i, j, ii, jj, n, allnum, jnum, itype, jtype;
int i, j, ii, jj, n, inum, jnum, itype, jtype;
double xtmp, ytmp, ztmp, delx, dely, delz, rsq;
int *ilist, *jlist, *numneigh, **firstneigh;
int *neighptr;
@ -649,7 +648,7 @@ void PairILPGrapheneHBN::ILP_neigh()
(int **) memory->smalloc(maxlocal * sizeof(int *), "ILPGrapheneHBN:firstneigh");
}
allnum = list->inum + list->gnum;
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
@ -659,7 +658,7 @@ void PairILPGrapheneHBN::ILP_neigh()
ipage->reset();
for (ii = 0; ii < allnum; ii++) {
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
n = 0;

View File

@ -39,7 +39,7 @@ class PairILPGrapheneHBN : public Pair {
static constexpr int NPARAMS_PER_LINE = 13;
enum { ILP_GrhBN, ILP_TMD, SAIP_METAL }; // for telling class variants apart in shared code
enum { ILP_GrhBN, ILP_TMD, SAIP_METAL, AIP_WATER_2DM }; // for telling class variants apart in shared code
protected:
int me;

View File

@ -15,7 +15,7 @@
e-mail: w.g.ouyang at gmail dot com
This is a full version of the potential described in
[Ouyang et al, J. Chem. Theory Comput. 17, 7237 (2021).]
[Ouyang et al., J. Chem. Theory Comput. 17, 7237 (2021).]
------------------------------------------------------------------------- */
#include "pair_ilp_tmd.h"
@ -35,10 +35,6 @@
using namespace LAMMPS_NS;
using namespace InterLayer;
#define MAXLINE 1024
#define DELTA 4
#define PGDELTA 1
static const char cite_ilp_tmd[] =
"ilp/tmd potential doi:10.1021/acs.jctc.1c00782\n"
"@Article{Ouyang2021\n"
@ -232,7 +228,7 @@ void PairILPTMD::calc_FRep(int eflag, int /* vflag */)
void PairILPTMD::ILP_neigh()
{
int i, j, l, ii, jj, ll, n, allnum, jnum, itype, jtype, ltype, imol, jmol, count;
int i, j, l, ii, jj, ll, n, inum, jnum, itype, jtype, ltype, imol, jmol, count;
double xtmp, ytmp, ztmp, delx, dely, delz, deljx, deljy, deljz, rsq, rsqlj;
int *ilist, *jlist, *numneigh, **firstneigh;
int *neighsort;
@ -249,7 +245,7 @@ void PairILPTMD::ILP_neigh()
ILP_firstneigh = (int **) memory->smalloc(maxlocal * sizeof(int *), "ILPTMD:firstneigh");
}
allnum = list->inum + list->gnum;
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
@ -259,7 +255,7 @@ void PairILPTMD::ILP_neigh()
ipage->reset();
for (ii = 0; ii < allnum; ii++) {
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
//initialize varibles
@ -288,14 +284,14 @@ void PairILPTMD::ILP_neigh()
delz = ztmp - x[j][2];
rsq = delx * delx + dely * dely + delz * delz;
// check if the atom i is TMD, i.e., Mo/S/W/Se
// check if the atom i is a TMD atom, i.e., Mo/S/W/Se
if (strcmp(elements[itype], "Mo") == 0 || strcmp(elements[itype], "W") == 0 ||
strcmp(elements[itype], "S") == 0 || strcmp(elements[itype], "Se") == 0 ||
strcmp(elements[itype], "Te") == 0) {
if (rsq != 0 && rsq < cutILPsq[itype][jtype] && imol == jmol && type[i] == type[j]) {
neighptr[n++] = j;
}
} else { // atom i is C, B, N or H.
} else { // atom i can be P, C, B, N or H.
if (rsq != 0 && rsq < cutILPsq[itype][jtype] && imol == jmol) { neighptr[n++] = j; }
}
} // loop over jj
@ -336,9 +332,6 @@ void PairILPTMD::ILP_neigh()
}
} // end of idenfying the first neighbor
} else if (n > Nnei) {
fprintf(screen, "Molecule ID = %d\n", imol);
fprintf(screen, "Atom Type = %d\n", type[i]);
fprintf(screen, "Neinum = %d\n", n);
error->one(FLERR,
"There are too many neighbors for TMD atoms, please check your configuration");
}
@ -367,11 +360,11 @@ void PairILPTMD::ILP_neigh()
ll++;
}
} // end of sorting the order of neighbors
} else { // for B/N/C/H atoms
} else { // for P/B/N/C/H atoms
if (n > 3)
error->one(
FLERR,
"There are too many neighbors for B/N/C/H atoms, please check your configuration");
"There are too many neighbors for P/B/N/C/H atoms, please check your configuration");
for (ll = 0; ll < n; ll++) { neighsort[ll] = neighptr[ll]; }
}
@ -390,12 +383,14 @@ void PairILPTMD::calc_normal()
{
int i, j, ii, jj, inum, jnum;
int cont, id, ip, m, k, itype;
double nn, xtp, ytp, ztp, delx, dely, delz, nn2;
int *ilist, *jlist;
int iH1,iH2,jH1,jH2;
double Nave[3], dni[3], dpvdri[3][3];
double nn, xtp, ytp, ztp, delx, dely, delz, nn2;
double **x = atom->x;
int *type = atom->type;
tagint *tag = atom->tag;
memory->destroy(dnn);
memory->destroy(vect);
@ -479,9 +474,60 @@ void PairILPTMD::calc_normal()
for (m = 0; m < Nnei; m++) { dnormal[i][id][m][ip] = 0.0; }
}
}
// for hydrogen in water molecule
if (strcmp(elements[itype], "Hw") == 0) {
if(cont == 0) {
jH1 = atom->map(tag[i] - 1);
jH2 = atom->map(tag[i] - 2);
iH1 = map[type[jH1]];
iH2 = map[type[jH2]];
if (strcmp(elements[iH1], "Ow") == 0 ) {
vect[0][0] = x[jH1][0] - xtp;
vect[0][1] = x[jH1][1] - ytp;
vect[0][2] = x[jH1][2] - ztp;
} else if (strcmp(elements[iH2], "Ow") == 0 ) {
vect[0][0] = x[jH2][0] - xtp;
vect[0][1] = x[jH2][1] - ytp;
vect[0][2] = x[jH2][2] - ztp;
} else {
error->one(FLERR, "The order of atoms in water molecule should be O H H !");
}
}
Nave[0] = vect[0][0];
Nave[1] = vect[0][1];
Nave[2] = vect[0][2];
// the magnitude of the normal vector
nn2 = Nave[0] * Nave[0] + Nave[1] * Nave[1] + Nave[2] * Nave[2];
nn = sqrt(nn2);
if (nn == 0) error->one(FLERR, "The magnitude of the normal vector is zero");
// the unit normal vector
normal[i][0] = Nave[0] / nn;
normal[i][1] = Nave[1] / nn;
normal[i][2] = Nave[2] / nn;
// Calculte dNave/dri, defined as dpvdri
for (id = 0; id < 3; id++) {
for (ip = 0; ip < 3; ip++) {
if (ip == id) { dpvdri[id][ip] = -1.0;}
else {dpvdri[id][ip] = 0.0;}
}
}
// derivatives of nn, dnn:3x1 vector
dni[0] = (Nave[0] * dpvdri[0][0] + Nave[1] * dpvdri[1][0] + Nave[2] * dpvdri[2][0]) / nn;
dni[1] = (Nave[0] * dpvdri[0][1] + Nave[1] * dpvdri[1][1] + Nave[2] * dpvdri[2][1]) / nn;
dni[2] = (Nave[0] * dpvdri[0][2] + Nave[1] * dpvdri[1][2] + Nave[2] * dpvdri[2][2]) / nn;
// derivatives of unit vector ni respect to ri, the result is 3x3 matrix
for (id = 0; id < 3; id++) {
for (ip = 0; ip < 3; ip++) {
dnormdri[i][id][ip] = dpvdri[id][ip] / nn - Nave[id] * dni[ip] / nn2;
dnormal[i][id][0][ip] = -dnormdri[i][id][ip];
}
}
}
}
//############################ For the edge atoms of TMD ################################
else if (cont < Nnei) {
else if (cont > 1 && cont < Nnei) {
if (strcmp(elements[itype], "Mo") == 0 || strcmp(elements[itype], "W") == 0 ||
strcmp(elements[itype], "S") == 0 || strcmp(elements[itype], "Se") == 0) {
// derivatives of Ni[l] respect to the cont neighbors
@ -557,8 +603,7 @@ void PairILPTMD::calc_normal()
for (m = 0; m < cont; m++) {
for (id = 0; id < 3; id++) {
dnn[m][id] = (Nave[0] * dNave[0][m][id] + Nave[1] * dNave[1][m][id] +
Nave[2] * dNave[2][m][id]) /
nn;
Nave[2] * dNave[2][m][id]) / nn;
}
}
// dnormal[i][id][m][ip]: the derivative of normal[i][id] respect to r[m][ip], id,ip=0,1,2.
@ -589,6 +634,93 @@ void PairILPTMD::calc_normal()
}
}
} // for TMD
//############################ For Oxygen in the water molecule #######################
else if (strcmp(elements[itype], "Ow") == 0) {
if(cont == 0) {
jH1 = atom->map(tag[i] + 1);
jH2 = atom->map(tag[i] + 2);
iH1 = map[type[jH1]];
iH2 = map[type[jH2]];
if (strcmp(elements[iH1], "Hw") == 0 && strcmp(elements[iH2], "Hw") == 0) {
vect[0][0] = x[jH1][0] - xtp;
vect[0][1] = x[jH1][1] - ytp;
vect[0][2] = x[jH1][2] - ztp;
vect[1][0] = x[jH2][0] - xtp;
vect[1][1] = x[jH2][1] - ytp;
vect[1][2] = x[jH2][2] - ztp;
cont = 2;
} else {
error->one(FLERR, "The order of atoms in water molecule should be O H H !");
}
}
if (cont == 2) {
Nave[0] = (vect[0][0] + vect[1][0])/cont;
Nave[1] = (vect[0][1] + vect[1][1])/cont;
Nave[2] = (vect[0][2] + vect[1][2])/cont;
// the magnitude of the normal vector
nn2 = Nave[0] * Nave[0] + Nave[1] * Nave[1] + Nave[2] * Nave[2];
nn = sqrt(nn2);
if (nn == 0) error->one(FLERR, "The magnitude of the normal vector is zero");
// the unit normal vector
normal[i][0] = Nave[0] / nn;
normal[i][1] = Nave[1] / nn;
normal[i][2] = Nave[2] / nn;
// derivatives of non-normalized normal vector, dNave:3xcontx3 array
// dNave[id][m][ip]: the derivatve of the id component of Nave
// respect to the ip component of atom m
for (id = 0; id < 3; id++) {
for (ip = 0; ip < 3; ip++) {
for (m = 0; m < cont; m++) {
if (ip == id) { dNave[id][m][ip] = 0.5;}
else {dNave[id][m][ip] = 0.0;}
}
}
}
// derivatives of nn, dnn:contx3 vector
// dnn[m][id]: the derivative of nn respect to r[m][id], m=0,...Nnei-1; id=0,1,2
// r[m][id]: the id's component of atom m
for (m = 0; m < cont; m++) {
for (id = 0; id < 3; id++) {
dnn[m][id] = (Nave[0] * dNave[0][m][id] + Nave[1] * dNave[1][m][id] +
Nave[2] * dNave[2][m][id]) / nn;
}
}
// dnormal[i][id][m][ip]: the derivative of normal[i][id] respect to r[m][ip], id,ip=0,1,2.
// for atom m, which is a neighbor atom of atom i, m = 0,...,Nnei-1
for (m = 0; m < cont; m++) {
for (id = 0; id < 3; id++) {
for (ip = 0; ip < 3; ip++) {
dnormal[i][id][m][ip] = dNave[id][m][ip] / nn - Nave[id] * dnn[m][ip] / nn2;
}
}
}
// Calculte dNave/dri, defined as dpvdri
for (id = 0; id < 3; id++) {
for (ip = 0; ip < 3; ip++) {
dpvdri[id][ip] = 0.0;
for (k = 0; k < cont; k++) { dpvdri[id][ip] -= dNave[id][k][ip]; }
}
}
// derivatives of nn, dnn:3x1 vector
dni[0] = (Nave[0] * dpvdri[0][0] + Nave[1] * dpvdri[1][0] + Nave[2] * dpvdri[2][0]) / nn;
dni[1] = (Nave[0] * dpvdri[0][1] + Nave[1] * dpvdri[1][1] + Nave[2] * dpvdri[2][1]) / nn;
dni[2] = (Nave[0] * dpvdri[0][2] + Nave[1] * dpvdri[1][2] + Nave[2] * dpvdri[2][2]) / nn;
// derivatives of unit vector ni respect to ri, the result is 3x3 matrix
for (id = 0; id < 3; id++) {
for (ip = 0; ip < 3; ip++) {
dnormdri[i][id][ip] = dpvdri[id][ip] / nn - Nave[id] * dni[ip] / nn2;
}
}
}
else if (cont >= 3) {
error->one(FLERR,
"There are too many neighbors for calculating normals of water molecules");
}
}
//############################ For the edge & bulk atoms of GrhBN ################################
else {
if (cont == 2) {
@ -657,8 +789,7 @@ void PairILPTMD::calc_normal()
for (m = 0; m < cont; m++) {
for (id = 0; id < 3; id++) {
dnn[m][id] = (Nave[0] * dNave[0][m][id] + Nave[1] * dNave[1][m][id] +
Nave[2] * dNave[2][m][id]) /
nn;
Nave[2] * dNave[2][m][id]) / nn;
}
}
// dnormal[i][id][m][ip]: the derivative of normal[i][id] respect to r[m][ip], id,ip=0,1,2.

View File

@ -12,7 +12,7 @@
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Wengen Ouyang (Tel Aviv University)
Contributing author: Wengen Ouyang (Wuhan University)
e-mail: w.g.ouyang at gmail dot com
based on previous versions by Jaap Kroes

View File

@ -15,7 +15,7 @@
e-mail: w.g.ouyang at gmail dot com
This is a full version of the potential described in
[Ouyang et al, J. Chem. Theory Comput. 17, 7215-7223 (2021)]
[Ouyang et al., J. Chem. Theory Comput. 17, 7215-7223 (2021)]
------------------------------------------------------------------------- */
#include "pair_saip_metal.h"

View File

@ -14,7 +14,7 @@
/* ----------------------------------------------------------------------
Contributing author: Aidan Thompson (SNL) - original Tersoff implementation
Wengen Ouyang (TAU) - Shift addition
Wengen Ouyang (WHU) - Shift addition
------------------------------------------------------------------------- */
#include "pair_tersoff.h"

View File

@ -0,0 +1,63 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://www.lammps.org/, Sandia National Laboratories
LAMMPS development team: developers@lammps.org
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
This is an optimized version of aip/water/2dm based on the contribution of:
author: Wengen Ouyang (Wuhan University)
e-mail: w.g.ouyang at gmail dot com
Optimizations are done by:
author1: Xiaohui Duan (National Supercomputing Center in Wuxi, China)
e-mail: sunrise_duan at 126 dot com
author2: Ping Gao (National Supercomputing Center in Wuxi, China)
e-mail: qdgaoping at gmail dot com
Optimizations are described in:
Gao, Ping and Duan, Xiaohui, et al.:
LMFF: Efficient and Scalable Layered Materials Force Field on Heterogeneous Many-Core Processors
DOI: 10.1145/3458817.3476137
Potential is described by:
[Feng and Ouyang et al, J. Phys. Chem. C 127, 8704-8713 (2023).]
*/
#include "pair_aip_water_2dm_opt.h"
#include "atom.h"
#include "memory.h"
#include <cstring>
using namespace LAMMPS_NS;
PairAIPWATER2DMOpt::PairAIPWATER2DMOpt(LAMMPS *lmp) :
PairILPGrapheneHBN(lmp), PairILPTMD(lmp), PairAIPWATER2DM(lmp), PairILPGrapheneHBNOpt(lmp)
{
}
void PairAIPWATER2DMOpt::coeff(int narg, char **args)
{
PairILPTMD::coeff(narg, args);
memory->create(special_type, atom->ntypes + 1, "PairAIPWATER2DMOpt:check_sublayer");
for (int i = 1; i <= atom->ntypes; i++) {
int itype = map[i];
if (strcmp(elements[itype], "Mo") == 0 || strcmp(elements[itype], "W") == 0 ||
strcmp(elements[itype], "S") == 0 || strcmp(elements[itype], "Se") == 0 ||
strcmp(elements[itype], "Te") == 0) {
special_type[i] = TMD_METAL;
} else if (strcmp(elements[itype], "Hw") == 0 || strcmp(elements[itype], "Ow") == 0) {
special_type[i] = WATER;
} else {
special_type[i] = NOT_SPECIAL;
}
}
}

View File

@ -0,0 +1,38 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://www.lammps.org/, Sandia National Laboratories
LAMMPS development team: developers@lammps.org
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifdef PAIR_CLASS
// clang-format off
PairStyle(aip/water/2dm/opt,PairAIPWATER2DMOpt);
// clang-format on
#else
#ifndef LMP_PAIR_AIP_WATER_2DM_OPT_H
#define LMP_PAIR_AIP_WATER_2DM_OPT_H
#include "pair_ilp_graphene_hbn_opt.h"
#include "pair_aip_water_2dm.h"
namespace LAMMPS_NS {
class PairAIPWATER2DMOpt : public PairAIPWATER2DM, public PairILPGrapheneHBNOpt {
public:
PairAIPWATER2DMOpt(class LAMMPS *);
void coeff(int narg, char **args) override;
protected:
};
} // namespace LAMMPS_NS
#endif
#endif

View File

@ -11,12 +11,24 @@
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Optimization author1: Ping Gao (National Supercomputing Center in Wuxi, China)
This is an optimized version of ilp/graphene/hbn based on the contirubtion of:
author: Wengen Ouyang (Wuhan University, China)
e-mail: w.g.ouyang at gmail dot com
Optimizations are done by:
author1: Ping Gao (National Supercomputing Center in Wuxi, China) implements the base ILP potential.
e-mail: qdgaoping at gmail dot com
Optimization author2: Xiaohui Duan (National Supercomputing Center in Wuxi, China)
author2: Xiaohui Duan (Shandong University, China) adjusts the framework to adopt SAIP, TMD, WATER2DM, etc.
e-mail: sunrise_duan at 126 dot com
Provides some bugfixes and performance optimizations in this potential.
Optimizations are described in:
Gao, Ping and Duan, Xiaohui, et al:
LMFF: Efficient and Scalable Layered Materials Force Field on Heterogeneous Many-Core Processors
DOI: 10.1145/3458817.3476137
Potential is described by:
[Ouyang et al., J. Chem. Theory Comput. 16(1), 666-676 (2020)]
*/
#include "pair_ilp_graphene_hbn_opt.h"
@ -28,19 +40,21 @@
#include "interlayer_taper.h"
#include "memory.h"
#include "neigh_list.h"
#include "neigh_request.h"
#include "neighbor.h"
#include "pointers.h"
#include <cmath>
#include <utility>
#include <cstring>
using namespace LAMMPS_NS;
using namespace InterLayer;
static const char cite_ilp_cur[] =
"ilp/graphene/hbn/opt potential: doi:10.1145/3458817.3476137\n"
"ilp/graphene/hbn/opt potential doi:10.1145/3458817.3476137\n"
"@inproceedings{gao2021lmff\n"
" author = {Gao, Ping and Duan, Xiaohui and others},\n"
" title = {{LMFF}: Efficient and Scalable Layered Materials Force Field on Heterogeneous "
" author = {Gao, Ping and Duan, Xiaohui and Others},\n"
" title = {LMFF: Efficient and Scalable Layered Materials Force Field on Heterogeneous "
"Many-Core Processors},\n"
" year = {2021},\n"
" isbn = {9781450384421},\n"
@ -50,14 +64,12 @@ static const char cite_ilp_cur[] =
" doi = {10.1145/3458817.3476137},\n"
" booktitle = {Proceedings of the International Conference for High Performance Computing, "
"Networking, Storage and Analysis},\n"
" pages = {42},\n"
" articleno = {42},\n"
" numpages = {14},\n"
" location = {St.~Louis, Missouri},\n"
" location = {St. Louis, Missouri},\n"
" series = {SC'21},\n"
"}\n\n";
static bool check_vdw(tagint itag, tagint jtag, double *xi, double *xj);
/* ---------------------------------------------------------------------- */
PairILPGrapheneHBNOpt::PairILPGrapheneHBNOpt(LAMMPS *lmp) :
@ -168,6 +180,36 @@ void PairILPGrapheneHBNOpt::compute(int eflag, int vflag)
}
}
}
} else if (variant == AIP_WATER_2DM) {
if (eflag_global || eflag_atom) {
if (vflag_either) {
if (tap_flag) {
eval<6, 1, 1, 1, AIP_WATER_2DM>();
} else {
eval<6, 1, 1, 0, AIP_WATER_2DM>();
}
} else {
if (tap_flag) {
eval<6, 1, 0, 1, AIP_WATER_2DM>();
} else {
eval<6, 1, 0, 0, AIP_WATER_2DM>();
}
}
} else {
if (vflag_either) {
if (tap_flag) {
eval<6, 0, 1, 1, AIP_WATER_2DM>();
} else {
eval<6, 0, 1, 0, AIP_WATER_2DM>();
}
} else {
if (tap_flag) {
eval<6, 0, 0, 1, AIP_WATER_2DM>();
} else {
eval<6, 0, 0, 0, AIP_WATER_2DM>();
}
}
}
} else if (variant == SAIP_METAL) {
if (eflag_global || eflag_atom) {
if (vflag_either) {
@ -255,7 +297,7 @@ void PairILPGrapheneHBNOpt::eval()
rsq = delx * delx + dely * dely + delz * delz;
if (rsq != 0 && rsq < cutILPsq[itype_map][jtype]) {
if (VARIANT == ILP_TMD && special_type[itype] && itype != type[j]) continue;
if ((VARIANT == ILP_TMD || VARIANT == AIP_WATER_2DM) && special_type[itype] == TMD_METAL && itype != type[j]) continue;
if (ILP_nneigh >= MAX_NNEIGH) {
error->one(FLERR, "There are too many neighbors for calculating normals");
}
@ -269,7 +311,8 @@ void PairILPGrapheneHBNOpt::eval()
dproddni[2] = 0.0;
double norm[3], dnormdxi[3][3], dnormdxk[MAX_NNEIGH][3][3];
calc_normal<MAX_NNEIGH>(i, ILP_neigh, ILP_nneigh, norm, dnormdxi, dnormdxk);
calc_atom_normal<MAX_NNEIGH>(i, itype, ILP_neigh, ILP_nneigh, norm, dnormdxi, dnormdxk);
for (jj = 0; jj < jnum_inter; jj++) {
j = jlist_inter[jj];
@ -298,7 +341,7 @@ void PairILPGrapheneHBNOpt::eval()
Tap = 1.0;
dTap = 0.0;
}
if (VARIANT != SAIP_METAL || !special_type[itype]) {
if (VARIANT != SAIP_METAL || special_type[itype] != SAIP_BNCH) {
// Calculate the transverse distance
prodnorm1 = norm[0] * delx + norm[1] * dely + norm[2] * delz;
rhosq1 = rsq - prodnorm1 * prodnorm1; // rho_ij
@ -310,7 +353,7 @@ void PairILPGrapheneHBNOpt::eval()
frho1 = exp1 * p.C;
Erep = 0.5 * p.epsilon + frho1;
if (VARIANT == SAIP_METAL && special_type[jtype]) { Erep += 0.5 * p.epsilon + p.C; }
if (VARIANT == SAIP_METAL && special_type[jtype] == SAIP_BNCH) { Erep += 0.5 * p.epsilon + p.C; }
Vilp = exp0 * Erep;
// derivatives
@ -428,6 +471,18 @@ inline void deriv_normal(double dndr[3][3], double *del, double *n, double rnnor
dndr[1][2] = (del[1] * n[0] * n[1] + del[0] * (n[0] * n[0] + n[2] * n[2])) * rnnorm;
dndr[2][2] = (del[1] * n[0] * n[2] - del[0] * n[1] * n[2]) * rnnorm;
}
inline void deriv_hat(double dnhatdn[3][3], double *n, double rnnorm, double factor){
double cfactor = rnnorm * factor;
dnhatdn[0][0] = (n[1]*n[1]+n[2]*n[2])*cfactor;
dnhatdn[1][0] = -n[1]*n[0]*cfactor;
dnhatdn[2][0] = -n[2]*n[0]*cfactor;
dnhatdn[0][1] = -n[0]*n[1]*cfactor;
dnhatdn[1][1] = (n[0]*n[0]+n[2]*n[2])*cfactor;
dnhatdn[2][1] = -n[2]*n[1]*cfactor;
dnhatdn[0][2] = -n[0]*n[2]*cfactor;
dnhatdn[1][2] = -n[1]*n[2]*cfactor;
dnhatdn[2][2] = (n[0]*n[0]+n[1]*n[1])*cfactor;
}
inline double normalize_factor(double *n)
{
double nnorm = sqrt(n[0] * n[0] + n[1] * n[1] + n[2] * n[2]);
@ -441,7 +496,7 @@ inline double normalize_factor(double *n)
Yet another normal calculation method for simpiler code.
*/
template <int MAX_NNEIGH>
void PairILPGrapheneHBNOpt::calc_normal(int i, int *ILP_neigh, int nneigh, double *n,
void PairILPGrapheneHBNOpt::calc_atom_normal(int i, int itype, int *ILP_neigh, int nneigh, double *n,
double (*dnormdri)[3], double (*dnormdrk)[3][3])
{
double **x = atom->x;
@ -475,6 +530,32 @@ void PairILPGrapheneHBNOpt::calc_normal(int i, int *ILP_neigh, int nneigh, doubl
vet[jj][2] = x[j][2] - x[i][2];
}
//specialize for AIP_WATER_2DM for hydrogen has special normal vector rule
if (variant == AIP_WATER_2DM && special_type[itype] == WATER) {
if (nneigh == 1){
n[0] = vet[0][0];
n[1] = vet[0][1];
n[2] = vet[0][2];
double rnnorm = normalize_factor(n);
deriv_hat(dnormdri, n, rnnorm, -1.0);
deriv_hat(dnormdrk[0], n, rnnorm, 1.0);
} else if (nneigh == 2){
n[0] = (vet[0][0] + vet[1][0])*0.5;
n[1] = (vet[0][1] + vet[1][1])*0.5;
n[2] = (vet[0][2] + vet[1][2])*0.5;
double rnnorm = normalize_factor(n);
deriv_hat(dnormdri, n, rnnorm, -1.0);
deriv_hat(dnormdrk[0], n, rnnorm, 0.5);
deriv_hat(dnormdrk[1], n, rnnorm, 0.5);
} else {
error->one(FLERR, "malformed water");
}
return;
}
if (nneigh <= 1) {
n[0] = 0.0;
n[1] = 0.0;

View File

@ -35,7 +35,7 @@ class PairILPGrapheneHBNOpt : virtual public PairILPGrapheneHBN {
protected:
void update_internal_list();
template <int MAX_NNEIGH>
void calc_normal(int i, int *ILP_neigh, int nneigh, double *normal, double (*dnormdri)[3],
void calc_atom_normal(int i, int itype, int *ILP_neigh, int nneigh, double *normal, double (*dnormdri)[3],
double (*dnormdrk)[3][3]);
template <int MAX_NNEIGH, int EFLAG, int VFLAG_EITHER, int TAP_FLAG, int VARIANT = ILP_GrhBN>
void eval();
@ -44,6 +44,14 @@ class PairILPGrapheneHBNOpt : virtual public PairILPGrapheneHBN {
int *special_type;
int *num_intra, *num_inter, *num_vdw;
int inum_max, jnum_max;
enum special_type_const {
NOT_SPECIAL = 0,
TMD_METAL,
SAIP_BNCH,
WATER,
};
};
} // namespace LAMMPS_NS

View File

@ -23,7 +23,7 @@
e-mail: qdgaoping at gmail dot com
Optimizations are described in:
Gao, Ping and Duan, Xiaohui, et al:
Gao, Ping and Duan, Xiaohui, et al.:
LMFF: Efficient and Scalable Layered Materials Force Field on Heterogeneous Many-Core Processors
DOI: 10.1145/3458817.3476137