forked from lijiext/lammps
Merge pull request #511 from akohlmey/add-compute-cnp
Integrate compute cnp/atom contributed by Paulo Branicio (USC)
This commit is contained in:
sjplimp
GitHub
commit
9a7207e34c
Binary file not shown.
|
After Width: | Height: | Size: 13 KiB |
|
|
@ -0,0 +1,14 @@
|
|||
\documentclass[12pt,article]{article}
|
||||
|
||||
\usepackage{indentfirst}
|
||||
\usepackage{amsmath}
|
||||
|
||||
\begin{document}
|
||||
|
||||
\begin{eqnarray*}
|
||||
r_{c}^{fcc} & = & \frac{1}{2} \left(\frac{\sqrt{2}}{2} + 1\right) \mathrm{a} \simeq 0.8536 \:\mathrm{a} \\
|
||||
r_{c}^{bcc} & = & \frac{1}{2}(\sqrt{2} + 1) \mathrm{a} \simeq 1.207 \:\mathrm{a} \\
|
||||
r_{c}^{hcp} & = & \frac{1}{2}\left(1+\sqrt{\frac{4+2x^{2}}{3}}\right) \mathrm{a}
|
||||
\end{eqnarray*}
|
||||
|
||||
\end{document}
|
||||
Binary file not shown.
|
After Width: | Height: | Size: 2.5 KiB |
|
|
@ -0,0 +1,12 @@
|
|||
\documentclass[12pt,article]{article}
|
||||
|
||||
\usepackage{indentfirst}
|
||||
\usepackage{amsmath}
|
||||
|
||||
\begin{document}
|
||||
|
||||
$$
|
||||
Rc + Rs > 2*{\rm cutoff}
|
||||
$$
|
||||
|
||||
\end{document}
|
||||
Binary file not shown.
|
After Width: | Height: | Size: 23 KiB |
|
|
@ -0,0 +1,9 @@
|
|||
\documentclass[12pt]{article}
|
||||
|
||||
\begin{document}
|
||||
|
||||
$$
|
||||
Q_{i} = \frac{1}{n_i}\sum_{j = 1}^{n_i} | \sum_{k = 1}^{n_{ij}} \vec{R}_{ik} + \vec{R}_{jk} |^2
|
||||
$$
|
||||
|
||||
\end{document}
|
||||
|
|
@ -831,6 +831,7 @@ package"_Section_start.html#start_3.
|
|||
|
||||
"ackland/atom"_compute_ackland_atom.html,
|
||||
"basal/atom"_compute_basal_atom.html,
|
||||
"cnp/atom"_compute_cnp_atom.html,
|
||||
"dpd"_compute_dpd.html,
|
||||
"dpd/atom"_compute_dpd_atom.html,
|
||||
"fep"_compute_fep.html,
|
||||
|
|
|
|||
|
|
@ -26,7 +26,7 @@ Define a computation that calculates the CNA (Common Neighbor
|
|||
Analysis) pattern for each atom in the group. In solid-state systems
|
||||
the CNA pattern is a useful measure of the local crystal structure
|
||||
around an atom. The CNA methodology is described in "(Faken)"_#Faken
|
||||
and "(Tsuzuki)"_#Tsuzuki.
|
||||
and "(Tsuzuki)"_#Tsuzuki1.
|
||||
|
||||
Currently, there are five kinds of CNA patterns LAMMPS recognizes:
|
||||
|
||||
|
|
@ -93,5 +93,5 @@ above.
|
|||
:link(Faken)
|
||||
[(Faken)] Faken, Jonsson, Comput Mater Sci, 2, 279 (1994).
|
||||
|
||||
:link(Tsuzuki)
|
||||
:link(Tsuzuki1)
|
||||
[(Tsuzuki)] Tsuzuki, Branicio, Rino, Comput Phys Comm, 177, 518 (2007).
|
||||
|
|
|
|||
|
|
@ -0,0 +1,111 @@
|
|||
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
|
||||
|
||||
:link(lws,http://lammps.sandia.gov)
|
||||
:link(ld,Manual.html)
|
||||
:link(lc,Section_commands.html#comm)
|
||||
|
||||
:line
|
||||
|
||||
compute cnp/atom command :h3
|
||||
|
||||
[Syntax:]
|
||||
|
||||
compute ID group-ID cnp/atom cutoff :pre
|
||||
|
||||
ID, group-ID are documented in "compute"_compute.html command
|
||||
cnp/atom = style name of this compute command
|
||||
cutoff = cutoff distance for nearest neighbors (distance units) :ul
|
||||
|
||||
[Examples:]
|
||||
|
||||
compute 1 all cnp/atom 3.08 :pre
|
||||
|
||||
[Description:]
|
||||
|
||||
Define a computation that calculates the Common Neighborhood
|
||||
Parameter (CNP) for each atom in the group. In solid-state systems
|
||||
the CNP is a useful measure of the local crystal structure
|
||||
around an atom and can be used to characterize whether the
|
||||
atom is part of a perfect lattice, a local defect (e.g. a dislocation
|
||||
or stacking fault), or at a surface.
|
||||
|
||||
The value of the CNP parameter will be 0.0 for atoms not in the
|
||||
specified compute group. Note that normally a CNP calculation should
|
||||
only be performed on single component systems.
|
||||
|
||||
This parameter is computed using the following formula from
|
||||
"(Tsuzuki)"_#Tsuzuki2
|
||||
|
||||
:c,image(Eqs/cnp_eq.jpg)
|
||||
|
||||
where the index {j} goes over the {n}i nearest neighbors of atom
|
||||
{i}, and the index {k} goes over the {n}ij common nearest neighbors
|
||||
between atom {i} and atom {j}. Rik and Rjk are the vectors connecting atom
|
||||
{k} to atoms {i} and {j}. The quantity in the double sum is computed
|
||||
for each atom.
|
||||
|
||||
The CNP calculation is sensitive to the specified cutoff value.
|
||||
You should ensure that the appropriate nearest neighbors of an atom are
|
||||
found within the cutoff distance for the presumed crystal structure.
|
||||
E.g. 12 nearest neighbor for perfect FCC and HCP crystals, 14 nearest
|
||||
neighbors for perfect BCC crystals. These formulas can be used to
|
||||
obtain a good cutoff distance:
|
||||
|
||||
:c,image(Eqs/cnp_cutoff.jpg)
|
||||
|
||||
where a is the lattice constant for the crystal structure concerned
|
||||
and in the HCP case, x = (c/a) / 1.633, where 1.633 is the ideal c/a
|
||||
for HCP crystals.
|
||||
|
||||
Also note that since the CNP calculation in LAMMPS uses the neighbors
|
||||
of an owned atom to find the nearest neighbors of a ghost atom, the
|
||||
following relation should also be satisfied:
|
||||
|
||||
:c,image(Eqs/cnp_cutoff2.jpg)
|
||||
|
||||
where Rc is the cutoff distance of the potential, Rs is the skin
|
||||
distance as specified by the "neighbor"_neighbor.html command, and
|
||||
cutoff is the argument used with the compute cnp/atom command. LAMMPS
|
||||
will issue a warning if this is not the case.
|
||||
|
||||
The neighbor list needed to compute this quantity is constructed each
|
||||
time the calculation is performed (e.g. each time a snapshot of atoms
|
||||
is dumped). Thus it can be inefficient to compute/dump this quantity
|
||||
too frequently or to have multiple compute/dump commands, each with a
|
||||
{cnp/atom} style.
|
||||
|
||||
[Output info:]
|
||||
|
||||
This compute calculates a per-atom vector, which can be accessed by
|
||||
any command that uses per-atom values from a compute as input. See
|
||||
"Section 6.15"_Section_howto.html#howto_15 for an overview of
|
||||
LAMMPS output options.
|
||||
|
||||
The per-atom vector values will be real positive numbers. Some typical CNP
|
||||
values:
|
||||
|
||||
FCC lattice = 0.0
|
||||
BCC lattice = 0.0
|
||||
HCP lattice = 4.4 :pre
|
||||
|
||||
FCC (111) surface ~ 13.0
|
||||
FCC (100) surface ~ 26.5
|
||||
FCC dislocation core ~ 11 :pre
|
||||
|
||||
[Restrictions:]
|
||||
|
||||
This compute is part of the USER-MISC package. It is only enabled if
|
||||
LAMMPS was built with that package. See the "Making
|
||||
LAMMPS"_Section_start.html#start_3 section for more info.
|
||||
|
||||
[Related commands:]
|
||||
|
||||
"compute cna/atom"_compute_cna_atom.html
|
||||
"compute centro/atom"_compute_centro_atom.html
|
||||
|
||||
[Default:] none
|
||||
|
||||
:line
|
||||
|
||||
:link(Tsuzuki2)
|
||||
[(Tsuzuki)] Tsuzuki, Branicio, Rino, Comput Phys Comm, 177, 518 (2007).
|
||||
|
|
@ -17,6 +17,7 @@ Computes :h1
|
|||
compute_chunk_atom
|
||||
compute_cluster_atom
|
||||
compute_cna_atom
|
||||
compute_cnp_atom
|
||||
compute_com
|
||||
compute_com_chunk
|
||||
compute_contact_atom
|
||||
|
|
|
|||
|
|
@ -301,6 +301,7 @@ compute_centro_atom.html
|
|||
compute_chunk_atom.html
|
||||
compute_cluster_atom.html
|
||||
compute_cna_atom.html
|
||||
compute_cnp_atom.html
|
||||
compute_com.html
|
||||
compute_com_chunk.html
|
||||
compute_contact_atom.html
|
||||
|
|
|
|||
File diff suppressed because it is too large
Load Diff
|
|
@ -0,0 +1,51 @@
|
|||
# Generation and relaxation of a partial dislocation in Cu perfect FCC crystal
|
||||
|
||||
# Initialization
|
||||
units metal
|
||||
boundary p p p
|
||||
atom_style atomic
|
||||
|
||||
# create simulation box and system
|
||||
lattice fcc 3.615 origin 0.01 0.01 0.01 orient x -1 -1 2 orient y 1 1 1 orient z -1 1 0
|
||||
region mdbox block 0 3 0.0 14.0 0 84 units lattice
|
||||
region system block 0 3 1.1 13.1 0 84 units lattice
|
||||
create_box 2 mdbox
|
||||
create_atoms 1 region system
|
||||
|
||||
# Define atoms mass and force field
|
||||
mass * 63.54
|
||||
pair_style eam/alloy
|
||||
pair_coeff * * Cu_Mishin1.eam Cu Cu
|
||||
|
||||
# Delete a plane of atoms along the z direction to generate a partial dislocation
|
||||
region dislocation_atoms block 0 3 7 14 41.9 42.1 units lattice
|
||||
delete_atoms region dislocation_atoms
|
||||
region quarter_up block 0 3 7 11 0 84 units lattice
|
||||
group middle region quarter_up
|
||||
|
||||
# specify simulation parameters
|
||||
timestep 0.004
|
||||
|
||||
# Relax configuration using conjugate gradient
|
||||
#min_style cg
|
||||
#minimize 1.0e-4 1.0e-6 100 1000
|
||||
|
||||
# Setup calculations
|
||||
compute 1 all cnp/atom 3.086
|
||||
compute 2 all cna/atom 3.086
|
||||
compute 3 all centro/atom fcc
|
||||
compute 4 all coord/atom cutoff 3.086
|
||||
dump 1 all custom 100 dump.lammpstrj id type xu yu zu c_1 c_2 c_3 c_4
|
||||
|
||||
### Set up thermo display
|
||||
thermo 10
|
||||
thermo_style custom step atoms temp press pe ke etotal
|
||||
|
||||
# Relax the system performing a langevin dynamics (freeze motion along y 111 direction)
|
||||
fix 1 all nve
|
||||
fix 2 all langevin 50 1 0.1 699483
|
||||
fix 3 all setforce NULL 0.0 NULL
|
||||
fix 4 middle setforce 0.0 0.0 0.0
|
||||
run 100
|
||||
unfix 4
|
||||
run 200
|
||||
|
|
@ -0,0 +1,185 @@
|
|||
LAMMPS (19 May 2017)
|
||||
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:90)
|
||||
using 1 OpenMP thread(s) per MPI task
|
||||
# Generation and relaxation of a partial dislocation in Cu perfect FCC crystal
|
||||
|
||||
# Initialization
|
||||
units metal
|
||||
boundary p p p
|
||||
atom_style atomic
|
||||
|
||||
# create simulation box and system
|
||||
lattice fcc 3.615 origin 0.01 0.01 0.01 orient x -1 -1 2 orient y 1 1 1 orient z -1 1 0
|
||||
Lattice spacing in x,y,z = 5.90327 6.26136 5.11238
|
||||
region mdbox block 0 3 0.0 14.0 0 84 units lattice
|
||||
region system block 0 3 1.1 13.1 0 84 units lattice
|
||||
create_box 2 mdbox
|
||||
Created orthogonal box = (0 0 0) to (17.7098 87.6591 429.44)
|
||||
1 by 1 by 4 MPI processor grid
|
||||
create_atoms 1 region system
|
||||
Created 48384 atoms
|
||||
|
||||
# Define atoms mass and force field
|
||||
mass * 63.54
|
||||
pair_style eam/alloy
|
||||
pair_coeff * * Cu_Mishin1.eam Cu Cu
|
||||
|
||||
# Delete a plane of atoms along the z direction to generate a partial dislocation
|
||||
region dislocation_atoms block 0 3 7 14 41.9 42.1 units lattice
|
||||
delete_atoms region dislocation_atoms
|
||||
Deleted 76 atoms, new total = 48308
|
||||
region quarter_up block 0 3 7 11 0 84 units lattice
|
||||
group middle region quarter_up
|
||||
16080 atoms in group middle
|
||||
|
||||
# specify simulation parameters
|
||||
timestep 0.004
|
||||
|
||||
# Relax configuration using conjugate gradient
|
||||
#min_style cg
|
||||
#minimize 1.0e-4 1.0e-6 100 1000
|
||||
|
||||
# Setup calculations
|
||||
compute 1 all cnp/atom 3.086
|
||||
compute 2 all cna/atom 3.086
|
||||
compute 3 all centro/atom fcc
|
||||
compute 4 all coord/atom cutoff 3.086
|
||||
dump 1 all custom 100 dump.lammpstrj id type xu yu zu c_1 c_2 c_3 c_4
|
||||
|
||||
### Set up thermo display
|
||||
thermo 10
|
||||
thermo_style custom step atoms temp press pe ke etotal
|
||||
|
||||
# Relax the system performing a langevin dynamics (freeze motion along y 111 direction)
|
||||
fix 1 all nve
|
||||
fix 2 all langevin 50 1 0.1 699483
|
||||
fix 3 all setforce NULL 0.0 NULL
|
||||
fix 4 middle setforce 0.0 0.0 0.0
|
||||
run 100
|
||||
Neighbor list info ...
|
||||
update every 1 steps, delay 10 steps, check yes
|
||||
max neighbors/atom: 2000, page size: 100000
|
||||
master list distance cutoff = 7.50679
|
||||
ghost atom cutoff = 7.50679
|
||||
binsize = 3.75339, bins = 5 24 115
|
||||
5 neighbor lists, perpetual/occasional/extra = 1 4 0
|
||||
(1) pair eam/alloy, perpetual
|
||||
attributes: half, newton on
|
||||
pair build: half/bin/atomonly/newton
|
||||
stencil: half/bin/3d/newton
|
||||
bin: standard
|
||||
(2) compute cnp/atom, occasional
|
||||
attributes: full, newton on
|
||||
pair build: full/bin/atomonly
|
||||
stencil: full/bin/3d
|
||||
bin: standard
|
||||
(3) compute cna/atom, occasional
|
||||
attributes: full, newton on
|
||||
pair build: full/bin/atomonly
|
||||
stencil: full/bin/3d
|
||||
bin: standard
|
||||
(4) compute centro/atom, occasional
|
||||
attributes: full, newton on
|
||||
pair build: full/bin/atomonly
|
||||
stencil: full/bin/3d
|
||||
bin: standard
|
||||
(5) compute coord/atom, occasional
|
||||
attributes: full, newton on
|
||||
pair build: full/bin/atomonly
|
||||
stencil: full/bin/3d
|
||||
bin: standard
|
||||
Per MPI rank memory allocation (min/avg/max) = 45.41 | 45.41 | 45.41 Mbytes
|
||||
Step Atoms Temp Press PotEng KinEng TotEng
|
||||
0 48308 0 -3388.0911 -169746.07 0 -169746.07
|
||||
10 48308 7.35092 -3091.0864 -169715.96 45.900393 -169670.05
|
||||
20 48308 9.9162268 -2822.7045 -169678.51 61.918604 -169616.59
|
||||
30 48308 12.351316 -2726.7195 -169666.35 77.123716 -169589.23
|
||||
40 48308 13.302856 -2703.586 -169662.9 83.06529 -169579.83
|
||||
50 48308 12.782228 -2706.8662 -169662.36 79.814401 -169582.55
|
||||
60 48308 12.198179 -2772.4206 -169670.02 76.167503 -169593.86
|
||||
70 48308 10.663322 -2841.3384 -169677.48 66.583595 -169610.9
|
||||
80 48308 9.1169804 -2932.3896 -169687.85 56.927974 -169630.92
|
||||
90 48308 7.2905076 -3029.9433 -169699.09 45.523167 -169653.56
|
||||
100 48308 5.4063635 -3139.4496 -169711.65 33.758252 -169677.89
|
||||
Loop time of 10.9003 on 4 procs for 100 steps with 48308 atoms
|
||||
|
||||
Performance: 3.171 ns/day, 7.570 hours/ns, 9.174 timesteps/s
|
||||
31.8% CPU use with 4 MPI tasks x 1 OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 9.8764 | 9.9587 | 10.021 | 1.6 | 91.36
|
||||
Neigh | 0 | 0 | 0 | 0.0 | 0.00
|
||||
Comm | 0.1232 | 0.18385 | 0.26683 | 12.1 | 1.69
|
||||
Output | 0.45385 | 0.45451 | 0.45634 | 0.2 | 4.17
|
||||
Modify | 0.25026 | 0.2537 | 0.25744 | 0.5 | 2.33
|
||||
Other | | 0.04949 | | | 0.45
|
||||
|
||||
Nlocal: 12077 ave 12096 max 12020 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 3
|
||||
Nghost: 14204 ave 14261 max 14109 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 0 2
|
||||
Neighs: 814050 ave 818584 max 809212 min
|
||||
Histogram: 1 0 0 0 0 2 0 0 0 1
|
||||
FullNghs: 1.6281e+06 ave 1.63296e+06 max 1.61808e+06 min
|
||||
Histogram: 1 0 0 0 0 0 1 0 0 2
|
||||
|
||||
Total # of neighbors = 6512400
|
||||
Ave neighs/atom = 134.81
|
||||
Neighbor list builds = 0
|
||||
Dangerous builds = 0
|
||||
unfix 4
|
||||
run 200
|
||||
Per MPI rank memory allocation (min/avg/max) = 45.41 | 45.41 | 45.41 Mbytes
|
||||
Step Atoms Temp Press PotEng KinEng TotEng
|
||||
100 48308 5.4063635 -3139.4496 -169711.65 33.758252 -169677.89
|
||||
110 48308 15.260795 -2793.119 -169677.24 95.290993 -169581.95
|
||||
120 48308 18.548656 -2433.1584 -169624.79 115.82096 -169508.97
|
||||
130 48308 22.15831 -2276.626 -169604.28 138.36025 -169465.92
|
||||
140 48308 24.393841 -2208.1771 -169596.16 152.31929 -169443.84
|
||||
150 48308 24.797558 -2173.3145 -169591.43 154.84016 -169436.59
|
||||
160 48308 24.73371 -2188.909 -169593.08 154.44148 -169438.64
|
||||
170 48308 24.128467 -2220.3404 -169596.96 150.66225 -169446.29
|
||||
180 48308 22.975708 -2275.1244 -169602.72 143.46422 -169459.26
|
||||
190 48308 21.936324 -2348.3762 -169610.59 136.97413 -169473.61
|
||||
200 48308 20.516249 -2432.8447 -169619.98 128.10694 -169491.87
|
||||
210 48308 19.000566 -2510.2915 -169628.58 118.64276 -169509.93
|
||||
220 48308 17.490407 -2597.299 -169638.24 109.21307 -169529.03
|
||||
230 48308 16.062482 -2684.1203 -169648.31 100.29687 -169548.01
|
||||
240 48308 14.360342 -2768.2313 -169657.7 89.668411 -169568.03
|
||||
250 48308 12.802315 -2852.6965 -169666.99 79.939831 -169587.05
|
||||
260 48308 11.258205 -2944.4533 -169677.52 70.298142 -169607.23
|
||||
270 48308 9.6159129 -3038.6304 -169688.06 60.043393 -169628.02
|
||||
280 48308 7.972425 -3129.0826 -169698.03 49.781176 -169648.25
|
||||
290 48308 6.3752377 -3219.2054 -169708.23 39.808067 -169668.42
|
||||
300 48308 4.7374688 -3306.1468 -169718.27 29.58156 -169688.69
|
||||
Loop time of 23.0164 on 4 procs for 200 steps with 48308 atoms
|
||||
|
||||
Performance: 3.003 ns/day, 7.992 hours/ns, 8.689 timesteps/s
|
||||
31.8% CPU use with 4 MPI tasks x 1 OpenMP threads
|
||||
|
||||
MPI task timing breakdown:
|
||||
Section | min time | avg time | max time |%varavg| %total
|
||||
---------------------------------------------------------------
|
||||
Pair | 20.221 | 20.423 | 20.57 | 3.1 | 88.73
|
||||
Neigh | 0 | 0 | 0 | 0.0 | 0.00
|
||||
Comm | 0.27748 | 0.42603 | 0.62832 | 21.4 | 1.85
|
||||
Output | 1.5454 | 1.5473 | 1.5529 | 0.3 | 6.72
|
||||
Modify | 0.48886 | 0.49773 | 0.50842 | 1.1 | 2.16
|
||||
Other | | 0.1221 | | | 0.53
|
||||
|
||||
Nlocal: 12077 ave 12096 max 12020 min
|
||||
Histogram: 1 0 0 0 0 0 0 0 0 3
|
||||
Nghost: 14204 ave 14261 max 14109 min
|
||||
Histogram: 1 0 0 0 0 1 0 0 0 2
|
||||
Neighs: 814094 ave 818584 max 809212 min
|
||||
Histogram: 1 0 0 0 0 2 0 0 0 1
|
||||
FullNghs: 1.62852e+06 ave 1.63296e+06 max 1.61892e+06 min
|
||||
Histogram: 1 0 0 0 0 0 0 1 0 2
|
||||
|
||||
Total # of neighbors = 6514094
|
||||
Ave neighs/atom = 134.845
|
||||
Neighbor list builds = 0
|
||||
Dangerous builds = 0
|
||||
Total wall time: 0:00:35
|
||||
|
|
@ -177,8 +177,8 @@
|
|||
/compute_basal_atom.h
|
||||
/compute_body_local.cpp
|
||||
/compute_body_local.h
|
||||
/compute_cna_atom2.cpp
|
||||
/compute_cna_atom2.h
|
||||
/compute_cnp_atom.cpp
|
||||
/compute_cnp_atom.h
|
||||
/compute_damage_atom.cpp
|
||||
/compute_damage_atom.h
|
||||
/compute_dilatation_atom.cpp
|
||||
|
|
|
|||
|
|
@ -28,6 +28,7 @@ bond_style harmonic/shift, Carsten Svaneborg, science at zqex.dk, 8 Aug 11
|
|||
bond_style harmonic/shift/cut, Carsten Svaneborg, science at zqex.dk, 8 Aug 11
|
||||
compute ackland/atom, Gerolf Ziegenhain, gerolf at ziegenhain.com, 4 Oct 2007
|
||||
compute basal/atom, Christopher Barrett, cdb333 at cavs.msstate.edu, 3 Mar 2013
|
||||
compute cnp/atom, Paulo Branicio (USC), branicio at usc.edu, 31 May 2017
|
||||
compute temp/rotate, Laurent Joly (U Lyon), ljoly.ulyon at gmail.com, 8 Aug 11
|
||||
compute PRESSURE/GREM, David Stelter, dstelter@bu.edu, 22 Nov 16
|
||||
dihedral_style cosine/shift/exp, Carsten Svaneborg, science at zqex.dk, 8 Aug 11
|
||||
|
|
|
|||
|
|
@ -0,0 +1,331 @@
|
|||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
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.
|
||||
|
||||
Common Neighbor Parameter as proposed in:
|
||||
Tsuzuki, Branicio, Rino, Comput Phys Comm, 177, 518 (2007)
|
||||
Cite: http://dx.doi.org/10.1063/1.2197987
|
||||
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Paulo Branicio (University of Southern California)
|
||||
branicio@usc.edu
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <string.h>
|
||||
#include <stdlib.h>
|
||||
#include <math.h>
|
||||
|
||||
#include "compute_cnp_atom.h"
|
||||
#include "atom.h"
|
||||
#include "update.h"
|
||||
#include "modify.h"
|
||||
#include "neighbor.h"
|
||||
#include "neigh_list.h"
|
||||
#include "neigh_request.h"
|
||||
#include "force.h"
|
||||
#include "pair.h"
|
||||
#include "comm.h"
|
||||
#include "memory.h"
|
||||
#include "error.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
//define maximum values
|
||||
#define MAXNEAR 24
|
||||
#define MAXCOMMON 12
|
||||
|
||||
enum{NCOMMON};
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
ComputeCNPAtom::ComputeCNPAtom(LAMMPS *lmp, int narg, char **arg) :
|
||||
Compute(lmp, narg, arg),
|
||||
nearest(NULL), nnearest(NULL), cnpv(NULL)
|
||||
{
|
||||
if (narg != 4) error->all(FLERR,"Illegal compute cnp/atom command");
|
||||
|
||||
peratom_flag = 1;
|
||||
size_peratom_cols = 0;
|
||||
|
||||
double cutoff = force->numeric(FLERR,arg[3]);
|
||||
if (cutoff < 0.0) error->all(FLERR,"Illegal compute cnp/atom command");
|
||||
cutsq = cutoff*cutoff;
|
||||
|
||||
// apply check for single type atoms in compute group
|
||||
int lasttype = -1;
|
||||
int n = -1;
|
||||
for (int i=0; i < atom->nlocal; ++i) {
|
||||
if (atom->mask[i] & groupbit) {
|
||||
if (lasttype != atom->type[i]) {
|
||||
lasttype = atom->type[i];
|
||||
++n;
|
||||
}
|
||||
}
|
||||
}
|
||||
int all_n = 0;
|
||||
MPI_Allreduce(&n,&all_n,1,MPI_INT,MPI_MAX,world);
|
||||
if (all_n > 0)
|
||||
error->warning(FLERR,"Compute cnp/atom requested on multi-type system");
|
||||
|
||||
nmax = 0;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
ComputeCNPAtom::~ComputeCNPAtom()
|
||||
{
|
||||
memory->destroy(nearest);
|
||||
memory->destroy(nnearest);
|
||||
memory->destroy(cnpv);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void ComputeCNPAtom::init()
|
||||
{
|
||||
if (force->pair == NULL)
|
||||
error->all(FLERR,"Compute cnp/atom requires a pair style be defined");
|
||||
|
||||
if (sqrt(cutsq) > force->pair->cutforce)
|
||||
error->all(FLERR,"Compute cnp/atom cutoff is longer than pairwise cutoff");
|
||||
|
||||
if (2.0*sqrt(cutsq) > force->pair->cutforce + neighbor->skin &&
|
||||
comm->me == 0)
|
||||
error->warning(FLERR,"Compute cnp/atom cutoff may be too large to find "
|
||||
"ghost atom neighbors");
|
||||
|
||||
int count = 0;
|
||||
for (int i = 0; i < modify->ncompute; i++)
|
||||
if (strcmp(modify->compute[i]->style,"cnp/atom") == 0) count++;
|
||||
if (count > 1 && comm->me == 0)
|
||||
error->warning(FLERR,"More than one compute cnp/atom defined");
|
||||
|
||||
// need an occasional full neighbor list
|
||||
int irequest = neighbor->request(this,instance_me);
|
||||
neighbor->requests[irequest]->pair = 0;
|
||||
neighbor->requests[irequest]->compute = 1;
|
||||
neighbor->requests[irequest]->half = 0;
|
||||
neighbor->requests[irequest]->full = 1;
|
||||
neighbor->requests[irequest]->occasional = 1;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void ComputeCNPAtom::init_list(int id, NeighList *ptr)
|
||||
{
|
||||
list = ptr;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void ComputeCNPAtom::compute_peratom()
|
||||
{
|
||||
int i,j,k,ii,jj,kk,m,n,inum,jnum,inear,jnear;
|
||||
int firstflag,ncommon;
|
||||
int *ilist,*jlist,*numneigh,**firstneigh;
|
||||
int onenearest[MAXNEAR];
|
||||
int common[MAXCOMMON];
|
||||
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
|
||||
double xjtmp,yjtmp,zjtmp,rjkx,rjky,rjkz;
|
||||
|
||||
invoked_peratom = update->ntimestep;
|
||||
|
||||
// grow arrays if necessary
|
||||
|
||||
if (atom->nmax > nmax) {
|
||||
memory->destroy(nearest);
|
||||
memory->destroy(nnearest);
|
||||
memory->destroy(cnpv);
|
||||
nmax = atom->nmax;
|
||||
memory->create(nearest,nmax,MAXNEAR,"cnp:nearest");
|
||||
memory->create(nnearest,nmax,"cnp:nnearest");
|
||||
memory->create(cnpv,nmax,"cnp:cnp_cnpv");
|
||||
vector_atom = cnpv;
|
||||
}
|
||||
|
||||
// invoke full neighbor list (will copy or build if necessary)
|
||||
|
||||
neighbor->build_one(list);
|
||||
|
||||
inum = list->inum;
|
||||
ilist = list->ilist;
|
||||
numneigh = list->numneigh;
|
||||
firstneigh = list->firstneigh;
|
||||
|
||||
// find the neigbors of each atom within cutoff using full neighbor list
|
||||
// nearest[] = atom indices of nearest neighbors, up to MAXNEAR
|
||||
// do this for all atoms, not just compute group
|
||||
// since CNP calculation requires neighbors of neighbors
|
||||
|
||||
double **x = atom->x;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
|
||||
int nerror = 0;
|
||||
for (ii = 0; ii < inum; ii++) {
|
||||
i = ilist[ii];
|
||||
xtmp = x[i][0];
|
||||
ytmp = x[i][1];
|
||||
ztmp = x[i][2];
|
||||
jlist = firstneigh[i];
|
||||
jnum = numneigh[i];
|
||||
|
||||
n = 0;
|
||||
for (jj = 0; jj < jnum; jj++) {
|
||||
j = jlist[jj];
|
||||
j &= NEIGHMASK;
|
||||
|
||||
delx = xtmp - x[j][0];
|
||||
dely = ytmp - x[j][1];
|
||||
delz = ztmp - x[j][2];
|
||||
rsq = delx*delx + dely*dely + delz*delz;
|
||||
if (rsq < cutsq) {
|
||||
if (n < MAXNEAR) nearest[i][n++] = j;
|
||||
else {
|
||||
nerror++;
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
nnearest[i] = n;
|
||||
}
|
||||
|
||||
// warning message
|
||||
|
||||
int nerrorall;
|
||||
MPI_Allreduce(&nerror,&nerrorall,1,MPI_INT,MPI_SUM,world);
|
||||
if (nerrorall && comm->me == 0) {
|
||||
char str[128];
|
||||
sprintf(str,"Too many neighbors in CNP for %d atoms",nerrorall);
|
||||
error->warning(FLERR,str,0);
|
||||
}
|
||||
|
||||
// compute CNP value for each atom in group
|
||||
// only performed if # of nearest neighbors = 12 or 14 (fcc,hcp)
|
||||
|
||||
nerror = 0;
|
||||
for (ii = 0; ii < inum; ii++) {
|
||||
i = ilist[ii];
|
||||
xtmp = x[i][0];
|
||||
ytmp = x[i][1];
|
||||
ztmp = x[i][2];
|
||||
// reset cnpv
|
||||
cnpv[i] = 0.0;
|
||||
|
||||
// skip computation of cnpv for atoms outside the compute group
|
||||
|
||||
if (!(mask[i] & groupbit)) continue;
|
||||
|
||||
// loop over nearest neighbors of I to build cnp data structure
|
||||
// cnp[k][NCOMMON] = # of common neighbors of I with each of its neighbors
|
||||
for (m = 0; m < nnearest[i]; m++) {
|
||||
j = nearest[i][m];
|
||||
xjtmp = x[j][0];
|
||||
yjtmp = x[j][1];
|
||||
zjtmp = x[j][2];
|
||||
|
||||
// common = list of neighbors common to atom I and atom J
|
||||
// if J is an owned atom, use its near neighbor list to find them
|
||||
// if J is a ghost atom, use full neighbor list of I to find them
|
||||
// in latter case, must exclude J from I's neighbor list
|
||||
|
||||
// find common neighbors of i and j using near neighbor list
|
||||
if (j < nlocal) {
|
||||
firstflag = 1;
|
||||
ncommon = 0;
|
||||
for (inear = 0; inear < nnearest[i]; inear++)
|
||||
for (jnear = 0; jnear < nnearest[j]; jnear++)
|
||||
if (nearest[i][inear] == nearest[j][jnear]) {
|
||||
if (ncommon < MAXCOMMON) common[ncommon++] = nearest[i][inear];
|
||||
else if (firstflag) {
|
||||
nerror++;
|
||||
firstflag = 0;
|
||||
}
|
||||
}
|
||||
|
||||
// find common neighbors of i and j using full neighbor list
|
||||
} else {
|
||||
jlist = firstneigh[i];
|
||||
jnum = numneigh[i];
|
||||
|
||||
n = 0;
|
||||
for (kk = 0; kk < jnum; kk++) {
|
||||
k = jlist[kk];
|
||||
k &= NEIGHMASK;
|
||||
if (k == j) continue;
|
||||
|
||||
delx = xjtmp - x[k][0];
|
||||
dely = yjtmp - x[k][1];
|
||||
delz = zjtmp - x[k][2];
|
||||
rsq = delx*delx + dely*dely + delz*delz;
|
||||
if (rsq < cutsq) {
|
||||
if (n < MAXNEAR) onenearest[n++] = k;
|
||||
else break;
|
||||
}
|
||||
}
|
||||
|
||||
firstflag = 1;
|
||||
ncommon = 0;
|
||||
for (inear = 0; inear < nnearest[i]; inear++)
|
||||
for (jnear = 0; (jnear < n) && (n < MAXNEAR); jnear++)
|
||||
if (nearest[i][inear] == onenearest[jnear]) {
|
||||
if (ncommon < MAXCOMMON) common[ncommon++] = nearest[i][inear];
|
||||
else if (firstflag) {
|
||||
nerror++;
|
||||
firstflag = 0;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
// Calculate and update sum |Rik+Rjk|ˆ2
|
||||
rjkx = 0.0;
|
||||
rjky = 0.0;
|
||||
rjkz = 0.0;
|
||||
for (kk = 0; kk < ncommon; kk++) {
|
||||
k = common[kk];
|
||||
rjkx += 2.0*x[k][0] - xjtmp - xtmp;
|
||||
rjky += 2.0*x[k][1] - yjtmp - ytmp;
|
||||
rjkz += 2.0*x[k][2] - zjtmp - ztmp;
|
||||
}
|
||||
// update cnpv with summed (valuejk)
|
||||
cnpv[i] += rjkx*rjkx + rjky*rjky + rjkz*rjkz;
|
||||
|
||||
// end of loop over j atoms
|
||||
}
|
||||
|
||||
// normalize cnp by the number of nearest neighbors
|
||||
cnpv[i] = cnpv[i] / nnearest[i];
|
||||
|
||||
// end of loop over i atoms
|
||||
}
|
||||
|
||||
// warning message
|
||||
MPI_Allreduce(&nerror,&nerrorall,1,MPI_INT,MPI_SUM,world);
|
||||
if (nerrorall && comm->me == 0) {
|
||||
char str[128];
|
||||
sprintf(str,"Too many common neighbors in CNP %d times",nerrorall);
|
||||
error->warning(FLERR,str);
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
memory usage of local atom-based array
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
double ComputeCNPAtom::memory_usage()
|
||||
{
|
||||
double bytes = nmax * sizeof(int);
|
||||
bytes += nmax * MAXNEAR * sizeof(int);
|
||||
bytes += nmax * sizeof(double);
|
||||
return bytes;
|
||||
}
|
||||
|
|
@ -0,0 +1,92 @@
|
|||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
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 COMPUTE_CLASS
|
||||
|
||||
ComputeStyle(cnp/atom,ComputeCNPAtom)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_COMPUTE_CNP_ATOM_H
|
||||
#define LMP_COMPUTE_CNP_ATOM_H
|
||||
|
||||
#include "compute.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class ComputeCNPAtom : public Compute {
|
||||
public:
|
||||
ComputeCNPAtom(class LAMMPS *, int, char **);
|
||||
~ComputeCNPAtom();
|
||||
void init();
|
||||
void init_list(int, class NeighList *);
|
||||
void compute_peratom();
|
||||
double memory_usage();
|
||||
|
||||
private:
|
||||
//revise
|
||||
int nmax;
|
||||
double cutsq;
|
||||
class NeighList *list;
|
||||
int **nearest;
|
||||
int *nnearest;
|
||||
double *cnpv;
|
||||
// int nmax;
|
||||
// double cutsq;
|
||||
// class NeighList *list;
|
||||
// int **nearest;
|
||||
// int *nnearest;
|
||||
// double *pattern;
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Illegal ... command
|
||||
|
||||
Self-explanatory. Check the input script syntax and compare to the
|
||||
documentation for the command. You can use -echo screen as a
|
||||
command-line option when running LAMMPS to see the offending line.
|
||||
|
||||
E: Compute cnp/atom requires a pair style be defined
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Compute cnp/atom cutoff is longer than pairwise cutoff
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
W: Compute cnp/atom cutoff may be too large to find ghost atom neighbors
|
||||
|
||||
The neighbor cutoff used may not encompass enough ghost atoms
|
||||
to perform this operation correctly.
|
||||
|
||||
W: More than one compute cnp/atom defined
|
||||
|
||||
It is not efficient to use compute cnp/atom more than once.
|
||||
|
||||
W: Too many neighbors in CNP for %d atoms
|
||||
|
||||
More than the maximum # of neighbors was found multiple times. This
|
||||
was unexpected.
|
||||
|
||||
W: Too many common neighbors in CNP %d times
|
||||
|
||||
More than the maximum # of neighbors was found multiple times. This
|
||||
was unexpected.
|
||||
|
||||
*/
|
||||
Loading…
Reference in New Issue