Merge branch 'remove-reax-meam' into collected-post-stable-patches

This commit is contained in:
Axel Kohlmeyer 2018-12-12 00:14:58 -05:00
commit 02bdda0d05
161 changed files with 87 additions and 21143 deletions

View File

@ -171,7 +171,7 @@ set(LAMMPS_DEPS)
set(LAMMPS_API_DEFINES)
set(DEFAULT_PACKAGES ASPHERE BODY CLASS2 COLLOID COMPRESS DIPOLE GRANULAR
KSPACE MANYBODY MC MEAM MESSAGE MISC MOLECULE PERI REAX REPLICA RIGID SHOCK
KSPACE MANYBODY MC MESSAGE MISC MOLECULE PERI REPLICA RIGID SHOCK
SPIN SNAP SRD KIM PYTHON MSCG MPIIO VORONOI POEMS LATTE USER-ATC USER-AWPMD
USER-BOCS USER-CGDNA USER-MESO USER-CGSDK USER-COLVARS USER-DIFFRACTION
USER-DPD USER-DRUDE USER-EFF USER-FEP USER-H5MD USER-LB USER-MANIFOLD
@ -191,11 +191,11 @@ endforeach()
######################################################
# packages with special compiler needs or external libs
######################################################
if(PKG_REAX OR PKG_MEAM OR PKG_USER-QUIP OR PKG_USER-QMMM OR PKG_LATTE OR PKG_USER-SCAFACOS)
if(PKG_USER-QUIP OR PKG_USER-QMMM OR PKG_LATTE OR PKG_USER-SCAFACOS)
enable_language(Fortran)
endif()
if(PKG_MEAM OR PKG_USER-H5MD OR PKG_USER-QMMM OR PKG_USER-SCAFACOS)
if(PKG_USER-H5MD OR PKG_USER-QMMM OR PKG_USER-SCAFACOS)
enable_language(C)
endif()
@ -826,7 +826,7 @@ endforeach()
##############################################
# add lib sources of (simple) enabled packages
############################################
foreach(SIMPLE_LIB REAX MEAM POEMS USER-ATC USER-AWPMD USER-COLVARS USER-H5MD
foreach(SIMPLE_LIB POEMS USER-ATC USER-AWPMD USER-COLVARS USER-H5MD
USER-QMMM)
if(PKG_${SIMPLE_LIB})
string(REGEX REPLACE "^USER-" "" PKG_LIB "${SIMPLE_LIB}")

View File

@ -36,7 +36,6 @@ This is the list of packages that may require additional steps.
"OPT"_#opt,
"POEMS"_#poems,
"PYTHON"_#python,
"REAX"_#reax,
"VORONOI"_#voronoi,
"USER-ATC"_#user-atc,
"USER-AWPMD"_#user-awpmd,
@ -536,45 +535,6 @@ build fails.
:line
REAX package :h4,link(reax)
NOTE: the use of the REAX package and its "pair_style
reax"_pair_reax.html command is discouraged, as it is no longer
maintained. Please use the USER-REAXC package and its "pair_style
reax/c"_pair_reaxc.html command instead, and possibly its KOKKOS
enabled variant (pair_style reax/c/kk), which has a more robust memory
management. See the "pair_style reax/c"_pair_reaxc.html doc page for
details.
[CMake build]:
No additional settings are needed besides "-D PKG_REAX=yes".
[Traditional make]:
Before building LAMMPS, you must build the REAX library in lib/reax.
You can do this manually if you prefer; follow the instructions in
lib/reax/README. You can also do it in one step from the lammps/src
dir, using a command like these, which simply invoke the
lib/reax/Install.py script with the specified args:
make lib-reax # print help message
make lib-reax args="-m serial" # build with GNU Fortran compiler (settings as with "make serial")
make lib-reax args="-m mpi" # build with default MPI Fortran compiler (settings as with "make mpi")
make lib-reax args="-m ifort" # build with Intel ifort compiler :pre
The build should produce two files: lib/reax/libreax.a and
lib/reax/Makefile.lammps. The latter is copied from an existing
Makefile.lammps.* and has settings needed to link C++ (LAMMPS) with
Fortran (REAX library). Typically the two compilers used for LAMMPS
and the REAX library need to be consistent (e.g. both Intel or both
GNU compilers). If necessary, you can edit/create a new
lib/reax/Makefile.machine file for your system, which should define an
EXTRAMAKE variable to specify a corresponding Makefile.lammps.machine
file.
:line
VORONOI package :h4,link(voronoi)
To build with this package, you must download and build the "Voro++

View File

@ -47,7 +47,6 @@ packages:
"OPT"_Build_extras.html#opt,
"POEMS"_Build_extras.html#poems,
"PYTHON"_Build_extras.html#python,
"REAX"_Build_extras.html#reax,
"VORONOI"_Build_extras.html#voronoi,
"USER-ATC"_Build_extras.html#user-atc,
"USER-AWPMD"_Build_extras.html#user-awpmd,

View File

@ -169,8 +169,7 @@ OPT.
"qmmm"_fix_qmmm.html,
"qtb"_fix_qtb.html,
"rattle"_fix_shake.html,
"reax/bonds"_fix_reax_bonds.html,
"reax/c/bonds (k)"_fix_reax_bonds.html,
"reax/c/bonds (k)"_fix_reaxc_bonds.html,
"reax/c/species (k)"_fix_reaxc_species.html,
"recenter"_fix_recenter.html,
"restrain"_fix_restrain.html,

View File

@ -160,8 +160,7 @@ OPT.
"lubricateU/poly"_pair_lubricateU.html,
"mdpd"_pair_meso.html,
"mdpd/rhosum"_pair_meso.html,
"meam"_pair_meam.html,
"meam/c"_pair_meam.html,
"meam/c"_pair_meamc.html,
"meam/spline (o)"_pair_meam_spline.html,
"meam/sw/spline"_pair_meam_sw_spline.html,
"mgpt"_pair_mgpt.html,
@ -194,7 +193,6 @@ OPT.
"polymorphic"_pair_polymorphic.html,
"python"_pair_python.html,
"quip"_pair_quip.html,
"reax"_pair_reax.html,
"reax/c (ko)"_pair_reaxc.html,
"rebo (io)"_pair_airebo.html,
"resquared (go)"_pair_resquared.html,

View File

@ -23,8 +23,8 @@ install the Windows MPI package (MPICH2 from Argonne National Labs),
needed to run in parallel.
The LAMMPS binaries contain all optional packages included in the
source distribution except: KIM, REAX, KOKKOS, USER-INTEL,
and USER-QMMM. The serial version also does not include the MPIIO and
source distribution except: KIM, KOKKOS, USER-INTEL, and USER-QMMM.
The serial version also does not include the MPIIO and
USER-LB packages. GPU support is provided for OpenCL.
The installer site also has instructions on how to run LAMMPS under

View File

@ -45,7 +45,6 @@ as contained in the file name.
"LATTE"_#PKG-LATTE,
"MANYBODY"_#PKG-MANYBODY,
"MC"_#PKG-MC,
"MEAM"_#PKG-MEAM,
"MESSAGE"_#PKG-MESSAGE,
"MISC"_#PKG-MISC,
"MOLECULE"_#PKG-MOLECULE,
@ -56,7 +55,6 @@ as contained in the file name.
"POEMS"_#PKG-POEMS,
"PYTHON"_#PKG-PYTHON,
"QEQ"_#PKG-QEQ,
"REAX"_#PKG-REAX,
"REPLICA"_#PKG-REPLICA2,
"RIGID"_#PKG-RIGID,
"SHOCK"_#PKG-SHOCK,
@ -533,37 +531,6 @@ http://lammps.sandia.gov/movies.html#gcmc :ul
:line
MEAM package :link(PKG-MEAM),h4
[Contents:]
A pair style for the modified embedded atom (MEAM) potential.
Please note that the use of the MEAM package is discouraged as
it has been superseded by the "USER-MEAMC"_#PKG-USER-MEAMC package,
which is a direct translation of the MEAM package to C++.
USER-MEAMC contains additional optimizations making it run faster
than MEAM on most machines, while providing the identical features
and user interface.
[Author:] Greg Wagner (Northwestern U) while at Sandia.
[Install:]
This package has "specific installation
instructions"_Build_extras.html#gpu on the "Build
extras"_Build_extras.html doc page.
[Supporting info:]
src/MEAM: filenames -> commands
src/meam/README
lib/meam/README
"pair_style meam"_pair_meam.html
examples/meam :ul
:line
MESSAGE package :link(PKG-MESSAGE),h4
[Contents:]
@ -834,38 +801,6 @@ examples/streitz :ul
:line
REAX package :link(PKG-REAX),h4
[Contents:]
NOTE: the use of the REAX package is discouraged, as it is no longer
maintained. Please use the "USER-REAXC"_#PKG-USER-REAXC package instead,
and possibly the KOKKOS enabled variant of that, which has a more robust
memory management.
A pair style which wraps a Fortran library which implements the ReaxFF
potential, which is a universal reactive force field. Also included is
a "fix reax/bonds"_fix_reax_bonds.html command for monitoring molecules
as bonds are created and destroyed.
[Author:] Aidan Thompson (Sandia).
[Install:]
This package has "specific installation
instructions"_Build_extras.html#gpu on the "Build
extras"_Build_extras.html doc page.
[Supporting info:]
src/REAX: filenames -> commands
lib/reax/README
"pair_style reax"_pair_reax.html
"fix reax/bonds"_fix_reax_bonds.html
examples/reax :ul
:line
REPLICA package :link(PKG-REPLICA2),h4
[Contents:]
@ -1552,10 +1487,9 @@ USER-MEAMC package :link(PKG-USER-MEAMC),h4
[Contents:]
A pair style for the modified embedded atom (MEAM) potential
translated from the Fortran version in the "MEAM"_#PKG-MEAM package
to plain C++. In contrast to the MEAM package, no library
needs to be compiled and the pair style can be instantiated
multiple times.
translated from the Fortran version in the (obsolete) "MEAM" package
to plain C++. The USER-MEAMC fully replaces the MEAM package, which
has been removed from LAMMPS after the 12 December 2018 version.
[Author:] Sebastian Huetter, (Otto-von-Guericke University Magdeburg)
based on the Fortran version of Greg Wagner (Northwestern U) while at
@ -1565,8 +1499,8 @@ Sandia.
src/USER-MEAMC: filenames -> commands
src/USER-MEAMC/README
"pair_style meam/c"_pair_meam.html
examples/meam :ul
"pair_style meam/c"_pair_meamc.html
examples/meamc :ul
:line
@ -1894,9 +1828,8 @@ USER-REAXC package :link(PKG-USER-REAXC),h4
[Contents:]
A pair style which implements the ReaxFF potential in C/C++ (in
contrast to the "REAX package"_#PKG-REAX and its Fortran library). ReaxFF
is universal reactive force field. See the src/USER-REAXC/README file
A pair style which implements the ReaxFF potential in C/C++. ReaxFF
is a universal reactive force field. See the src/USER-REAXC/README file
for more info on differences between the two packages. Also two fixes
for monitoring molecules as bonds are created and destroyed.
@ -1907,7 +1840,7 @@ for monitoring molecules as bonds are created and destroyed.
src/USER-REAXC: filenames -> commands
src/USER-REAXC/README
"pair_style reax/c"_pair_reaxc.html
"fix reax/c/bonds"_fix_reax_bonds.html
"fix reax/c/bonds"_fix_reaxc_bonds.html
"fix reax/c/species"_fix_reaxc_species.html
examples/reax :ul

View File

@ -46,7 +46,6 @@ Package, Description, Doc page, Example, Library
"LATTE"_Packages_details.html#PKG-LATTE, quantum DFTB forces via LATTE, "fix latte"_fix_latte.html, latte, ext
"MANYBODY"_Packages_details.html#PKG-MANYBODY, many-body potentials, "pair_style tersoff"_pair_tersoff.html, shear, no
"MC"_Packages_details.html#PKG-MC, Monte Carlo options, "fix gcmc"_fix_gcmc.html, n/a, no
"MEAM"_Packages_details.html#PKG-MEAM, modified EAM potential, "pair_style meam"_pair_meam.html, meam, int
"MESSAGE"_Packages_details.html#PKG-MESSAGE, client/server messaging, "message"_message.html, message, int
"MISC"_Packages_details.html#PKG-MISC, miscellaneous single-file commands, n/a, no, no
"MOLECULE"_Packages_details.html#PKG-MOLECULE, molecular system force fields, "Howto bioFF"_Howto_bioFF.html, peptide, no
@ -57,7 +56,6 @@ Package, Description, Doc page, Example, Library
"POEMS"_Packages_details.html#PKG-POEMS, coupled rigid body motion, "fix poems"_fix_poems.html, rigid, int
"PYTHON"_Packages_details.html#PKG-PYTHON, embed Python code in an input script, "python"_python.html, python, sys
"QEQ"_Packages_details.html#PKG-QEQ, QEq charge equilibration, "fix qeq"_fix_qeq.html, qeq, no
"REAX"_Packages_details.html#PKG-REAX, ReaxFF potential (Fortran), "pair_style reax"_pair_reax.html, reax, int
"REPLICA"_Packages_details.html#PKG-REPLICA2, multi-replica methods, "Howto replica"_Howto_replica.html, tad, no
"RIGID"_Packages_details.html#PKG-RIGID, rigid bodies and constraints, "fix rigid"_fix_rigid.html, rigid, no
"SHOCK"_Packages_details.html#PKG-SHOCK, shock loading methods, "fix msst"_fix_msst.html, n/a, no

View File

@ -53,7 +53,7 @@ Package, Description, Doc page, Example, Library
"USER-INTEL"_Packages_details.html#PKG-USER-INTEL, optimized Intel CPU and KNL styles,"Speed intel"_Speed_intel.html, "Benchmarks"_http://lammps.sandia.gov/bench.html, no
"USER-LB"_Packages_details.html#PKG-USER-LB, Lattice Boltzmann fluid,"fix lb/fluid"_fix_lb_fluid.html, USER/lb, no
"USER-MANIFOLD"_Packages_details.html#PKG-USER-MANIFOLD, motion on 2d surfaces,"fix manifoldforce"_fix_manifoldforce.html, USER/manifold, no
"USER-MEAMC"_Packages_details.html#PKG-USER-MEAMC, modified EAM potential (C++), "pair_style meam/c"_pair_meam.html, meam, no
"USER-MEAMC"_Packages_details.html#PKG-USER-MEAMC, modified EAM potential (C++), "pair_style meam/c"_pair_meamc.html, meamc, no
"USER-MESO"_Packages_details.html#PKG-USER-MESO, mesoscale DPD models, "pair_style edpd"_pair_meso.html, USER/meso, no
"USER-MGPT"_Packages_details.html#PKG-USER-MGPT, fast MGPT multi-ion potentials, "pair_style mgpt"_pair_mgpt.html, USER/mgpt, no
"USER-MISC"_Packages_details.html#PKG-USER-MISC, single-file contributions, USER-MISC/README, USER/misc, no

View File

@ -486,8 +486,8 @@ README for more info on Pizza.py and how to use these scripts.
reax tool :h4,link(reax_tool)
The reax sub-directory contains stand-alone codes that can
post-process the output of the "fix reax/bonds"_fix_reax_bonds.html
command from a LAMMPS simulation using "ReaxFF"_pair_reax.html. See
post-process the output of the "fix reax/c/bonds"_fix_reaxc_bonds.html
command from a LAMMPS simulation using "ReaxFF"_pair_reaxc.html. See
the README.txt file for more info.
These tools were written by Aidan Thompson at Sandia.

View File

@ -24,7 +24,7 @@ nsub = {n}-instance of a sub-style, if a pair style is used multiple times in a
compute 1 all pair gauss
compute 1 all pair lj/cut/coul/cut ecoul
compute 1 all pair tersoff 2 epair
compute 1 all pair reax :pre
compute 1 all pair reax/c :pre
[Description:]
@ -60,8 +60,8 @@ corrections, even if they are enabled via the
"pair_modify"_pair_modify.html command.
Some pair styles tally additional quantities, e.g. a breakdown of
potential energy into a dozen or so components is tallied by the
"pair_style reax"_pair_reax.html command. These values (1 or more)
potential energy into 14 components is tallied by the "pair_style
reax/c"_pair_reaxc.html command. These values (1 or more)
are stored as a global vector by this compute. See the doc page for
"individual pair styles"_pair_style.html for info on these values.

View File

@ -312,9 +312,8 @@ accelerated styles exist.
"qmmm"_fix_qmmm.html -
"qtb"_fix_qtb.html -
"rattle"_fix_shake.html - RATTLE constraints on bonds and/or angles
"reax/bonds"_fix_reax_bonds.html - write out ReaxFF bond information
"reax/c/bonds"_fix_reax_bonds.html -
"reax/c/species"_fix_reaxc_species.html -
"reax/c/bonds"_fix_reaxc_bonds.html - write out ReaxFF bond information
"reax/c/species"_fix_reaxc_species.html - write out ReaxFF molecule information
"recenter"_fix_recenter.html - constrain the center-of-mass position of a group of atoms
"restrain"_fix_restrain.html - constrain a bond, angle, dihedral
"rhok"_fix_rhok.html -

View File

@ -6,13 +6,12 @@
:line
fix reax/bonds command :h3
fix reax/c/bonds command :h3
fix reax/c/bonds/kk command :h3
[Syntax:]
fix ID group-ID reax/bonds Nevery filename :pre
fix ID group-ID reaxc/bonds Nevery filename :pre
ID, group-ID are documented in "fix"_fix.html command
reax/bonds = style name of this fix command
@ -21,16 +20,14 @@ filename = name of output file :ul
[Examples:]
fix 1 all reax/bonds 100 bonds.tatb
fix 1 all reax/c/bonds 100 bonds.reaxc :pre
[Description:]
Write out the bond information computed by the ReaxFF potential
specified by "pair_style reax"_pair_reax.html or "pair_style
reax/c"_pair_reaxc.html in the exact same format as the original
stand-alone ReaxFF code of Adri van Duin. The bond information is
written to {filename} on timesteps that are multiples of {Nevery},
Write out the bond information computed by the ReaxFF potential specified
by "pair_style reax/c"_pair_reaxc.html in the exact same format as the
original stand-alone ReaxFF code of Adri van Duin. The bond information
is written to {filename} on timesteps that are multiples of {Nevery},
including timestep 0. For time-averaged chemical species analysis,
please see the "fix reaxc/c/species"_fix_reaxc_species.html command.
@ -94,12 +91,8 @@ more instructions on how to use the accelerated styles effectively.
[Restrictions:]
The fix reax/bonds command requires that the "pair_style
reax"_pair_reax.html be invoked. This fix is part of the REAX
package. It is only enabled if LAMMPS was built with that package,
which also requires the REAX library be built and linked with LAMMPS.
The fix reax/c/bonds command requires that the "pair_style
reax/c"_pair_reaxc.html be invoked. This fix is part of the
reax/c"_pair_reaxc.html is invoked. This fix is part of the
USER-REAXC package. It is only enabled if LAMMPS was built with that
package. See the "Build package"_Build_package.html doc page for more
info.
@ -109,7 +102,6 @@ To write gzipped bond files, you must compile LAMMPS with the
[Related commands:]
"pair_style reax"_pair_reax.html, "pair_style
reax/c"_pair_reaxc.html, "fix reax/c/species"_fix_reaxc_species.html
"pair_style reax/c"_pair_reaxc.html, "fix reax/c/species"_fix_reaxc_species.html
[Default:] none

View File

@ -161,7 +161,7 @@ more instructions on how to use the accelerated styles effectively.
[Restrictions:]
The fix species currently only works with "pair_style
The "fix reax/c/species" currently only works with "pair_style
reax/c"_pair_reaxc.html and it requires that the "pair_style
reax/c"_pair_reaxc.html be invoked. This fix is part of the
USER-REAXC package. It is only enabled if LAMMPS was built with that
@ -177,8 +177,7 @@ It should be possible to extend it to other reactive pair_styles (such as
[Related commands:]
"pair_style reax/c"_pair_reaxc.html, "fix
reax/bonds"_fix_reax_bonds.html
"pair_style reax/c"_pair_reaxc.html, "fix reax/c/bonds"_fix_reaxc_bonds.html
[Default:]

View File

@ -135,7 +135,7 @@ Fixes :h1
fix_qeq_reax
fix_qmmm
fix_qtb
fix_reax_bonds
fix_reaxc_bonds
fix_reaxc_species
fix_recenter
fix_restrain

View File

@ -356,7 +356,7 @@ fix_qeq_comb.html
fix_qeq_reax.html
fix_qmmm.html
fix_qtb.html
fix_reax_bonds.html
fix_reaxc_bonds.html
fix_reaxc_species.html
fix_recenter.html
fix_restrain.html
@ -598,7 +598,7 @@ pair_lj_soft.html
pair_lubricate.html
pair_lubricateU.html
pair_mdf.html
pair_meam.html
pair_meamc.html
pair_meam_spline.html
pair_meam_sw_spline.html
pair_meso.html
@ -617,7 +617,6 @@ pair_peri.html
pair_polymorphic.html
pair_python.html
pair_quip.html
pair_reax.html
pair_reaxc.html
pair_resquared.html
pair_sdk.html

View File

@ -70,15 +70,10 @@ other pairwise potential for several different atom type pairs in your
model, then you should just list the sub-style once and use the
pair_coeff command to assign parameters for the different type pairs.
NOTE: There are two exceptions to this option to list an individual
pair style multiple times. The first is for pair styles implemented
as Fortran libraries: "pair_style meam"_pair_meam.html and "pair_style
reax"_pair_reax.html ("pair_style reax/c"_pair_reaxc.html is OK).
This is because unlike a C++ class, they can not be instantiated
multiple times, due to the manner in which they were coded in Fortran.
The second is for GPU-enabled pair styles in the GPU package. This is
b/c the GPU package also currently assumes that only one instance of a
pair style is being used.
NOTE: There is one exception to this option to list an individual
pair style multiple times: GPU-enabled pair styles in the GPU package.
This is because the GPU package currently assumes that only one
instance of a pair style is being used.
In the pair_coeff commands, the name of a pair style must be added
after the I,J type specification, with the remaining coefficients

View File

@ -152,7 +152,7 @@ info.
[Related commands:]
"pair_coeff"_pair_coeff.html, "pair_style meam"_pair_meam.html
"pair_coeff"_pair_coeff.html, "pair_style meam/c"_pair_meamc.html
[Default:] none

View File

@ -116,7 +116,7 @@ info.
[Related commands:]
"pair_coeff"_pair_coeff.html, "pair_style meam"_pair_meam.html,
"pair_coeff"_pair_coeff.html, "pair_style meam/c"_pair_meamc.html,
"pair_style meam/spline"_pair_meam_spline.html
[Default:] none

View File

@ -6,18 +6,17 @@
:line
pair_style meam command :h3
pair_style meam/c command :h3
[Syntax:]
pair_style style :pre
style = {meam} or {meam/c}
style = {meam/c}
[Examples:]
pair_style meam
pair_style meam/c
pair_coeff * * ../potentials/library.meam Si ../potentials/si.meam Si
pair_coeff * * ../potentials/library.meam Ni Al NULL Ni Al Ni Ni :pre
@ -27,14 +26,16 @@ NOTE: The behavior of the MEAM potential for alloy systems has changed
as of November 2010; see description below of the mixture_ref_t
parameter
Style {meam} computes pairwise interactions for a variety of materials
Style {meam/c} computes pairwise interactions for a variety of materials
using modified embedded-atom method (MEAM) potentials
"(Baskes)"_#Baskes. Conceptually, it is an extension to the original
"EAM potentials"_pair_eam.html which adds angular forces. It is
thus suitable for modeling metals and alloys with fcc, bcc, hcp and
diamond cubic structures, as well as covalently bonded materials like
silicon and carbon. Style {meam/c} is a translation of the {meam} code
from (mostly) Fortran to C++. It is functionally equivalent to {meam}.
silicon and carbon. Style {meam/c} is a translation of the (now obsolete)
{meam} code from Fortran to C++. It is functionally equivalent to {meam}
but more efficient, and thus {meam} has been removed from LAMMPS after
the 12 December 2018 release.
In the MEAM formulation, the total energy E of a system of atoms is
given by:
@ -352,13 +353,8 @@ This pair style can only be used via the {pair} keyword of the
[Restrictions:]
The {meam} style is part of the MEAM package. It is only enabled if
LAMMPS was built with that package, which also requires the MEAM
library be built and linked with LAMMPS. The {meam/c} style is
provided in the USER-MEAMC package. It is only enabled if LAMMPS was
built with that package. In contrast to the {meam} style, {meam/c}
does not require a separate library to be compiled and it can be
instantiated multiple times in a "hybrid"_pair_hybrid.html pair style.
The {meam/c} style is provided in the USER-MEAMC package. It is
only enabled if LAMMPS was built with that package.
See the "Build package"_Build_package.html doc page for more info.
[Related commands:]

View File

@ -1,216 +0,0 @@
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
:link(lws,http://lammps.sandia.gov)
:link(ld,Manual.html)
:link(lc,Commands_all.html)
:line
pair_style reax command :h3
[Syntax:]
pair_style reax hbcut hbnewflag tripflag precision :pre
hbcut = hydrogen-bond cutoff (optional) (distance units)
hbnewflag = use old or new hbond function style (0 or 1) (optional)
tripflag = apply stabilization to all triple bonds (0 or 1) (optional)
precision = precision for charge equilibration (optional) :ul
[Examples:]
pair_style reax
pair_style reax 10.0 0 1 1.0e-5
pair_coeff * * ffield.reax 3 1 2 2
pair_coeff * * ffield.reax 3 NULL NULL 3 :pre
[Description:]
Style {reax} computes the ReaxFF potential of van Duin, Goddard and
co-workers. ReaxFF uses distance-dependent bond-order functions to
represent the contributions of chemical bonding to the potential
energy. There is more than one version of ReaxFF. The version
implemented in LAMMPS uses the functional forms documented in the
supplemental information of the following paper:
"(Chenoweth)"_#Chenoweth_20081. The version integrated into LAMMPS matches
the most up-to-date version of ReaxFF as of summer 2010.
WARNING: pair style reax is now deprecated and will soon be retired. Users
should switch to "pair_style reax/c"_pair_reaxc.html. The {reax} style
differs from the {reax/c} style in the lo-level implementation details.
The {reax} style is a
Fortran library, linked to LAMMPS. The {reax/c} style was initially
implemented as stand-alone C code and is now integrated into LAMMPS as
a package.
LAMMPS requires that a file called ffield.reax be provided, containing
the ReaxFF parameters for each atom type, bond type, etc. The format
is identical to the ffield file used by van Duin and co-workers. The
filename is required as an argument in the pair_coeff command. Any
value other than "ffield.reax" will be rejected (see below).
LAMMPS provides several different versions of ffield.reax in its
potentials dir, each called potentials/ffield.reax.label. These are
documented in potentials/README.reax. The default ffield.reax
contains parameterizations for the following elements: C, H, O, N.
NOTE: We do not distribute a wide variety of ReaxFF force field files
with LAMMPS. Adri van Duin's group at PSU is the central repository
for this kind of data as they are continuously deriving and updating
parameterizations for different classes of materials. You can submit
a contact request at the Materials Computation Center (MCC) website
"https://www.mri.psu.edu/materials-computation-center/connect-mcc"_https://www.mri.psu.edu/materials-computation-center/connect-mcc,
describing the material(s) you are interested in modeling with ReaxFF.
They can tell
you what is currently available or what it would take to create a
suitable ReaxFF parameterization.
The format of these files is identical to that used originally by van
Duin. We have tested the accuracy of {pair_style reax} potential
against the original ReaxFF code for the systems mentioned above. You
can use other ffield files for specific chemical systems that may be
available elsewhere (but note that their accuracy may not have been
tested).
The {hbcut}, {hbnewflag}, {tripflag}, and {precision} settings are
optional arguments. If none are provided, default settings are used:
{hbcut} = 6 (which is Angstroms in real units), {hbnewflag} = 1 (use
new hbond function style), {tripflag} = 1 (apply stabilization to all
triple bonds), and {precision} = 1.0e-6 (one part in 10^6). If you
wish to override any of these defaults, then all of the settings must
be specified.
Two examples using {pair_style reax} are provided in the examples/reax
sub-directory, along with corresponding examples for
"pair_style reax/c"_pair_reaxc.html. Note that while the energy and force
calculated by both of these pair styles match very closely, the
contributions due to the valence angles differ slightly due to
the fact that with {pair_style reax/c} the default value of {thb_cutoff_sq}
is 0.00001, while for {pair_style reax} it is hard-coded to be 0.001.
Use of this pair style requires that a charge be defined for every
atom since the {reax} pair style performs a charge equilibration (QEq)
calculation. See the "atom_style"_atom_style.html and
"read_data"_read_data.html commands for details on how to specify
charges.
The thermo variable {evdwl} stores the sum of all the ReaxFF potential
energy contributions, with the exception of the Coulombic and charge
equilibration contributions which are stored in the thermo variable
{ecoul}. The output of these quantities is controlled by the
"thermo"_thermo.html command.
This pair style tallies a breakdown of the total ReaxFF potential
energy into sub-categories, which can be accessed via the "compute
pair"_compute_pair.html command as a vector of values of length 14.
The 14 values correspond to the following sub-categories (the variable
names in italics match those used in the ReaxFF FORTRAN library):
{eb} = bond energy
{ea} = atom energy
{elp} = lone-pair energy
{emol} = molecule energy (always 0.0)
{ev} = valence angle energy
{epen} = double-bond valence angle penalty
{ecoa} = valence angle conjugation energy
{ehb} = hydrogen bond energy
{et} = torsion energy
{eco} = conjugation energy
{ew} = van der Waals energy
{ep} = Coulomb energy
{efi} = electric field energy (always 0.0)
{eqeq} = charge equilibration energy :ol
To print these quantities to the log file (with descriptive column
headings) the following commands could be included in an input script:
compute reax all pair reax
variable eb equal c_reax\[1\]
variable ea equal c_reax\[2\]
...
variable eqeq equal c_reax\[14\]
thermo_style custom step temp epair v_eb v_ea ... v_eqeq :pre
Only a single pair_coeff command is used with the {reax} style which
specifies a ReaxFF potential file with parameters for all needed
elements. These are mapped to LAMMPS atom types by specifying N
additional arguments after the filename in the pair_coeff command,
where N is the number of LAMMPS atom types:
filename
N indices = mapping of ReaxFF elements to atom types :ul
The specification of the filename and the mapping of LAMMPS atom types
recognized by the ReaxFF is done differently than for other LAMMPS
potentials, due to the non-portable difficulty of passing character
strings (e.g. filename, element names) between C++ and Fortran.
The filename has to be "ffield.reax" and it has to exist in the
directory you are running LAMMPS in. This means you cannot prepend a
path to the file in the potentials dir. Rather, you should copy that
file into the directory you are running from. If you wish to use
another ReaxFF potential file, then name it "ffield.reax" and put it
in the directory you run from.
In the ReaxFF potential file, near the top, after the general
parameters, is the atomic parameters section that contains element
names, each with a couple dozen numeric parameters. If there are M
elements specified in the {ffield} file, think of these as numbered 1
to M. Each of the N indices you specify for the N atom types of LAMMPS
atoms must be an integer from 1 to M. Atoms with LAMMPS type 1 will
be mapped to whatever element you specify as the first index value,
etc. If a mapping value is specified as NULL, the mapping is not
performed. This can be used when a ReaxFF potential is used as part
of the {hybrid} pair style. The NULL values are placeholders for atom
types that will be used with other potentials.
NOTE: Currently the reax pair style cannot be used as part of the
{hybrid} pair style. Some additional changes still need to be made to
enable this.
As an example, say your LAMMPS simulation has 4 atom types and the
elements are ordered as C, H, O, N in the {ffield} file. If you want
the LAMMPS atom type 1 and 2 to be C, type 3 to be N, and type 4 to be
H, you would use the following pair_coeff command:
pair_coeff * * ffield.reax 1 1 4 2 :pre
:line
[Mixing, shift, table, tail correction, restart, rRESPA info]:
This pair style does not support the "pair_modify"_pair_modify.html
mix, shift, table, and tail options.
This pair style does not write its information to "binary restart
files"_restart.html, 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.
This pair style can only be used via the {pair} keyword of the
"run_style respa"_run_style.html command. It does not support the
{inner}, {middle}, {outer} keywords.
[Restrictions:]
The ReaxFF potential files provided with LAMMPS in the potentials
directory are parameterized for real "units"_units.html. You can use
the ReaxFF potential with any LAMMPS units, but you would need to
create your own potential file with coefficients listed in the
appropriate units if your simulation doesn't use "real" units.
[Related commands:]
"pair_coeff"_pair_coeff.html, "pair_style reax/c"_pair_reaxc.html,
"fix_reax_bonds"_fix_reax_bonds.html
[Default:]
The keyword defaults are {hbcut} = 6, {hbnewflag} = 1, {tripflag} = 1,
{precision} = 1.0e-6.
:line
:link(Chenoweth_20081)
[(Chenoweth_2008)] Chenoweth, van Duin and Goddard,
Journal of Physical Chemistry A, 112, 1040-1053 (2008).

View File

@ -37,7 +37,7 @@ pair_coeff * * ffield.reax C H O N :pre
Style {reax/c} computes the ReaxFF potential of van Duin, Goddard and
co-workers. ReaxFF uses distance-dependent bond-order functions to
represent the contributions of chemical bonding to the potential
energy. There is more than one version of ReaxFF. The version
energy. There is more than one version of ReaxFF. The version
implemented in LAMMPS uses the functional forms documented in the
supplemental information of the following paper: "(Chenoweth et al.,
2008)"_#Chenoweth_20082. The version integrated into LAMMPS matches
@ -56,11 +56,10 @@ consideration when using the {reax/c/kk} style is the choice of either
half or full neighbor lists. This setting can be changed using the
Kokkos "package"_package.html command.
The {reax/c} style differs from the "pair_style reax"_pair_reax.html
command in the lo-level implementation details. The {reax} style is a
Fortran library, linked to LAMMPS. The {reax/c} style was initially
implemented as stand-alone C code and is now integrated into LAMMPS as
a package.
The {reax/c} style differs from the (obsolete) "pair_style reax"
command in the implementation details. The {reax} style was a
Fortran library, linked to LAMMPS. The {reax} style has been removed
from LAMMPS after the 12 December 2018 version.
LAMMPS provides several different versions of ffield.reax in its
potentials dir, each called potentials/ffield.reax.label. These are
@ -98,9 +97,8 @@ correspond to those used by Adri van Duin's stand-alone serial
code. If these are changed by setting control variables in the control
file, the results from LAMMPS and the serial code will not agree.
Two examples using {pair_style reax/c} are provided in the examples/reax
sub-directory, along with corresponding examples for
"pair_style reax"_pair_reax.html.
Examples using {pair_style reax/c} are provided in the examples/reax
sub-directory.
Use of this pair style requires that a charge be defined for every
atom. See the "atom_style"_atom_style.html and
@ -193,8 +191,7 @@ where N is the number of LAMMPS atom types:
filename
N indices = ReaxFF elements :ul
The filename is the ReaxFF potential file. Unlike for the {reax}
pair style, any filename can be used.
The filename is the ReaxFF potential file.
In the ReaxFF potential file, near the top, after the general
parameters, is the atomic parameters section that contains element
@ -337,9 +334,8 @@ appropriate units if your simulation doesn't use "real" units.
[Related commands:]
"pair_coeff"_pair_coeff.html, "fix qeq/reax"_fix_qeq_reax.html, "fix
reax/c/bonds"_fix_reax_bonds.html, "fix
reax/c/species"_fix_reaxc_species.html, "pair_style
reax"_pair_reax.html
reax/c/bonds"_fix_reaxc_bonds.html, "fix
reax/c/species"_fix_reaxc_species.html
[Default:]

View File

@ -226,8 +226,7 @@ accelerated styles exist.
"lubricateU/poly"_pair_lubricateU.html - hydrodynamic lubrication forces for Fast Lubrication with polydispersity
"mdpd"_pair_meso.html - mDPD particle interactions
"mdpd/rhosum"_pair_meso.html - mDPD particle interactions for mass density
"meam"_pair_meam.html - modified embedded atom method (MEAM) in Fortran
"meam/c"_pair_meam.html - modified embedded atom method (MEAM) in C
"meam/c"_pair_meamc.html - modified embedded atom method (MEAM) in C
"meam/spline"_pair_meam_spline.html - splined version of MEAM
"meam/sw/spline"_pair_meam_sw_spline.html - splined version of MEAM with a Stillinger-Weber term
"mgpt"_pair_mgpt.html - simplified model generalized pseudopotential theory (MGPT) potential
@ -260,7 +259,6 @@ accelerated styles exist.
"polymorphic"_pair_polymorphic.html - polymorphic 3-body potential
"python"_pair_python.html -
"quip"_pair_quip.html -
"reax"_pair_reax.html - ReaxFF potential in Fortran
"reax/c"_pair_reaxc.html - ReaxFF potential in C
"rebo"_pair_airebo.html - 2nd generation REBO potential of Brenner
"resquared"_pair_resquared.html - Everaers RE-Squared ellipsoidal potential

View File

@ -63,7 +63,7 @@ Pair Styles :h1
pair_lubricate
pair_lubricateU
pair_mdf
pair_meam
pair_meamc
pair_meam_spline
pair_meam_sw_spline
pair_meso
@ -82,7 +82,6 @@ Pair Styles :h1
pair_polymorphic
pair_python
pair_quip
pair_reax
pair_reaxc
pair_resquared
pair_sdk

View File

@ -1526,6 +1526,7 @@ mdf
mdpd
mDPD
meam
meamc
MEAMC
meamf
meanDist

View File

@ -1,30 +0,0 @@
# Test of MEAM potential for SiC system
units metal
boundary p p p
atom_style atomic
read_data data.meam
pair_style meam
pair_coeff * * library.meam Si C SiC.meam Si C
neighbor 0.3 bin
neigh_modify delay 10
fix 1 all nve
thermo 10
timestep 0.001
#dump 1 all atom 50 dump.meam
#dump 2 all image 10 image.*.jpg element element &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3 element Si C
#dump 3 all movie 10 movie.mpg element element &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3 element Si C
run 100

View File

@ -1,96 +0,0 @@
LAMMPS (27 Nov 2018)
using 1 OpenMP thread(s) per MPI task
# Test of MEAM potential for SiC system
units metal
boundary p p p
atom_style atomic
read_data data.meam
orthogonal box = (-6 -6 -6) to (5.97232 5.97232 5.97232)
1 by 1 by 1 MPI processor grid
reading atoms ...
128 atoms
pair_style meam
WARNING: THE pair_style meam COMMAND IS OBSOLETE AND WILL BE REMOVED VERY SOON. PLEASE USE meam/c (src/MEAM/pair_meam.cpp:51)
pair_coeff * * library.meam Si C SiC.meam Si C
Reading potential file library.meam with DATE: 2012-06-29
Reading potential file SiC.meam with DATE: 2007-06-11
neighbor 0.3 bin
neigh_modify delay 10
fix 1 all nve
thermo 10
timestep 0.001
#dump 1 all atom 50 dump.meam
#dump 2 all image 10 image.*.jpg element element # axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3 element Si C
#dump 3 all movie 10 movie.mpg element element # axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3 element Si C
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 = 4.3
ghost atom cutoff = 4.3
binsize = 2.15, bins = 6 6 6
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair meam, perpetual
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
(2) pair meam, perpetual, half/full from (1)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 8.103 | 8.103 | 8.103 Mbytes
Step Temp E_pair E_mol TotEng Press
0 0 -636.38121 0 -636.38121 -76571.819
10 1807.8862 -666.21959 0 -636.54126 -150571.49
20 1932.4467 -668.2581 0 -636.53498 -120223.52
30 1951.3652 -668.58139 0 -636.54771 -100508.4
40 2172.5974 -672.22715 0 -636.5617 -110753.34
50 2056.9149 -670.33108 0 -636.56468 -105418.07
60 1947.9564 -668.52788 0 -636.55015 -111413.04
70 1994.7712 -669.28849 0 -636.54225 -109645.76
80 2126.0903 -671.43755 0 -636.53557 -97475.831
90 2065.755 -670.4349 0 -636.52338 -95858.837
100 2051.4553 -670.20799 0 -636.53122 -107068.9
Loop time of 0.0545483 on 1 procs for 100 steps with 128 atoms
Performance: 158.392 ns/day, 0.152 hours/ns, 1833.239 timesteps/s
99.4% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.050821 | 0.050821 | 0.050821 | 0.0 | 93.17
Neigh | 0.0026484 | 0.0026484 | 0.0026484 | 0.0 | 4.86
Comm | 0.0006423 | 0.0006423 | 0.0006423 | 0.0 | 1.18
Output | 0.00011492 | 0.00011492 | 0.00011492 | 0.0 | 0.21
Modify | 0.00021195 | 0.00021195 | 0.00021195 | 0.0 | 0.39
Other | | 0.0001101 | | | 0.20
Nlocal: 128 ave 128 max 128 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 543 ave 543 max 543 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 1526 ave 1526 max 1526 min
Histogram: 1 0 0 0 0 0 0 0 0 0
FullNghs: 3052 ave 3052 max 3052 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 3052
Ave neighs/atom = 23.8438
Neighbor list builds = 10
Dangerous builds = 10
Total wall time: 0:00:00

View File

@ -1,96 +0,0 @@
LAMMPS (27 Nov 2018)
using 1 OpenMP thread(s) per MPI task
# Test of MEAM potential for SiC system
units metal
boundary p p p
atom_style atomic
read_data data.meam
orthogonal box = (-6 -6 -6) to (5.97232 5.97232 5.97232)
1 by 2 by 2 MPI processor grid
reading atoms ...
128 atoms
pair_style meam
WARNING: THE pair_style meam COMMAND IS OBSOLETE AND WILL BE REMOVED VERY SOON. PLEASE USE meam/c (src/MEAM/pair_meam.cpp:51)
pair_coeff * * library.meam Si C SiC.meam Si C
Reading potential file library.meam with DATE: 2012-06-29
Reading potential file SiC.meam with DATE: 2007-06-11
neighbor 0.3 bin
neigh_modify delay 10
fix 1 all nve
thermo 10
timestep 0.001
#dump 1 all atom 50 dump.meam
#dump 2 all image 10 image.*.jpg element element # axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3 element Si C
#dump 3 all movie 10 movie.mpg element element # axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3 element Si C
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 = 4.3
ghost atom cutoff = 4.3
binsize = 2.15, bins = 6 6 6
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair meam, perpetual
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
(2) pair meam, perpetual, half/full from (1)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 8.024 | 8.026 | 8.027 Mbytes
Step Temp E_pair E_mol TotEng Press
0 0 -636.38121 0 -636.38121 -76571.819
10 1807.8862 -666.21959 0 -636.54126 -150571.49
20 1932.4467 -668.2581 0 -636.53498 -120223.52
30 1951.3652 -668.58139 0 -636.54771 -100508.4
40 2172.5974 -672.22715 0 -636.5617 -110753.34
50 2056.9149 -670.33108 0 -636.56468 -105418.07
60 1947.9564 -668.52788 0 -636.55015 -111413.04
70 1994.7712 -669.28849 0 -636.54225 -109645.76
80 2126.0903 -671.43755 0 -636.53557 -97475.831
90 2065.755 -670.4349 0 -636.52338 -95858.837
100 2051.4553 -670.20799 0 -636.53122 -107068.9
Loop time of 0.023721 on 4 procs for 100 steps with 128 atoms
Performance: 364.234 ns/day, 0.066 hours/ns, 4215.667 timesteps/s
95.5% 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.019888 | 0.020242 | 0.020626 | 0.2 | 85.33
Neigh | 0.00071859 | 0.00076133 | 0.00082922 | 0.0 | 3.21
Comm | 0.0019681 | 0.0022618 | 0.002651 | 0.5 | 9.53
Output | 0.00018048 | 0.0002225 | 0.00034213 | 0.0 | 0.94
Modify | 8.2016e-05 | 8.6308e-05 | 9.203e-05 | 0.0 | 0.36
Other | | 0.000147 | | | 0.62
Nlocal: 32 ave 36 max 30 min
Histogram: 1 2 0 0 0 0 0 0 0 1
Nghost: 293.75 ave 305 max 285 min
Histogram: 2 0 0 0 0 0 0 1 0 1
Neighs: 381.5 ave 413 max 334 min
Histogram: 1 0 0 0 1 0 0 0 0 2
FullNghs: 763 ave 866 max 678 min
Histogram: 1 0 1 0 0 1 0 0 0 1
Total # of neighbors = 3052
Ave neighs/atom = 23.8438
Neighbor list builds = 10
Dangerous builds = 10
Total wall time: 0:00:00

View File

@ -1,209 +0,0 @@
LAMMPS (27 Nov 2018)
using 1 OpenMP thread(s) per MPI task
# 3d metal shear simulation
units metal
boundary s s p
atom_style atomic
lattice fcc 3.52
Lattice spacing in x,y,z = 3.52 3.52 3.52
region box block 0 16.0 0 10.0 0 2.828427
create_box 3 box
Created orthogonal box = (0 0 0) to (56.32 35.2 9.95606)
1 by 1 by 1 MPI processor grid
lattice fcc 3.52 orient x 1 0 0 orient y 0 1 1 orient z 0 -1 1 origin 0.5 0 0
Lattice spacing in x,y,z = 3.52 4.97803 4.97803
create_atoms 1 box
Created 1912 atoms
Time spent = 0.0007267 secs
pair_style meam
WARNING: THE pair_style meam COMMAND IS OBSOLETE AND WILL BE REMOVED VERY SOON. PLEASE USE meam/c (src/MEAM/pair_meam.cpp:51)
pair_coeff * * library.meam Ni4 Ni.meam Ni4 Ni4 Ni4
Reading potential file library.meam with DATE: 2012-06-29
Reading potential file Ni.meam with DATE: 2007-06-11
neighbor 0.3 bin
neigh_modify delay 5
region lower block INF INF INF 0.9 INF INF
region upper block INF INF 6.1 INF INF INF
group lower region lower
264 atoms in group lower
group upper region upper
264 atoms in group upper
group boundary union lower upper
528 atoms in group boundary
group mobile subtract all boundary
1384 atoms in group mobile
set group lower type 2
264 settings made for type
set group upper type 3
264 settings made for type
# void
#region void cylinder z 8 5 2.5 INF INF
#delete_atoms region void
# temp controllers
compute new3d mobile temp
compute new2d mobile temp/partial 0 1 1
# equilibrate
velocity mobile create 300.0 5812775 temp new3d
fix 1 all nve
fix 2 boundary setforce 0.0 0.0 0.0
fix 3 mobile temp/rescale 10 300.0 300.0 10.0 1.0
fix_modify 3 temp new3d
thermo 25
thermo_modify temp new3d
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:488)
timestep 0.001
run 100
Neighbor list info ...
update every 1 steps, delay 5 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 4.3
ghost atom cutoff = 4.3
binsize = 2.15, bins = 27 17 5
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair meam, perpetual
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
(2) pair meam, perpetual, half/full from (1)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 9.788 | 9.788 | 9.788 Mbytes
Step Temp E_pair E_mol TotEng Press Volume
0 300 -8232.7767 0 -8179.1466 1386.6643 19547.02
25 222.78953 -8188.1215 0 -8148.2941 9095.9008 19547.02
50 300 -8149.7654 0 -8096.1353 10633.141 19684.382
75 304.80657 -8163.4557 0 -8108.9665 7045.457 19759.745
100 300 -8173.6884 0 -8120.0584 5952.521 19886.589
Loop time of 0.894544 on 1 procs for 100 steps with 1912 atoms
Performance: 9.659 ns/day, 2.485 hours/ns, 111.789 timesteps/s
99.4% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.87012 | 0.87012 | 0.87012 | 0.0 | 97.27
Neigh | 0.01798 | 0.01798 | 0.01798 | 0.0 | 2.01
Comm | 0.0016143 | 0.0016143 | 0.0016143 | 0.0 | 0.18
Output | 0.00011492 | 0.00011492 | 0.00011492 | 0.0 | 0.01
Modify | 0.0035381 | 0.0035381 | 0.0035381 | 0.0 | 0.40
Other | | 0.001176 | | | 0.13
Nlocal: 1912 ave 1912 max 1912 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 1672 ave 1672 max 1672 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 23806 ave 23806 max 23806 min
Histogram: 1 0 0 0 0 0 0 0 0 0
FullNghs: 47612 ave 47612 max 47612 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 47612
Ave neighs/atom = 24.9017
Neighbor list builds = 5
Dangerous builds = 0
# shear
velocity upper set 1.0 0 0
velocity mobile ramp vx 0.0 1.0 y 1.4 8.6 sum yes
unfix 3
fix 3 mobile temp/rescale 10 300.0 300.0 10.0 1.0
fix_modify 3 temp new2d
#dump 1 all atom 500 dump.meam.shear
#dump 2 all image 100 image.*.jpg type type # axes yes 0.8 0.02 view 0 0 zoom 1.5 up 0 1 0 adiam 2.0
#dump_modify 2 pad 4
#dump 3 all movie 100 movie.mpg type type # axes yes 0.8 0.02 view 0 0 zoom 1.5 up 0 1 0 adiam 2.0
#dump_modify 3 pad 4
thermo 100
thermo_modify temp new2d
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:488)
reset_timestep 0
run 3000
Per MPI rank memory allocation (min/avg/max) = 9.964 | 9.964 | 9.964 Mbytes
Step Temp E_pair E_mol TotEng Press Volume
0 300.39988 -8173.6884 0 -8137.8874 4992.9811 19894.297
100 292.06374 -8177.7096 0 -8142.9021 2568.3762 19871.53
200 306.69894 -8177.1357 0 -8140.584 874.24259 20047.24
300 295.68229 -8172.9213 0 -8137.6825 -1049.0836 20091.759
400 308.99958 -8169.6355 0 -8132.8096 -1785.9335 20121.698
500 303.85723 -8163.984 0 -8127.7709 -150.56268 20183.813
600 300 -8157.7632 0 -8122.0099 1492.5742 20279.887
700 300 -8148.1328 0 -8112.3794 3506.9234 20435.302
800 300 -8139.1821 0 -8103.4288 3628.3957 20509.519
900 305.03425 -8126.7734 0 -8090.4201 5316.2206 20638.992
1000 304.00321 -8112.1616 0 -8075.9311 7441.9638 20767.243
1100 304.14041 -8096.5041 0 -8060.2573 9646.6972 20888.167
1200 302.78454 -8080.5931 0 -8044.5079 11516.208 20995.917
1300 308.67064 -8061.6724 0 -8024.8857 11496.471 21130.013
1400 309.82994 -8046.27 0 -8009.3451 12926.819 21247.271
1500 300 -8035.0317 0 -7999.2784 15346.797 21370.637
1600 300 -8030.6636 0 -7994.9102 14803.43 21496.446
1700 300 -8024.4819 0 -7988.7286 13175.257 21611.262
1800 300 -8022.8531 0 -7987.0998 10315.63 21743.178
1900 300 -8028.4095 0 -7992.6561 6882.0635 21855.551
2000 300 -8036.9005 0 -8001.1472 3508.9237 21983.802
2100 300 -8037.8224 0 -8002.0691 2724.0594 22112.054
2200 306.93248 -8035.3297 0 -7998.7501 4400.6008 22228.091
2300 306.24125 -8036.748 0 -8000.2508 6075.0546 22352.678
2400 300 -8038.8534 0 -8003.1 8701.8498 22465.051
2500 308.34129 -8034.0796 0 -7997.3322 10977.68 22600.632
2600 299.70072 -8028.8815 0 -7993.1638 15468.97 22715.447
2700 298.78276 -8019.1655 0 -7983.5572 18076.132 22844.921
2800 305.57845 -8014.3542 0 -7977.936 17573.035 22962.179
2900 300 -8015.7677 0 -7980.0144 13461.463 23087.988
3000 300 -8010.5908 0 -7974.8375 9333.4855 23199.139
Loop time of 29.4592 on 1 procs for 3000 steps with 1912 atoms
Performance: 8.799 ns/day, 2.728 hours/ns, 101.836 timesteps/s
99.8% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 28.489 | 28.489 | 28.489 | 0.0 | 96.71
Neigh | 0.77356 | 0.77356 | 0.77356 | 0.0 | 2.63
Comm | 0.052517 | 0.052517 | 0.052517 | 0.0 | 0.18
Output | 0.00095224 | 0.00095224 | 0.00095224 | 0.0 | 0.00
Modify | 0.10813 | 0.10813 | 0.10813 | 0.0 | 0.37
Other | | 0.03525 | | | 0.12
Nlocal: 1912 ave 1912 max 1912 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 1668 ave 1668 max 1668 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 23391 ave 23391 max 23391 min
Histogram: 1 0 0 0 0 0 0 0 0 0
FullNghs: 46782 ave 46782 max 46782 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 46782
Ave neighs/atom = 24.4676
Neighbor list builds = 219
Dangerous builds = 0
Total wall time: 0:00:30

View File

@ -1,209 +0,0 @@
LAMMPS (27 Nov 2018)
using 1 OpenMP thread(s) per MPI task
# 3d metal shear simulation
units metal
boundary s s p
atom_style atomic
lattice fcc 3.52
Lattice spacing in x,y,z = 3.52 3.52 3.52
region box block 0 16.0 0 10.0 0 2.828427
create_box 3 box
Created orthogonal box = (0 0 0) to (56.32 35.2 9.95606)
2 by 2 by 1 MPI processor grid
lattice fcc 3.52 orient x 1 0 0 orient y 0 1 1 orient z 0 -1 1 origin 0.5 0 0
Lattice spacing in x,y,z = 3.52 4.97803 4.97803
create_atoms 1 box
Created 1912 atoms
Time spent = 0.000408649 secs
pair_style meam
WARNING: THE pair_style meam COMMAND IS OBSOLETE AND WILL BE REMOVED VERY SOON. PLEASE USE meam/c (src/MEAM/pair_meam.cpp:51)
pair_coeff * * library.meam Ni4 Ni.meam Ni4 Ni4 Ni4
Reading potential file library.meam with DATE: 2012-06-29
Reading potential file Ni.meam with DATE: 2007-06-11
neighbor 0.3 bin
neigh_modify delay 5
region lower block INF INF INF 0.9 INF INF
region upper block INF INF 6.1 INF INF INF
group lower region lower
264 atoms in group lower
group upper region upper
264 atoms in group upper
group boundary union lower upper
528 atoms in group boundary
group mobile subtract all boundary
1384 atoms in group mobile
set group lower type 2
264 settings made for type
set group upper type 3
264 settings made for type
# void
#region void cylinder z 8 5 2.5 INF INF
#delete_atoms region void
# temp controllers
compute new3d mobile temp
compute new2d mobile temp/partial 0 1 1
# equilibrate
velocity mobile create 300.0 5812775 temp new3d
fix 1 all nve
fix 2 boundary setforce 0.0 0.0 0.0
fix 3 mobile temp/rescale 10 300.0 300.0 10.0 1.0
fix_modify 3 temp new3d
thermo 25
thermo_modify temp new3d
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:488)
timestep 0.001
run 100
Neighbor list info ...
update every 1 steps, delay 5 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 4.3
ghost atom cutoff = 4.3
binsize = 2.15, bins = 27 17 5
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair meam, perpetual
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
(2) pair meam, perpetual, half/full from (1)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 8.954 | 8.957 | 8.959 Mbytes
Step Temp E_pair E_mol TotEng Press Volume
0 300 -8232.7767 0 -8179.1466 1386.6643 19547.02
25 221.59546 -8187.6813 0 -8148.0673 9100.4509 19547.02
50 300 -8150.0685 0 -8096.4384 10317.407 19685.743
75 307.76021 -8164.6669 0 -8109.6496 6289.7138 19757.814
100 300 -8176.5141 0 -8122.884 4162.2559 19873.327
Loop time of 0.263516 on 4 procs for 100 steps with 1912 atoms
Performance: 32.787 ns/day, 0.732 hours/ns, 379.483 timesteps/s
98.3% 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.24401 | 0.2486 | 0.25128 | 0.6 | 94.34
Neigh | 0.0046518 | 0.0047416 | 0.0048261 | 0.1 | 1.80
Comm | 0.0054512 | 0.0082117 | 0.012793 | 3.1 | 3.12
Output | 0.00010562 | 0.00013095 | 0.00019932 | 0.0 | 0.05
Modify | 0.0010016 | 0.0010275 | 0.0010595 | 0.1 | 0.39
Other | | 0.0008045 | | | 0.31
Nlocal: 478 ave 492 max 465 min
Histogram: 2 0 0 0 0 0 0 0 1 1
Nghost: 809 ave 822 max 795 min
Histogram: 1 1 0 0 0 0 0 0 0 2
Neighs: 5916 ave 6133 max 5658 min
Histogram: 1 0 0 1 0 0 0 0 1 1
FullNghs: 11832 ave 12277 max 11299 min
Histogram: 1 0 0 1 0 0 0 0 1 1
Total # of neighbors = 47328
Ave neighs/atom = 24.7531
Neighbor list builds = 5
Dangerous builds = 0
# shear
velocity upper set 1.0 0 0
velocity mobile ramp vx 0.0 1.0 y 1.4 8.6 sum yes
unfix 3
fix 3 mobile temp/rescale 10 300.0 300.0 10.0 1.0
fix_modify 3 temp new2d
#dump 1 all atom 500 dump.meam.shear
#dump 2 all image 100 image.*.jpg type type # axes yes 0.8 0.02 view 0 0 zoom 1.5 up 0 1 0 adiam 2.0
#dump_modify 2 pad 4
#dump 3 all movie 100 movie.mpg type type # axes yes 0.8 0.02 view 0 0 zoom 1.5 up 0 1 0 adiam 2.0
#dump_modify 3 pad 4
thermo 100
thermo_modify temp new2d
WARNING: Temperature for thermo pressure is not for group all (src/thermo.cpp:488)
reset_timestep 0
run 3000
Per MPI rank memory allocation (min/avg/max) = 8.999 | 9.002 | 9.005 Mbytes
Step Temp E_pair E_mol TotEng Press Volume
0 295.32113 -8176.5141 0 -8141.3183 3169.3113 19886.93
100 292.00251 -8176.5358 0 -8141.7356 -825.04802 19918.765
200 306.11682 -8176.7719 0 -8140.2895 -1370.6886 19948.877
300 300 -8172.6262 0 -8136.8729 -1735.9765 20085.714
400 306.88489 -8168.435 0 -8131.8611 -933.02058 20117.012
500 308.99003 -8166.2906 0 -8129.4658 -1049.3138 20198.256
600 304.23435 -8158.0946 0 -8121.8366 583.93595 20328.848
700 296.44479 -8149.7914 0 -8114.4618 1985.4155 20421.046
800 307.75738 -8139.1649 0 -8102.487 4319.078 20513.183
900 304.61422 -8126.9246 0 -8090.6214 6654.0962 20640.213
1000 300 -8113.8464 0 -8078.0931 7760.1237 20768.465
1100 300.17873 -8097.7469 0 -8061.9722 8438.1259 20874.731
1200 306.01441 -8083.3367 0 -8046.8665 10835.588 20994.432
1300 300 -8067.022 0 -8031.2687 11216.061 21126.348
1400 300 -8053.223 0 -8017.4696 10570.211 21253.378
1500 300 -8043.4848 0 -8007.7314 11360.762 21375.523
1600 300 -8034.6201 0 -7998.8667 11371.282 21498.889
1700 300 -8028.6797 0 -7992.9263 9596.8454 21613.705
1800 300 -8033.0802 0 -7997.3268 8767.8176 21743.178
1900 303.23288 -8035.1821 0 -7999.0434 8065.2879 21859.215
2000 300 -8025.0795 0 -7989.3262 9321.8098 21980.138
2100 300 -8041.3621 0 -8005.6088 6674.2623 22108.39
2200 300 -8039.7261 0 -8003.9727 7548.8847 22225.648
2300 300 -8052.3497 0 -8016.5964 8936.4935 22352.678
2400 300 -8049.395 0 -8013.6416 12633.909 22476.044
2500 308.48099 -8039.9448 0 -8003.1807 16242.081 22593.303
2600 300 -8032.1953 0 -7996.442 18386.669 22722.776
2700 303.49413 -8027.6563 0 -7991.4865 14415.581 22829.042
2800 304.13476 -8017.3394 0 -7981.0933 7457.1076 22953.629
2900 300 -8010.3658 0 -7974.6124 2815.5155 23074.552
3000 309.49253 -7999.74 0 -7962.8553 756.7511 23210.132
Loop time of 8.57528 on 4 procs for 3000 steps with 1912 atoms
Performance: 30.226 ns/day, 0.794 hours/ns, 349.843 timesteps/s
98.5% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 8.0046 | 8.0451 | 8.1075 | 1.5 | 93.82
Neigh | 0.20142 | 0.20699 | 0.21469 | 1.2 | 2.41
Comm | 0.1972 | 0.2657 | 0.312 | 9.3 | 3.10
Output | 0.00087762 | 0.0015897 | 0.0037148 | 3.1 | 0.02
Modify | 0.030267 | 0.031624 | 0.032929 | 0.7 | 0.37
Other | | 0.02427 | | | 0.28
Nlocal: 478 ave 507 max 447 min
Histogram: 1 1 0 0 0 0 0 0 1 1
Nghost: 799.75 ave 842 max 763 min
Histogram: 2 0 0 0 0 0 0 1 0 1
Neighs: 5806.5 ave 6097 max 5507 min
Histogram: 1 0 1 0 0 0 0 1 0 1
FullNghs: 11613 ave 12159 max 11039 min
Histogram: 1 0 1 0 0 0 0 1 0 1
Total # of neighbors = 46452
Ave neighs/atom = 24.295
Neighbor list builds = 224
Dangerous builds = 0
Total wall time: 0:00:08

View File

@ -1,12 +1,8 @@
This directory contains input files for two short ReaxFF simulations
(RDX and TATB crystals) using the ReaxFF parameterization developed
for nitramines. The parameter file ffield.reax is the same as that in
subdirectory RDX (see below). Input files for both pair_style reax and
pair_style reax/c are provided.
In addition, each subdirectory contains validated parameter files for
a particular published version of ReaxFF. In all cases, the examples
use pair_style reax/c.
subdirectory RDX (see below). In addition, each subdirectory contains
validated parameter files for a particular published version of ReaxFF.
Disclaimer: Using these force fields for systems they have not been
explicitly trained against may produce unrealistic results. Please

View File

@ -1,52 +0,0 @@
# ReaxFF potential for RDX system
units real
atom_style charge
read_data data.rdx
# reax args: hbcut hbnewflag tripflag precision
pair_style reax 6.0 1 1 1.0e-6
pair_coeff * * ffield.reax 1 2 3 4
compute reax all pair reax
variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]
neighbor 2.5 bin
neigh_modify every 10 delay 0 check no
fix 1 all nve
thermo 10
thermo_style custom step temp epair etotal press &
v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb &
v_et v_eco v_ew v_ep v_efi v_eqeq
timestep 1.0
#dump 1 all custom 10 dump.reax.rdx id type q xs ys zs
#dump 2 all image 25 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
run 100

View File

@ -1,53 +0,0 @@
# ReaxFF potential for TATB system
units real
atom_style charge
read_data data.tatb
# reax args: hbcut hbnewflag tripflag precision
pair_style reax 6.0 1 1 1.0e-6
pair_coeff * * ffield.reax 1 2 3 4
compute reax all pair reax
variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]
neighbor 2.5 bin
neigh_modify delay 0 every 5 check no
fix 1 all nve
thermo 5
thermo_style custom step temp epair etotal press &
v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa &
v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
timestep 0.0625
#dump 1 all custom 100 dump.reax.tatb id type q x y z
#dump 2 all image 5 image.*.jpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 5 movie.mpg type type &
# axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
fix 2 all reax/bonds 25 bonds.reax.tatb
run 25

View File

@ -1,107 +0,0 @@
LAMMPS (8 Mar 2018)
using 1 OpenMP thread(s) per MPI task
# ReaxFF potential for RDX system
units real
atom_style charge
read_data data.rdx
orthogonal box = (35 35 35) to (48 48 48)
1 by 1 by 1 MPI processor grid
reading atoms ...
21 atoms
# reax args: hbcut hbnewflag tripflag precision
pair_style reax 6.0 1 1 1.0e-6
WARNING: The pair_style reax command is unsupported. Please switch to pair_style reax/c instead (../pair_reax.cpp:49)
pair_coeff * * ffield.reax 1 2 3 4
compute reax all pair reax
variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]
neighbor 2.5 bin
neigh_modify every 10 delay 0 check no
fix 1 all nve
thermo 10
thermo_style custom step temp epair etotal press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
timestep 1.0
#dump 1 all custom 10 dump.reax.rdx id type q xs ys zs
#dump 2 all image 25 image.*.jpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
run 100
Neighbor list info ...
update every 10 steps, delay 0 steps, check no
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12.5
ghost atom cutoff = 12.5
binsize = 6.25, bins = 3 3 3
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair reax, perpetual
attributes: half, newton off
pair build: half/bin/newtoff
stencil: half/bin/3d/newtoff
bin: standard
Per MPI rank memory allocation (min/avg/max) = 3.278 | 3.278 | 3.278 Mbytes
Step Temp E_pair TotEng Press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
0 0 -1885.1269 -1885.1269 27233.074 -2958.4712 79.527715 0.31082031 0 97.771125 25.846176 -0.18034154 0 16.709078 -9.1620736 938.43732 -244.79973 0 168.8842
10 1281.7558 -1989.1322 -1912.7188 -19609.913 -2733.8828 -15.775275 0.20055725 0 55.02023 3.1070523 -77.710916 0 14.963568 -5.8082203 843.41939 -180.17724 0 107.5115
20 516.83079 -1941.677 -1910.8655 -12525.412 -2801.8626 7.410797 0.073134186 0 81.986983 0.2281551 -57.494871 0 30.656735 -10.102557 877.78695 -158.93385 0 88.574159
30 467.26411 -1940.978 -1913.1215 -35957.489 -2755.021 -6.9179958 0.049322453 0 78.853173 0.13604393 -51.653635 0 19.862871 -9.7098575 853.79334 -151.232 0 80.86177
40 647.45528 -1951.1994 -1912.6006 -5883.713 -2798.3556 17.334814 0.15102862 0 63.235117 0.18070924 -54.598957 0 17.325007 -12.052278 883.0167 -164.21335 0 96.777424
50 716.38088 -1949.4735 -1906.7656 5473.1969 -2800.9309 9.2056861 0.15413274 0 85.371466 3.2986127 -78.253597 0 34.861774 -8.553123 882.01431 -193.85254 0 117.21068
60 1175.2705 -1975.961 -1905.8958 -1939.4966 -2726.5816 -11.651996 0.24296786 0 48.320654 7.1799691 -75.363638 0 16.520127 -4.8869441 844.75401 -194.23297 0 119.73841
70 1156.701 -1975.3497 -1906.3916 24628.304 -2880.5225 25.652501 0.26894311 0 83.724852 7.1049152 -68.70096 0 24.750735 -8.6338267 911.20079 -183.40562 0 113.21047
80 840.23677 -1955.4769 -1905.3851 -17731.334 -2755.7299 -8.0167723 0.1386797 0 86.147417 2.2387319 -76.945843 0 23.595869 -7.260968 853.63487 -167.88288 0 94.603961
90 365.79122 -1926.4061 -1904.599 898.38479 -2842.1832 47.368107 0.23109002 0 92.288071 0.38031213 -61.361485 0 18.476336 -12.25546 900.24233 -186.48046 0 116.88827
100 801.32158 -1953.418 -1905.6462 -2417.6887 -2802.7247 4.6676477 0.18046575 0 76.729987 5.4177322 -77.102566 0 24.997175 -7.7554074 898.67337 -196.89114 0 120.38946
Loop time of 0.463306 on 1 procs for 100 steps with 21 atoms
Performance: 18.649 ns/day, 1.287 hours/ns, 215.840 timesteps/s
99.6% CPU use with 1 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.46143 | 0.46143 | 0.46143 | 0.0 | 99.60
Neigh | 0.00087953 | 0.00087953 | 0.00087953 | 0.0 | 0.19
Comm | 0.00042653 | 0.00042653 | 0.00042653 | 0.0 | 0.09
Output | 0.00034237 | 0.00034237 | 0.00034237 | 0.0 | 0.07
Modify | 0.00010109 | 0.00010109 | 0.00010109 | 0.0 | 0.02
Other | | 0.000124 | | | 0.03
Nlocal: 21 ave 21 max 21 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 546 ave 546 max 546 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 1106 ave 1106 max 1106 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 1106
Ave neighs/atom = 52.6667
Neighbor list builds = 10
Dangerous builds not checked
Total wall time: 0:00:00

View File

@ -1,107 +0,0 @@
LAMMPS (8 Mar 2018)
using 1 OpenMP thread(s) per MPI task
# ReaxFF potential for RDX system
units real
atom_style charge
read_data data.rdx
orthogonal box = (35 35 35) to (48 48 48)
1 by 2 by 2 MPI processor grid
reading atoms ...
21 atoms
# reax args: hbcut hbnewflag tripflag precision
pair_style reax 6.0 1 1 1.0e-6
WARNING: The pair_style reax command is unsupported. Please switch to pair_style reax/c instead (../pair_reax.cpp:49)
pair_coeff * * ffield.reax 1 2 3 4
compute reax all pair reax
variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]
neighbor 2.5 bin
neigh_modify every 10 delay 0 check no
fix 1 all nve
thermo 10
thermo_style custom step temp epair etotal press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
timestep 1.0
#dump 1 all custom 10 dump.reax.rdx id type q xs ys zs
#dump 2 all image 25 image.*.jpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 25 movie.mpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
run 100
Neighbor list info ...
update every 10 steps, delay 0 steps, check no
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12.5
ghost atom cutoff = 12.5
binsize = 6.25, bins = 3 3 3
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair reax, perpetual
attributes: half, newton off
pair build: half/bin/newtoff
stencil: half/bin/3d/newtoff
bin: standard
Per MPI rank memory allocation (min/avg/max) = 3.262 | 3.36 | 3.647 Mbytes
Step Temp E_pair TotEng Press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
0 0 -1885.1268 -1885.1268 27233.074 -2958.4712 79.527715 0.31082031 0 97.771125 25.846176 -0.18034154 0 16.709078 -9.1620736 938.43732 -244.79972 0 168.88428
10 1281.7558 -1989.1322 -1912.7187 -19609.913 -2733.8828 -15.775275 0.20055725 0 55.020231 3.1070523 -77.710916 0 14.963568 -5.8082203 843.41939 -180.17724 0 107.51152
20 516.83079 -1941.677 -1910.8655 -12525.412 -2801.8626 7.410797 0.073134187 0 81.986983 0.2281551 -57.494871 0 30.656735 -10.102557 877.78695 -158.93385 0 88.574168
30 467.26411 -1940.978 -1913.1215 -35957.489 -2755.021 -6.9179959 0.049322449 0 78.853173 0.13604392 -51.653635 0 19.862871 -9.7098575 853.79334 -151.232 0 80.861765
40 647.45479 -1951.1995 -1912.6007 -5883.7199 -2798.3556 17.334805 0.15102868 0 63.235116 0.18070946 -54.59897 0 17.32501 -12.052277 883.0166 -164.21339 0 96.777473
50 716.37927 -1949.466 -1906.7582 5473.2486 -2800.9309 9.2056758 0.15413278 0 85.37143 3.2986099 -78.253596 0 34.861773 -8.5531243 882.01424 -193.85223 0 117.21791
60 1175.2698 -1975.9612 -1905.896 -1939.5206 -2726.5818 -11.651942 0.24296793 0 48.320679 7.1799538 -75.36365 0 16.520134 -4.8869515 844.75405 -194.23289 0 119.7383
70 1156.6963 -1975.3494 -1906.3915 24628.423 -2880.5221 25.65242 0.26894312 0 83.724787 7.1049615 -68.700925 0 24.750729 -8.6338123 911.2006 -183.40591 0 113.21091
80 840.238 -1955.4788 -1905.387 -17731.371 -2755.7301 -8.0167357 0.13868007 0 86.147246 2.2387405 -76.945868 0 23.595868 -7.2609697 853.6349 -167.88312 0 94.602512
90 365.78645 -1926.4072 -1904.6004 898.36945 -2842.1831 47.368307 0.23108998 0 92.288039 0.38031101 -61.361464 0 18.476388 -12.255481 900.24216 -186.48066 0 116.88716
100 801.31322 -1953.4165 -1905.6452 -2417.2041 -2802.7247 4.6678077 0.18046498 0 76.730367 5.4176812 -77.102592 0 24.9973 -7.7554425 898.6732 -196.89097 0 120.39043
Loop time of 0.404551 on 4 procs for 100 steps with 21 atoms
Performance: 21.357 ns/day, 1.124 hours/ns, 247.188 timesteps/s
97.4% 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.2191 | 0.28038 | 0.39839 | 13.2 | 69.31
Neigh | 5.8651e-05 | 0.00025928 | 0.00062203 | 0.0 | 0.06
Comm | 0.0046599 | 0.12307 | 0.1845 | 19.9 | 30.42
Output | 0.00055337 | 0.00062728 | 0.00071192 | 0.0 | 0.16
Modify | 5.3167e-05 | 7.844e-05 | 0.00010109 | 0.0 | 0.02
Other | | 0.0001363 | | | 0.03
Nlocal: 5.25 ave 15 max 0 min
Histogram: 1 0 2 0 0 0 0 0 0 1
Nghost: 355.5 ave 432 max 282 min
Histogram: 1 0 0 0 1 1 0 0 0 1
Neighs: 301.25 ave 827 max 0 min
Histogram: 1 0 2 0 0 0 0 0 0 1
Total # of neighbors = 1205
Ave neighs/atom = 57.381
Neighbor list builds = 10
Dangerous builds not checked
Total wall time: 0:00:00

View File

@ -1,103 +0,0 @@
LAMMPS (8 Mar 2018)
using 1 OpenMP thread(s) per MPI task
# ReaxFF potential for TATB system
units real
atom_style charge
read_data data.tatb
triclinic box = (0 0 0) to (13.624 17.1149 15.1826) with tilt (-5.75316 -6.32547 7.42573)
1 by 1 by 1 MPI processor grid
reading atoms ...
384 atoms
# reax args: hbcut hbnewflag tripflag precision
pair_style reax 6.0 1 1 1.0e-6
WARNING: The pair_style reax command is unsupported. Please switch to pair_style reax/c instead (../pair_reax.cpp:49)
pair_coeff * * ffield.reax 1 2 3 4
compute reax all pair reax
variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]
neighbor 2.5 bin
neigh_modify delay 0 every 5 check no
fix 1 all nve
thermo 5
thermo_style custom step temp epair etotal press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
timestep 0.0625
#dump 1 all custom 100 dump.reax.tatb id type q x y z
#dump 2 all image 5 image.*.jpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 5 movie.mpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
fix 2 all reax/bonds 25 bonds.reax.tatb
run 25
Neighbor list info ...
update every 5 steps, delay 0 steps, check no
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12.5
ghost atom cutoff = 12.5
binsize = 6.25, bins = 5 4 3
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair reax, perpetual
attributes: half, newton off
pair build: half/bin/newtoff
stencil: half/bin/3d/newtoff
bin: standard
Per MPI rank memory allocation (min/avg/max) = 7.764 | 7.764 | 7.764 Mbytes
Step Temp E_pair TotEng Press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
0 0 -44767.08 -44767.08 7294.6353 -61120.591 486.4378 4.7236377 0 1568.024 20.788929 -279.51642 -1556.4696 252.57147 -655.84699 18862.412 -8740.6378 0 6391.0231
5 0.63682806 -44767.737 -44767.01 8391.5964 -61118.763 486.82916 4.723415 0 1567.835 20.768662 -278.20804 -1557.6962 252.64683 -655.74117 18859.328 -8738.2728 0 6388.8127
10 2.4306958 -44769.409 -44766.634 11717.376 -61113.142 487.89093 4.7227063 0 1567.2936 20.705084 -274.37509 -1560.8546 252.87219 -655.43578 18850.19 -8731.0693 0 6381.7942
15 5.0590493 -44772.631 -44766.855 17125.067 -61103.34 489.28007 4.7214008 0 1566.4744 20.590604 -268.28962 -1566.5961 252.97781 -654.93836 18835.335 -8719.3013 0 6370.4551
20 8.067859 -44775.936 -44766.725 24620.627 -61088.791 490.42346 4.7193467 0 1565.5541 20.415031 -260.38512 -1574.1001 253.39805 -654.26837 18815.312 -8703.3748 0 6355.1614
25 10.975538 -44777.233 -44764.702 34381.173 -61068.889 490.53149 4.7164093 0 1566.5715 20.169755 -251.23109 -1582.8552 253.88696 -653.46042 18790.855 -8683.8691 0 6336.3409
Loop time of 7.80129 on 1 procs for 25 steps with 384 atoms
Performance: 0.017 ns/day, 1386.896 hours/ns, 3.205 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 | 7.7384 | 7.7384 | 7.7384 | 0.0 | 99.19
Neigh | 0.058615 | 0.058615 | 0.058615 | 0.0 | 0.75
Comm | 0.0022428 | 0.0022428 | 0.0022428 | 0.0 | 0.03
Output | 0.00033212 | 0.00033212 | 0.00033212 | 0.0 | 0.00
Modify | 0.0013618 | 0.0013618 | 0.0013618 | 0.0 | 0.02
Other | | 0.0003309 | | | 0.00
Nlocal: 384 ave 384 max 384 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 7559 ave 7559 max 7559 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 286828 ave 286828 max 286828 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 286828
Ave neighs/atom = 746.948
Neighbor list builds = 5
Dangerous builds not checked
Total wall time: 0:00:08

View File

@ -1,103 +0,0 @@
LAMMPS (8 Mar 2018)
using 1 OpenMP thread(s) per MPI task
# ReaxFF potential for TATB system
units real
atom_style charge
read_data data.tatb
triclinic box = (0 0 0) to (13.624 17.1149 15.1826) with tilt (-5.75316 -6.32547 7.42573)
1 by 2 by 2 MPI processor grid
reading atoms ...
384 atoms
# reax args: hbcut hbnewflag tripflag precision
pair_style reax 6.0 1 1 1.0e-6
WARNING: The pair_style reax command is unsupported. Please switch to pair_style reax/c instead (../pair_reax.cpp:49)
pair_coeff * * ffield.reax 1 2 3 4
compute reax all pair reax
variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]
neighbor 2.5 bin
neigh_modify delay 0 every 5 check no
fix 1 all nve
thermo 5
thermo_style custom step temp epair etotal press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
timestep 0.0625
#dump 1 all custom 100 dump.reax.tatb id type q x y z
#dump 2 all image 5 image.*.jpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 2 pad 3
#dump 3 all movie 5 movie.mpg type type # axes yes 0.8 0.02 view 60 -30
#dump_modify 3 pad 3
fix 2 all reax/bonds 25 bonds.reax.tatb
run 25
Neighbor list info ...
update every 5 steps, delay 0 steps, check no
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12.5
ghost atom cutoff = 12.5
binsize = 6.25, bins = 5 4 3
1 neighbor lists, perpetual/occasional/extra = 1 0 0
(1) pair reax, perpetual
attributes: half, newton off
pair build: half/bin/newtoff
stencil: half/bin/3d/newtoff
bin: standard
Per MPI rank memory allocation (min/avg/max) = 4.402 | 4.402 | 4.402 Mbytes
Step Temp E_pair TotEng Press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq
0 0 -44767.08 -44767.08 7294.6353 -61120.591 486.4378 4.7236377 0 1568.024 20.788929 -279.51642 -1556.4696 252.57147 -655.84699 18862.412 -8740.6378 0 6391.0231
5 0.63682727 -44767.816 -44767.089 8391.1708 -61118.763 486.82916 4.723415 0 1567.835 20.768662 -278.20804 -1557.6962 252.64683 -655.74117 18859.328 -8738.3973 0 6388.8581
10 2.4306941 -44769.405 -44766.63 11717.306 -61113.142 487.89094 4.7227063 0 1567.2936 20.705084 -274.3751 -1560.8546 252.87219 -655.43578 18850.19 -8731.08 0 6381.8083
15 5.0590444 -44772.6 -44766.824 17125.207 -61103.34 489.28008 4.7214008 0 1566.4744 20.590604 -268.28963 -1566.5961 252.97781 -654.93836 18835.335 -8719.2653 0 6370.4505
20 8.0678523 -44775.983 -44766.772 24620.114 -61088.791 490.42348 4.7193467 0 1565.5541 20.415031 -260.38513 -1574.1001 253.39804 -654.26837 18815.312 -8703.5228 0 6355.2629
25 10.975532 -44777.234 -44764.704 34381.065 -61068.889 490.53151 4.7164093 0 1566.5715 20.169755 -251.23111 -1582.8552 253.88696 -653.46042 18790.855 -8683.898 0 6336.3682
Loop time of 3.74388 on 4 procs for 25 steps with 384 atoms
Performance: 0.036 ns/day, 665.579 hours/ns, 6.678 timesteps/s
98.7% CPU use with 4 MPI tasks x 1 OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 3.478 | 3.6025 | 3.7215 | 4.8 | 96.22
Neigh | 0.012731 | 0.01299 | 0.013174 | 0.2 | 0.35
Comm | 0.0073411 | 0.12653 | 0.25119 | 25.4 | 3.38
Output | 0.00050354 | 0.00081849 | 0.0011628 | 0.0 | 0.02
Modify | 0.00049281 | 0.00082356 | 0.001157 | 0.0 | 0.02
Other | | 0.0002663 | | | 0.01
Nlocal: 96 ave 96 max 96 min
Histogram: 4 0 0 0 0 0 0 0 0 0
Nghost: 5118 ave 5118 max 5118 min
Histogram: 4 0 0 0 0 0 0 0 0 0
Neighs: 79754 ave 79754 max 79754 min
Histogram: 4 0 0 0 0 0 0 0 0 0
Total # of neighbors = 319016
Ave neighs/atom = 830.771
Neighbor list builds = 5
Dangerous builds not checked
Total wall time: 0:00:03

View File

@ -33,8 +33,6 @@ kokkos Kokkos package for GPU and many-core acceleration
from Kokkos development team (Sandia)
linalg set of BLAS and LAPACK routines needed by USER-ATC package
from Axel Kohlmeyer (Temple U)
meam modified embedded atom method (MEAM) potential, MEAM package
from Greg Wagner (Sandia)
message client/server communication library via MPI, sockets, files
from Steve Plimpton (Sandia)
molfile hooks to VMD molfile plugins, used by the USER-MOLFILE package
@ -51,8 +49,6 @@ qmmm quantum mechanics/molecular mechanics coupling interface
from Axel Kohlmeyer (Temple U)
quip interface to QUIP/libAtoms framework, USER-QUIP package
from Albert Bartok-Partay and Gabor Csanyi (U Cambridge)
reax ReaxFF potential, REAX package
from Adri van Duin (Penn State) and Aidan Thompson (Sandia)
smd hooks to Eigen library, used by USER-SMD package
from Georg Ganzenmueller (Ernst Mach Institute, Germany)
voronoi hooks to the Voro++ library, used by compute voronoi/atom command

View File

@ -1,9 +0,0 @@
# dependencies. needed for parallel make
$(DIR)meam_data.o: meam_data.F
$(DIR)meam_cleanup.o: meam_cleanup.F $(DIR)meam_data.o
$(DIR)meam_dens_final.o: meam_dens_final.F $(DIR)meam_data.o
$(DIR)meam_dens_init.o: meam_dens_init.F $(DIR)meam_data.o
$(DIR)meam_force.o: meam_force.F $(DIR)meam_data.o
$(DIR)meam_setup_done.o: meam_setup_done.F $(DIR)meam_data.o
$(DIR)meam_setup_global.o: meam_setup_global.F $(DIR)meam_data.o
$(DIR)meam_setup_param.o: meam_setup_param.F $(DIR)meam_data.o

1
lib/meam/.gitignore vendored
View File

@ -1 +0,0 @@
*.mod

View File

@ -1 +0,0 @@
../Install.py

View File

@ -1,57 +0,0 @@
# *
# *_________________________________________________________________________*
# * MEAM: MODEFIED EMBEDDED ATOM METHOD *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * AUTHORS: Greg Wagner, Sandia National Laboratories *
# * CONTACT: gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.gfortran
# ------ FILES ------
SRC = meam_data.F meam_setup_done.F meam_setup_global.F meam_setup_param.F meam_dens_init.F meam_dens_final.F meam_force.F meam_cleanup.F
FILES = $(SRC) Makefile
# ------ DEFINITIONS ------
LIB = libmeam.a
OBJ = $(SRC:.F=.o) fm_exp.o
# ------ SETTINGS ------
F90 = g95
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F
$(F90) $(F90FLAGS) -c $<
%.o:%.c
$(CC) $(F90FLAGS) -c $<
include .depend
# ------ CLEAN ------
clean:
-rm *.o *.mod *~ $(LIB)
tar:
-tar -cvf ../MEAM.tar $(FILES)

View File

@ -1,61 +0,0 @@
# *
# *_________________________________________________________________________*
# * MEAM: MODEFIED EMBEDDED ATOM METHOD *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * AUTHORS: Greg Wagner, Sandia National Laboratories *
# * CONTACT: gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.gfortran
# ------ FILES ------
SRC = meam_data.F meam_setup_done.F meam_setup_global.F meam_setup_param.F meam_dens_init.F meam_dens_final.F meam_force.F meam_cleanup.F
FILES = $(SRC) Makefile
# ------ DEFINITIONS ------
LIB = libmeam.a
OBJ = $(SRC:.F=.o) fm_exp.o
# ------ SETTINGS ------
F90 = gfortran
CC = gcc
F90FLAGS = -O3 -fPIC -ffast-math -ftree-vectorize -fexpensive-optimizations -fno-second-underscore -g
#F90FLAGS = -O
ARCHIVE = ar
ARCHFLAG = -rc
LINK = g++
LINKFLAGS = -O
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F
$(F90) $(F90FLAGS) -c $<
%.o:%.c
$(CC) $(F90FLAGS) -c $<
include .depend
# ------ CLEAN ------
clean:
-rm *.o *.mod *~ $(LIB)
tar:
-tar -cvf ../MEAM.tar $(FILES)

View File

@ -1,57 +0,0 @@
# *
# *_________________________________________________________________________*
# * MEAM: MODEFIED EMBEDDED ATOM METHOD *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * AUTHORS: Greg Wagner, Sandia National Laboratories *
# * CONTACT: gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.ifort
# ------ FILES ------
SRC = meam_data.F meam_setup_done.F meam_setup_global.F meam_setup_param.F meam_dens_init.F meam_dens_final.F meam_force.F meam_cleanup.F
FILES = $(SRC) Makefile
# ------ DEFINITIONS ------
LIB = libmeam.a
OBJ = $(SRC:.F=.o) fm_exp.o
# ------ SETTINGS ------
F90 = ifort
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F
$(F90) $(F90FLAGS) -c $<
%.o:%.c
$(CC) $(F90FLAGS) -c $<
include .depend
# ------ CLEAN ------
clean:
-rm *.o *.mod *~ $(LIB)
tar:
-tar -cvf ../MEAM.tar $(FILES)

View File

@ -1,5 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
meam_SYSINC =
meam_SYSLIB =
meam_SYSPATH =

View File

@ -1,5 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
meam_SYSINC =
meam_SYSLIB = -lgfortran
meam_SYSPATH =

View File

@ -1,5 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
meam_SYSINC =
meam_SYSLIB = -lifcore -lsvml -lompstub -limf
meam_SYSPATH = -L/opt/intel-11.1.046/lib/intel64

View File

@ -1,5 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
meam_SYSINC =
meam_SYSLIB = -lifcore -lsvml -lompstub -limf
meam_SYSPATH = -L/opt/intel/fce/10.0.023/lib

View File

@ -1,61 +0,0 @@
# *
# *_________________________________________________________________________*
# * MEAM: MODEFIED EMBEDDED ATOM METHOD *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * AUTHORS: Greg Wagner, Sandia National Laboratories *
# * CONTACT: gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.empty
# ------ FILES ------
SRC = meam_data.F meam_setup_done.F meam_setup_global.F meam_setup_param.F meam_dens_init.F meam_dens_final.F meam_force.F meam_cleanup.F
FILES = $(SRC) Makefile
# ------ DEFINITIONS ------
LIB = libmeam.a
OBJ = $(SRC:.F=.o) fm_exp.o
# ------ SETTINGS ------
F90 = mpifort
CC = mpicc
F90FLAGS = -O3 -fPIC
#F90FLAGS = -O
ARCHIVE = ar
ARCHFLAG = -rc
LINK = mpicxx
LINKFLAGS = -O
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F
$(F90) $(F90FLAGS) -c $<
%.o:%.c
$(CC) $(F90FLAGS) -c $<
include .depend
# ------ CLEAN ------
clean:
-rm *.o *.mod *~ $(LIB)
tar:
-tar -cvf ../MEAM.tar $(FILES)

View File

@ -1,57 +0,0 @@
# *
# *_________________________________________________________________________*
# * MEAM: MODEFIED EMBEDDED ATOM METHOD *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * AUTHORS: Greg Wagner, Sandia National Laboratories *
# * CONTACT: gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.pgf90
# ------ FILES ------
SRC = meam_data.F meam_setup_done.F meam_setup_global.F meam_setup_param.F meam_dens_init.F meam_dens_final.F meam_force.F meam_cleanup.F
FILES = $(SRC) Makefile
# ------ DEFINITIONS ------
LIB = libmeam.a
OBJ = $(SRC:.F=.o) fm_exp.o
# ------ SETTINGS ------
F90 = pgf90
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F
$(F90) $(F90FLAGS) -c $<
%.o:%.c
$(CC) $(F90FLAGS) -c $<
include .depend
# ------ CLEAN ------
clean:
-rm *.o *.mod *~ $(LIB)
tar:
-tar -cvf ../MEAM.tar $(FILES)

View File

@ -1 +0,0 @@
Makefile.gfortran

View File

@ -1,59 +0,0 @@
# *
# *_________________________________________________________________________*
# * MEAM: MODEFIED EMBEDDED ATOM METHOD *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * AUTHORS: Greg Wagner, Sandia National Laboratories *
# * CONTACT: gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.glory
# ------ FILES ------
SRC = meam_data.F meam_setup_done.F meam_setup_global.F meam_setup_param.F meam_dens_init.F meam_dens_final.F meam_force.F meam_cleanup.F
FILES = $(SRC) Makefile
# ------ DEFINITIONS ------
LIB = libmeam.a
OBJ = $(SRC:.F=.o) fm_exp.o
# ------ SETTINGS ------
F90 = mpif90
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
LINK = g++
LINKFLAGS = -O
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F
$(F90) $(F90FLAGS) -c $<
%.o:%.c
$(CC) $(F90FLAGS) -c $<
include .depend
# ------ CLEAN ------
clean:
-rm *.o *.mod *~ $(LIB)
tar:
-tar -cvf ../MEAM.tar $(FILES)

View File

@ -1,51 +0,0 @@
MEAM (modified embedded atom method) library
Greg Wagner, Sandia National Labs
gjwagne at sandia.gov
Jan 2007
This library is in implementation of the MEAM potential, specifically
designed to work with LAMMPS.
-------------------------------------------------
This directory has source files to build a library that LAMMPS
links against when using the MEAM package.
This library must be built with a F90 compiler, before LAMMPS is
built, so LAMMPS can link against it.
You can type "make lib-meam" from the src directory to see help on how
to build this library via make commands, or you can do the same thing
by typing "python Install.py" from within this directory, or you can
do it manually by following the instructions below.
Build the library using one of the provided Makefile.* files or create
your own, specific to your compiler and system. For example:
make -f Makefile.gfortran
When you are done building this library, two files should
exist in this directory:
libmeam.a the library LAMMPS will link against
Makefile.lammps settings the LAMMPS Makefile will import
Makefile.lammps is created by the make command, by copying one of the
Makefile.lammps.* files. See the EXTRAMAKE setting at the top of the
Makefile.* files.
IMPORTANT: You must examine the final Makefile.lammps to insure it is
correct for your system, else the LAMMPS build will likely fail.
Makefile.lammps has settings for 3 variables:
user-meam_SYSINC = leave blank for this package
user-meam_SYSLIB = auxiliary F90 libs needed to link a F90 lib with
a C++ program (LAMMPS) via a C++ compiler
user-meam_SYSPATH = path(s) to where those libraries are
Because you have a F90 compiler on your system, you should have these
libraries. But you will have to figure out which ones are needed and
where they are. Examples of common configurations are in the
Makefile.lammps.* files.

View File

@ -1,133 +0,0 @@
/*
Copyright (c) 2012,2013 Axel Kohlmeyer <akohlmey@gmail.com>
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above copyright
notice, this list of conditions and the following disclaimer in the
documentation and/or other materials provided with the distribution.
* Neither the name of the <organization> nor the
names of its contributors may be used to endorse or promote products
derived from this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
ARE DISCLAIMED. IN NO EVENT SHALL <COPYRIGHT HOLDER> BE LIABLE FOR ANY
DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/
/* faster versions of 2**x, e**x, and 10**x in single and double precision.
*
* Based on the Cephes math library 2.8
*/
#include <math.h>
#include <stdint.h>
/* internal definitions for the fastermath library */
/* IEEE 754 double precision floating point data manipulation */
typedef union
{
double f;
uint64_t u;
struct {int32_t i0,i1;};
} udi_t;
#define FM_DOUBLE_BIAS 1023
#define FM_DOUBLE_EMASK 2146435072
#define FM_DOUBLE_MBITS 20
#define FM_DOUBLE_MMASK 1048575
#define FM_DOUBLE_EZERO 1072693248
/* generate 2**num in floating point by bitshifting */
#define FM_DOUBLE_INIT_EXP(var,num) \
var.i0 = 0; \
var.i1 = (((int) num) + FM_DOUBLE_BIAS) << 20
/* double precision constants */
#define FM_DOUBLE_LOG2OFE 1.4426950408889634074
#define FM_DOUBLE_LOGEOF2 6.9314718055994530942e-1
#define FM_DOUBLE_LOG2OF10 3.32192809488736234789
#define FM_DOUBLE_LOG10OF2 3.0102999566398119521e-1
#define FM_DOUBLE_LOG10OFE 4.3429448190325182765e-1
#define FM_DOUBLE_SQRT2 1.41421356237309504880
#define FM_DOUBLE_SQRTH 0.70710678118654752440
/* optimizer friendly implementation of exp2(x).
*
* strategy:
*
* split argument into an integer part and a fraction:
* ipart = floor(x+0.5);
* fpart = x - ipart;
*
* compute exp2(ipart) from setting the ieee754 exponent
* compute exp2(fpart) using a pade' approximation for x in [-0.5;0.5[
*
* the result becomes: exp2(x) = exp2(ipart) * exp2(fpart)
*/
static const double fm_exp2_q[] = {
/* 1.00000000000000000000e0, */
2.33184211722314911771e2,
4.36821166879210612817e3
};
static const double fm_exp2_p[] = {
2.30933477057345225087e-2,
2.02020656693165307700e1,
1.51390680115615096133e3
};
static double fm_exp2(double x)
{
double ipart, fpart, px, qx;
udi_t epart;
ipart = floor(x+0.5);
fpart = x - ipart;
FM_DOUBLE_INIT_EXP(epart,ipart);
x = fpart*fpart;
px = fm_exp2_p[0];
px = px*x + fm_exp2_p[1];
qx = x + fm_exp2_q[0];
px = px*x + fm_exp2_p[2];
qx = qx*x + fm_exp2_q[1];
px = px * fpart;
x = 1.0 + 2.0*(px/(qx-px));
return epart.f*x;
}
double fm_exp_(double *x)
{
#if defined(__BYTE_ORDER__)
#if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
return fm_exp2(FM_DOUBLE_LOG2OFE * (*x));
#endif
#endif
return exp(*x);
}
/*
* Local Variables:
* mode: c
* compile-command: "make -C .."
* c-basic-offset: 4
* fill-column: 76
* indent-tabs-mode: nil
* End:
*/

View File

@ -1,26 +0,0 @@
c Declaration in pair_meam.h:
c
c void meam_cleanup()
c
c Call from PairMEAM destructor
c
c meam_cleanup()
c
subroutine meam_cleanup
use meam_data
implicit none
integer dealloc_error
deallocate(phir,STAT=dealloc_error)
deallocate(phirar,STAT=dealloc_error)
deallocate(phirar1,STAT=dealloc_error)
deallocate(phirar2,STAT=dealloc_error)
deallocate(phirar3,STAT=dealloc_error)
deallocate(phirar4,STAT=dealloc_error)
deallocate(phirar5,STAT=dealloc_error)
deallocate(phirar6,STAT=dealloc_error)
return
end

View File

@ -1,87 +0,0 @@
module meam_data
integer, parameter :: maxelt = 5
real*8 , external :: fm_exp
c cutforce = force cutoff
c cutforcesq = force cutoff squared
real*8 cutforce,cutforcesq
c Ec_meam = cohesive energy
c re_meam = nearest-neighbor distance
c Omega_meam = atomic volume
c B_meam = bulk modulus
c Z_meam = number of first neighbors for reference structure
c ielt_meam = atomic number of element
c A_meam = adjustable parameter
c alpha_meam = sqrt(9*Omega*B/Ec)
c rho0_meam = density scaling parameter
c delta_meam = heat of formation for alloys
c beta[0-3]_meam = electron density constants
c t[0-3]_meam = coefficients on densities in Gamma computation
c rho_ref_meam = background density for reference structure
c ibar_meam(i) = selection parameter for Gamma function for elt i,
c lattce_meam(i,j) = lattce configuration for elt i or alloy (i,j)
c neltypes = maximum number of element type defined
c eltind = index number of pair (similar to Voigt notation; ij = ji)
c phir = pair potential function array
c phirar[1-6] = spline coeffs
c attrac_meam = attraction parameter in Rose energy
c repuls_meam = repulsion parameter in Rose energy
c nn2_meam = 1 if second nearest neighbors are to be computed, else 0
c zbl_meam = 1 if zbl potential for small r to be use, else 0
c emb_lin_neg = 1 if linear embedding function for rhob to be used, else 0
c bkgd_dyn = 1 if reference densities follows Dynamo, else 0
c Cmin_meam, Cmax_meam = min and max values in screening cutoff
c rc_meam = cutoff distance for meam
c delr_meam = cutoff region for meam
c ebound_meam = factor giving maximum boundary of sceen fcn ellipse
c augt1 = flag for whether t1 coefficient should be augmented
c ialloy = flag for newer alloy formulation (as in dynamo code)
c mix_ref_t = flag to recover "old" way of computing t in reference config
c erose_form = selection parameter for form of E_rose function
c gsmooth_factor = factor determining length of G smoothing region
c vind[23]D = Voight notation index maps for 2 and 3D
c v2D,v3D = array of factors to apply for Voight notation
c nr,dr = pair function discretization parameters
c nrar,rdrar = spline coeff array parameters
real*8 Ec_meam(maxelt,maxelt),re_meam(maxelt,maxelt)
real*8 Omega_meam(maxelt),Z_meam(maxelt)
real*8 A_meam(maxelt),alpha_meam(maxelt,maxelt),rho0_meam(maxelt)
real*8 delta_meam(maxelt,maxelt)
real*8 beta0_meam(maxelt),beta1_meam(maxelt)
real*8 beta2_meam(maxelt),beta3_meam(maxelt)
real*8 t0_meam(maxelt),t1_meam(maxelt)
real*8 t2_meam(maxelt),t3_meam(maxelt)
real*8 rho_ref_meam(maxelt)
integer ibar_meam(maxelt),ielt_meam(maxelt)
character*3 lattce_meam(maxelt,maxelt)
integer nn2_meam(maxelt,maxelt)
integer zbl_meam(maxelt,maxelt)
integer eltind(maxelt,maxelt)
integer neltypes
real*8, allocatable :: phir(:,:)
real*8, allocatable :: phirar(:,:),phirar1(:,:),phirar2(:,:),
$ phirar3(:,:),phirar4(:,:),phirar5(:,:),phirar6(:,:)
real*8 attrac_meam(maxelt,maxelt),repuls_meam(maxelt,maxelt)
real*8 Cmin_meam(maxelt,maxelt,maxelt)
real*8 Cmax_meam(maxelt,maxelt,maxelt)
real*8 rc_meam,delr_meam,ebound_meam(maxelt,maxelt)
integer augt1, ialloy, mix_ref_t, erose_form
integer emb_lin_neg, bkgd_dyn
real*8 gsmooth_factor
integer vind2D(3,3),vind3D(3,3,3)
integer v2D(6),v3D(10)
integer nr,nrar
real*8 dr,rdrar
end module

View File

@ -1,296 +0,0 @@
c Extern "C" declaration has the form:
c
c void meam_dens_final_(int *, int *, int *, int *, int *, double *, double *,
c int *, int *, int *,
c double *, double *, double *, double *, double *, double *,
c double *, double *, double *, double *, double *, double *,
c double *, double *, double *, double *, double *, int *);
c
c Call from pair_meam.cpp has the form:
c
c meam_dens_final_(&nlocal,&nmax,&eflag_either,&eflag_global,&eflag_atom,
c &eng_vdwl,eatom,ntype,type,fmap,
c &arho1[0][0],&arho2[0][0],arho2b,&arho3[0][0],
c &arho3b[0][0],&t_ave[0][0],&tsq_ave[0][0],gamma,dgamma1,
c dgamma2,dgamma3,rho,rho0,rho1,rho2,rho3,frhop,&errorflag);
c
subroutine meam_dens_final(nlocal, nmax,
$ eflag_either, eflag_global, eflag_atom, eng_vdwl, eatom,
$ ntype, type, fmap,
$ Arho1, Arho2, Arho2b, Arho3, Arho3b, t_ave, tsq_ave,
$ Gamma, dGamma1, dGamma2, dGamma3,
$ rho, rho0, rho1, rho2, rho3, fp, errorflag)
use meam_data
implicit none
integer nlocal, nmax, eflag_either, eflag_global, eflag_atom
integer ntype, type, fmap
real*8 eng_vdwl, eatom, Arho1, Arho2
real*8 Arho2b, Arho3, Arho3b
real*8 t_ave, tsq_ave
real*8 Gamma, dGamma1, dGamma2, dGamma3
real*8 rho, rho0, rho1, rho2, rho3
real*8 fp
integer errorflag
dimension eatom(nmax)
dimension type(nmax), fmap(ntype)
dimension Arho1(3,nmax), Arho2(6,nmax), Arho2b(nmax)
dimension Arho3(10,nmax), Arho3b(3,nmax), t_ave(3,nmax)
dimension tsq_ave(3,nmax)
dimension Gamma(nmax), dGamma1(nmax), dGamma2(nmax)
dimension dGamma3(nmax), rho(nmax), rho0(nmax)
dimension rho1(nmax), rho2(nmax), rho3(nmax)
dimension fp(nmax)
integer i, elti
integer m
real*8 rhob, G, dG, Gbar, dGbar, gam, shp(3), shpi(3), Z
real*8 B, denom, rho_bkgd
c Complete the calculation of density
do i = 1,nlocal
elti = fmap(type(i))
if (elti.gt.0) then
rho1(i) = 0.d0
rho2(i) = -1.d0/3.d0*Arho2b(i)*Arho2b(i)
rho3(i) = 0.d0
do m = 1,3
rho1(i) = rho1(i) + Arho1(m,i)*Arho1(m,i)
rho3(i) = rho3(i) - 3.d0/5.d0*Arho3b(m,i)*Arho3b(m,i)
enddo
do m = 1,6
rho2(i) = rho2(i) + v2D(m)*Arho2(m,i)*Arho2(m,i)
enddo
do m = 1,10
rho3(i) = rho3(i) + v3D(m)*Arho3(m,i)*Arho3(m,i)
enddo
if( rho0(i) .gt. 0.0 ) then
if (ialloy.eq.1) then
if (tsq_ave(1,i) .ne. 0.0d0) then
t_ave(1,i) = t_ave(1,i)/tsq_ave(1,i)
else
t_ave(1,i) = 0.0d0
endif
if (tsq_ave(2,i) .ne. 0.0d0) then
t_ave(2,i) = t_ave(2,i)/tsq_ave(2,i)
else
t_ave(2,i) = 0.0d0
endif
if (tsq_ave(3,i) .ne. 0.0d0) then
t_ave(3,i) = t_ave(3,i)/tsq_ave(3,i)
else
t_ave(3,i) = 0.0d0
endif
else if (ialloy.eq.2) then
t_ave(1,i) = t1_meam(elti)
t_ave(2,i) = t2_meam(elti)
t_ave(3,i) = t3_meam(elti)
else
t_ave(1,i) = t_ave(1,i)/rho0(i)
t_ave(2,i) = t_ave(2,i)/rho0(i)
t_ave(3,i) = t_ave(3,i)/rho0(i)
endif
endif
Gamma(i) = t_ave(1,i)*rho1(i)
$ + t_ave(2,i)*rho2(i) + t_ave(3,i)*rho3(i)
if( rho0(i) .gt. 0.0 ) then
Gamma(i) = Gamma(i)/(rho0(i)*rho0(i))
end if
Z = Z_meam(elti)
call G_gam(Gamma(i),ibar_meam(elti),
$ gsmooth_factor,G,errorflag)
if (errorflag.ne.0) return
call get_shpfcn(shp,lattce_meam(elti,elti))
if (ibar_meam(elti).le.0) then
Gbar = 1.d0
dGbar = 0.d0
else
if (mix_ref_t.eq.1) then
gam = (t_ave(1,i)*shp(1)+t_ave(2,i)*shp(2)
$ +t_ave(3,i)*shp(3))/(Z*Z)
else
gam = (t1_meam(elti)*shp(1)+t2_meam(elti)*shp(2)
$ +t3_meam(elti)*shp(3))/(Z*Z)
endif
call G_gam(gam,ibar_meam(elti),gsmooth_factor,
$ Gbar,errorflag)
endif
rho(i) = rho0(i) * G
if (mix_ref_t.eq.1) then
if (ibar_meam(elti).le.0) then
Gbar = 1.d0
dGbar = 0.d0
else
gam = (t_ave(1,i)*shp(1)+t_ave(2,i)*shp(2)
$ +t_ave(3,i)*shp(3))/(Z*Z)
call dG_gam(gam,ibar_meam(elti),gsmooth_factor,
$ Gbar,dGbar)
endif
rho_bkgd = rho0_meam(elti)*Z*Gbar
else
if (bkgd_dyn.eq.1) then
rho_bkgd = rho0_meam(elti)*Z
else
rho_bkgd = rho_ref_meam(elti)
endif
endif
rhob = rho(i)/rho_bkgd
denom = 1.d0/rho_bkgd
call dG_gam(Gamma(i),ibar_meam(elti),gsmooth_factor,G,dG)
dGamma1(i) = (G - 2*dG*Gamma(i))*denom
if( rho0(i) .ne. 0.d0 ) then
dGamma2(i) = (dG/rho0(i))*denom
else
dGamma2(i) = 0.d0
end if
c dGamma3 is nonzero only if we are using the "mixed" rule for
c computing t in the reference system (which is not correct, but
c included for backward compatibility
if (mix_ref_t.eq.1) then
dGamma3(i) = rho0(i)*G*dGbar/(Gbar*Z*Z)*denom
else
dGamma3(i) = 0.0
endif
B = A_meam(elti)*Ec_meam(elti,elti)
if( rhob .ne. 0.d0 ) then
if (emb_lin_neg.eq.1 .and. rhob.le.0) then
fp(i) = -B
else
fp(i) = B*(log(rhob)+1.d0)
endif
if (eflag_either.ne.0) then
if (eflag_global.ne.0) then
if (emb_lin_neg.eq.1 .and. rhob.le.0) then
eng_vdwl = eng_vdwl - B*rhob
else
eng_vdwl = eng_vdwl + B*rhob*log(rhob)
endif
endif
if (eflag_atom.ne.0) then
if (emb_lin_neg.eq.1 .and. rhob.le.0) then
eatom(i) = eatom(i) - B*rhob
else
eatom(i) = eatom(i) + B*rhob*log(rhob)
endif
endif
endif
else
if (emb_lin_neg.eq.1) then
fp(i) = -B
else
fp(i) = B
endif
endif
endif
enddo
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine G_gam(Gamma,ibar,gsmooth_factor,G,errorflag)
c Compute G(Gamma) based on selection flag ibar:
c 0 => G = sqrt(1+Gamma)
c 1 => G = exp(Gamma/2)
c 2 => not implemented
c 3 => G = 2/(1+exp(-Gamma))
c 4 => G = sqrt(1+Gamma)
c -5 => G = +-sqrt(abs(1+Gamma))
use meam_data , only: fm_exp
implicit none
real*8 Gamma,G
real*8 gsmooth_factor, gsmooth_switchpoint
integer ibar, errorflag
if (ibar.eq.0.or.ibar.eq.4) then
gsmooth_switchpoint = -gsmooth_factor / (gsmooth_factor+1)
if (Gamma.lt.gsmooth_switchpoint) then
c e.g. gsmooth_factor is 99, then:
c gsmooth_switchpoint = -0.99
c G = 0.01*(-0.99/Gamma)**99
G = 1/(gsmooth_factor+1)
$ *(gsmooth_switchpoint/Gamma)**gsmooth_factor
G = sqrt(G)
else
G = sqrt(1.d0+Gamma)
endif
else if (ibar.eq.1) then
G = fm_exp(Gamma/2.d0)
else if (ibar.eq.3) then
G = 2.d0/(1.d0+exp(-Gamma))
else if (ibar.eq.-5) then
if ((1.d0+Gamma).ge.0) then
G = sqrt(1.d0+Gamma)
else
G = -sqrt(-1.d0-Gamma)
endif
else
errorflag = 1
endif
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine dG_gam(Gamma,ibar,gsmooth_factor,G,dG)
c Compute G(Gamma) and dG(gamma) based on selection flag ibar:
c 0 => G = sqrt(1+Gamma)
c 1 => G = fm_exp(Gamma/2)
c 2 => not implemented
c 3 => G = 2/(1+fm_exp(-Gamma))
c 4 => G = sqrt(1+Gamma)
c -5 => G = +-sqrt(abs(1+Gamma))
use meam_data , only: fm_exp
real*8 Gamma,G,dG
real*8 gsmooth_factor, gsmooth_switchpoint
integer ibar
if (ibar.eq.0.or.ibar.eq.4) then
gsmooth_switchpoint = -gsmooth_factor / (gsmooth_factor+1)
if (Gamma.lt.gsmooth_switchpoint) then
c e.g. gsmooth_factor is 99, then:
c gsmooth_switchpoint = -0.99
c G = 0.01*(-0.99/Gamma)**99
G = 1/(gsmooth_factor+1)
$ *(gsmooth_switchpoint/Gamma)**gsmooth_factor
G = sqrt(G)
dG = -gsmooth_factor*G/(2.0*Gamma)
else
G = sqrt(1.d0+Gamma)
dG = 1.d0/(2.d0*G)
endif
else if (ibar.eq.1) then
G = fm_exp(Gamma/2.d0)
dG = G/2.d0
else if (ibar.eq.3) then
G = 2.d0/(1.d0+fm_exp(-Gamma))
dG = G*(2.d0-G)/2
else if (ibar.eq.-5) then
if ((1.d0+Gamma).ge.0) then
G = sqrt(1.d0+Gamma)
dG = 1.d0/(2.d0*G)
else
G = -sqrt(-1.d0-Gamma)
dG = -1.d0/(2.d0*G)
endif
endif
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc

View File

@ -1,564 +0,0 @@
c Extern "C" declaration has the form:
c
c void meam_dens_init_(int *, int *, int *, double *, int *, int *, int *, double *,
c int *, int *, int *, int *,
c double *, double *, double *, double *, double *, double *,
c double *, double *, double *, double *, double *, int *);
c
c
c Call from pair_meam.cpp has the form:
c
c meam_dens_init_(&i,&nmax,ntype,type,fmap,&x[0][0],
c &numneigh[i],firstneigh[i],&numneigh_full[i],firstneigh_full[i],
c &scrfcn[offset],&dscrfcn[offset],&fcpair[offset],
c rho0,&arho1[0][0],&arho2[0][0],arho2b,
c &arho3[0][0],&arho3b[0][0],&t_ave[0][0],&tsq_ave[0][0],&errorflag);
c
subroutine meam_dens_init(i, nmax,
$ ntype, type, fmap, x,
$ numneigh, firstneigh,
$ numneigh_full, firstneigh_full,
$ scrfcn, dscrfcn, fcpair, rho0, arho1, arho2, arho2b,
$ arho3, arho3b, t_ave, tsq_ave, errorflag)
use meam_data
implicit none
integer i, nmax, ntype, type, fmap
real*8 x
integer numneigh, firstneigh, numneigh_full, firstneigh_full
real*8 scrfcn, dscrfcn, fcpair
real*8 rho0, arho1, arho2
real*8 arho2b, arho3, arho3b, t_ave, tsq_ave
integer errorflag
integer j,jn
dimension x(3,nmax)
dimension type(nmax), fmap(ntype)
dimension firstneigh(numneigh), firstneigh_full(numneigh_full)
dimension scrfcn(numneigh), dscrfcn(numneigh), fcpair(numneigh)
dimension rho0(nmax), arho1(3,nmax), arho2(6,nmax)
dimension arho2b(nmax), arho3(10,nmax), arho3b(3,nmax)
dimension t_ave(3,nmax), tsq_ave(3,nmax)
errorflag = 0
c Compute screening function and derivatives
call getscreen(i, nmax, scrfcn, dscrfcn, fcpair, x,
$ numneigh, firstneigh,
$ numneigh_full, firstneigh_full,
$ ntype, type, fmap)
c Calculate intermediate density terms to be communicated
call calc_rho1(i, nmax, ntype, type, fmap, x,
$ numneigh, firstneigh,
$ scrfcn, fcpair, rho0, arho1, arho2, arho2b,
$ arho3, arho3b, t_ave, tsq_ave)
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine getscreen(i, nmax, scrfcn, dscrfcn, fcpair, x,
$ numneigh, firstneigh,
$ numneigh_full, firstneigh_full,
$ ntype, type, fmap)
use meam_data
implicit none
integer i, nmax
real*8 scrfcn, dscrfcn, fcpair, x
integer numneigh, firstneigh, numneigh_full, firstneigh_full
integer ntype, type, fmap
dimension scrfcn(numneigh), dscrfcn(numneigh)
dimension fcpair(numneigh), x(3,nmax)
dimension firstneigh(numneigh), firstneigh_full(numneigh_full)
dimension type(nmax), fmap(ntype)
integer jn,j,kn,k
integer elti,eltj,eltk
real*8 xitmp,yitmp,zitmp,delxij,delyij,delzij,rij2,rij
real*8 xjtmp,yjtmp,zjtmp,delxik,delyik,delzik,rik2,rik
real*8 xktmp,yktmp,zktmp,delxjk,delyjk,delzjk,rjk2,rjk
real*8 xik,xjk,sij,fcij,sfcij,dfcij,sikj,dfikj,cikj
real*8 Cmin,Cmax,delc,ebound,rbound,a,coef1,coef2
real*8 coef1a,coef1b,coef2a,coef2b
real*8 dcikj
real*8 dC1a,dC1b,dC2a,dC2b
real*8 rnorm,fc,dfc,drinv
drinv = 1.d0/delr_meam
elti = fmap(type(i))
if (elti.gt.0) then
xitmp = x(1,i)
yitmp = x(2,i)
zitmp = x(3,i)
do jn = 1,numneigh
j = firstneigh(jn)
eltj = fmap(type(j))
if (eltj.gt.0) then
c First compute screening function itself, sij
xjtmp = x(1,j)
yjtmp = x(2,j)
zjtmp = x(3,j)
delxij = xjtmp - xitmp
delyij = yjtmp - yitmp
delzij = zjtmp - zitmp
rij2 = delxij*delxij + delyij*delyij + delzij*delzij
rij = sqrt(rij2)
if (rij.gt.rc_meam) then
fcij = 0.0
dfcij = 0.d0
sij = 0.d0
else
rnorm = (rc_meam-rij)*drinv
call screen(i, j, nmax, x, rij2, sij,
$ numneigh_full, firstneigh_full, ntype, type, fmap)
call dfcut(rnorm,fc,dfc)
fcij = fc
dfcij = dfc*drinv
endif
c Now compute derivatives
dscrfcn(jn) = 0.d0
sfcij = sij*fcij
if (sfcij.eq.0.d0.or.sfcij.eq.1.d0) goto 100
rbound = ebound_meam(elti,eltj) * rij2
do kn = 1,numneigh_full
k = firstneigh_full(kn)
if (k.eq.j) goto 10
eltk = fmap(type(k))
if (eltk.eq.0) goto 10
xktmp = x(1,k)
yktmp = x(2,k)
zktmp = x(3,k)
delxjk = xktmp - xjtmp
delyjk = yktmp - yjtmp
delzjk = zktmp - zjtmp
rjk2 = delxjk*delxjk + delyjk*delyjk + delzjk*delzjk
if (rjk2.gt.rbound) goto 10
delxik = xktmp - xitmp
delyik = yktmp - yitmp
delzik = zktmp - zitmp
rik2 = delxik*delxik + delyik*delyik + delzik*delzik
if (rik2.gt.rbound) goto 10
xik = rik2/rij2
xjk = rjk2/rij2
a = 1 - (xik-xjk)*(xik-xjk)
c if a < 0, then ellipse equation doesn't describe this case and
c atom k can't possibly screen i-j
if (a.le.0.d0) goto 10
cikj = (2.d0*(xik+xjk) + a - 2.d0)/a
Cmax = Cmax_meam(elti,eltj,eltk)
Cmin = Cmin_meam(elti,eltj,eltk)
if (cikj.ge.Cmax) then
goto 10
c Note that cikj may be slightly negative (within numerical
c tolerance) if atoms are colinear, so don't reject that case here
c (other negative cikj cases were handled by the test on "a" above)
c Note that we never have 0<cikj<Cmin here, else sij=0 (rejected above)
else
delc = Cmax - Cmin
cikj = (cikj-Cmin)/delc
call dfcut(cikj,sikj,dfikj)
coef1 = dfikj/(delc*sikj)
call dCfunc(rij2,rik2,rjk2,dCikj)
dscrfcn(jn) = dscrfcn(jn) + coef1*dCikj
endif
10 continue
enddo
coef1 = sfcij
coef2 = sij*dfcij/rij
dscrfcn(jn) = dscrfcn(jn)*coef1 - coef2
100 continue
scrfcn(jn) = sij
fcpair(jn) = fcij
endif
enddo
endif
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine calc_rho1(i, nmax, ntype, type, fmap, x,
$ numneigh, firstneigh,
$ scrfcn, fcpair, rho0, arho1, arho2, arho2b,
$ arho3, arho3b, t_ave, tsq_ave)
use meam_data
implicit none
integer i, nmax, ntype, type, fmap
real*8 x
integer numneigh, firstneigh
real*8 scrfcn, fcpair, rho0, arho1, arho2
real*8 arho2b, arho3, arho3b, t_ave, tsq_ave
dimension type(nmax), fmap(ntype), x(3,nmax)
dimension firstneigh(numneigh)
dimension scrfcn(numneigh), fcpair(numneigh)
dimension rho0(nmax), arho1(3,nmax), arho2(6,nmax)
dimension arho2b(nmax), arho3(10,nmax), arho3b(3,nmax)
dimension t_ave(3,nmax), tsq_ave(3,nmax)
integer jn,j,m,n,p,elti,eltj
integer nv2,nv3
real*8 xtmp,ytmp,ztmp,delij(3),rij2,rij,sij
real*8 ai,aj,rhoa0j,rhoa1j,rhoa2j,rhoa3j,A1j,A2j,A3j
real*8 G,Gbar,gam,shp(3)
real*8 ro0i,ro0j
real*8 rhoa0i,rhoa1i,rhoa2i,rhoa3i,A1i,A2i,A3i
elti = fmap(type(i))
xtmp = x(1,i)
ytmp = x(2,i)
ztmp = x(3,i)
do jn = 1,numneigh
if (scrfcn(jn).ne.0.d0) then
j = firstneigh(jn)
sij = scrfcn(jn)*fcpair(jn)
delij(1) = x(1,j) - xtmp
delij(2) = x(2,j) - ytmp
delij(3) = x(3,j) - ztmp
rij2 = delij(1)*delij(1) + delij(2)*delij(2)
$ + delij(3)*delij(3)
if (rij2.lt.cutforcesq) then
eltj = fmap(type(j))
rij = sqrt(rij2)
ai = rij/re_meam(elti,elti) - 1.d0
aj = rij/re_meam(eltj,eltj) - 1.d0
ro0i = rho0_meam(elti)
ro0j = rho0_meam(eltj)
rhoa0j = ro0j*fm_exp(-beta0_meam(eltj)*aj)*sij
rhoa1j = ro0j*fm_exp(-beta1_meam(eltj)*aj)*sij
rhoa2j = ro0j*fm_exp(-beta2_meam(eltj)*aj)*sij
rhoa3j = ro0j*fm_exp(-beta3_meam(eltj)*aj)*sij
rhoa0i = ro0i*fm_exp(-beta0_meam(elti)*ai)*sij
rhoa1i = ro0i*fm_exp(-beta1_meam(elti)*ai)*sij
rhoa2i = ro0i*fm_exp(-beta2_meam(elti)*ai)*sij
rhoa3i = ro0i*fm_exp(-beta3_meam(elti)*ai)*sij
if (ialloy.eq.1) then
rhoa1j = rhoa1j * t1_meam(eltj)
rhoa2j = rhoa2j * t2_meam(eltj)
rhoa3j = rhoa3j * t3_meam(eltj)
rhoa1i = rhoa1i * t1_meam(elti)
rhoa2i = rhoa2i * t2_meam(elti)
rhoa3i = rhoa3i * t3_meam(elti)
endif
rho0(i) = rho0(i) + rhoa0j
rho0(j) = rho0(j) + rhoa0i
c For ialloy = 2, use single-element value (not average)
if (ialloy.ne.2) then
t_ave(1,i) = t_ave(1,i) + t1_meam(eltj)*rhoa0j
t_ave(2,i) = t_ave(2,i) + t2_meam(eltj)*rhoa0j
t_ave(3,i) = t_ave(3,i) + t3_meam(eltj)*rhoa0j
t_ave(1,j) = t_ave(1,j) + t1_meam(elti)*rhoa0i
t_ave(2,j) = t_ave(2,j) + t2_meam(elti)*rhoa0i
t_ave(3,j) = t_ave(3,j) + t3_meam(elti)*rhoa0i
endif
if (ialloy.eq.1) then
tsq_ave(1,i) = tsq_ave(1,i) +
$ t1_meam(eltj)*t1_meam(eltj)*rhoa0j
tsq_ave(2,i) = tsq_ave(2,i) +
$ t2_meam(eltj)*t2_meam(eltj)*rhoa0j
tsq_ave(3,i) = tsq_ave(3,i) +
$ t3_meam(eltj)*t3_meam(eltj)*rhoa0j
tsq_ave(1,j) = tsq_ave(1,j) +
$ t1_meam(elti)*t1_meam(elti)*rhoa0i
tsq_ave(2,j) = tsq_ave(2,j) +
$ t2_meam(elti)*t2_meam(elti)*rhoa0i
tsq_ave(3,j) = tsq_ave(3,j) +
$ t3_meam(elti)*t3_meam(elti)*rhoa0i
endif
Arho2b(i) = Arho2b(i) + rhoa2j
Arho2b(j) = Arho2b(j) + rhoa2i
A1j = rhoa1j/rij
A2j = rhoa2j/rij2
A3j = rhoa3j/(rij2*rij)
A1i = rhoa1i/rij
A2i = rhoa2i/rij2
A3i = rhoa3i/(rij2*rij)
nv2 = 1
nv3 = 1
do m = 1,3
Arho1(m,i) = Arho1(m,i) + A1j*delij(m)
Arho1(m,j) = Arho1(m,j) - A1i*delij(m)
Arho3b(m,i) = Arho3b(m,i) + rhoa3j*delij(m)/rij
Arho3b(m,j) = Arho3b(m,j) - rhoa3i*delij(m)/rij
do n = m,3
Arho2(nv2,i) = Arho2(nv2,i) + A2j*delij(m)*delij(n)
Arho2(nv2,j) = Arho2(nv2,j) + A2i*delij(m)*delij(n)
nv2 = nv2+1
do p = n,3
Arho3(nv3,i) = Arho3(nv3,i)
$ + A3j*delij(m)*delij(n)*delij(p)
Arho3(nv3,j) = Arho3(nv3,j)
$ - A3i*delij(m)*delij(n)*delij(p)
nv3 = nv3+1
enddo
enddo
enddo
endif
endif
enddo
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine screen(i, j, nmax, x, rijsq, sij,
$ numneigh_full, firstneigh_full, ntype, type, fmap)
c Screening function
c Inputs: i = atom 1 id (integer)
c j = atom 2 id (integer)
c rijsq = squared distance between i and j
c Outputs: sij = screening function
use meam_data
implicit none
integer i,j,nmax,k,nk,m
real*8 x,rijsq,sij
integer numneigh_full, firstneigh_full
integer ntype, type, fmap
dimension x(3,nmax), firstneigh_full(numneigh_full)
dimension type(nmax), fmap(ntype)
integer elti,eltj,eltk
real*8 delxik,delyik,delzik
real*8 delxjk,delyjk,delzjk
real*8 riksq,rjksq,xik,xjk,cikj,a,delc,sikj,fcij,rij
real*8 Cmax,Cmin,rbound
sij = 1.d0
elti = fmap(type(i))
eltj = fmap(type(j))
c if rjksq > ebound*rijsq, atom k is definitely outside the ellipse
rbound = ebound_meam(elti,eltj)*rijsq
do nk = 1,numneigh_full
k = firstneigh_full(nk)
eltk = fmap(type(k))
if (k.eq.j) goto 10
delxjk = x(1,k) - x(1,j)
delyjk = x(2,k) - x(2,j)
delzjk = x(3,k) - x(3,j)
rjksq = delxjk*delxjk + delyjk*delyjk + delzjk*delzjk
if (rjksq.gt.rbound) goto 10
delxik = x(1,k) - x(1,i)
delyik = x(2,k) - x(2,i)
delzik = x(3,k) - x(3,i)
riksq = delxik*delxik + delyik*delyik + delzik*delzik
if (riksq.gt.rbound) goto 10
xik = riksq/rijsq
xjk = rjksq/rijsq
a = 1 - (xik-xjk)*(xik-xjk)
c if a < 0, then ellipse equation doesn't describe this case and
c atom k can't possibly screen i-j
if (a.le.0.d0) goto 10
cikj = (2.d0*(xik+xjk) + a - 2.d0)/a
Cmax = Cmax_meam(elti,eltj,eltk)
Cmin = Cmin_meam(elti,eltj,eltk)
if (cikj.ge.Cmax) then
goto 10
c note that cikj may be slightly negative (within numerical
c tolerance) if atoms are colinear, so don't reject that case here
c (other negative cikj cases were handled by the test on "a" above)
else if (cikj.le.Cmin) then
sij = 0.d0
goto 20
else
delc = Cmax - Cmin
cikj = (cikj-Cmin)/delc
call fcut(cikj,sikj)
endif
sij = sij * sikj
10 continue
enddo
20 continue
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine dsij(i,j,k,jn,nmax,numneigh,rij2,dsij1,dsij2,
$ ntype,type,fmap,x,scrfcn,fcpair)
c Inputs: i,j,k = id's of 3 atom triplet
c jn = id of i-j pair
c rij2 = squared distance between i and j
c Outputs: dsij1 = deriv. of sij w.r.t. rik
c dsij2 = deriv. of sij w.r.t. rjk
use meam_data
implicit none
integer i,j,k,jn,nmax,numneigh
integer elti,eltj,eltk
real*8 rij2,rik2,rjk2,dsij1,dsij2
integer ntype, type, fmap
real*8 x, scrfcn, fcpair
dimension type(nmax), fmap(ntype)
dimension x(3,nmax), scrfcn(numneigh), fcpair(numneigh)
real*8 dxik,dyik,dzik
real*8 dxjk,dyjk,dzjk
real*8 rbound,delc,sij,xik,xjk,cikj,sikj,dfc,a
real*8 Cmax,Cmin,dCikj1,dCikj2
sij = scrfcn(jn)*fcpair(jn)
elti = fmap(type(i))
eltj = fmap(type(j))
eltk = fmap(type(k))
Cmax = Cmax_meam(elti,eltj,eltk)
Cmin = Cmin_meam(elti,eltj,eltk)
dsij1 = 0.d0
dsij2 = 0.d0
if ((sij.ne.0.d0).and.(sij.ne.1.d0)) then
rbound = rij2*ebound_meam(elti,eltj)
delc = Cmax-Cmin
dxjk = x(1,k) - x(1,j)
dyjk = x(2,k) - x(2,j)
dzjk = x(3,k) - x(3,j)
rjk2 = dxjk*dxjk + dyjk*dyjk + dzjk*dzjk
if (rjk2.le.rbound) then
dxik = x(1,k) - x(1,i)
dyik = x(2,k) - x(2,i)
dzik = x(3,k) - x(3,i)
rik2 = dxik*dxik + dyik*dyik + dzik*dzik
if (rik2.le.rbound) then
xik = rik2/rij2
xjk = rjk2/rij2
a = 1 - (xik-xjk)*(xik-xjk)
if (a.ne.0.d0) then
cikj = (2.d0*(xik+xjk) + a - 2.d0)/a
if (cikj.ge.Cmin.and.cikj.le.Cmax) then
cikj = (cikj-Cmin)/delc
call dfcut(cikj,sikj,dfc)
call dCfunc2(rij2,rik2,rjk2,dCikj1,dCikj2)
a = sij/delc*dfc/sikj
dsij1 = a*dCikj1
dsij2 = a*dCikj2
endif
endif
endif
endif
endif
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine fcut(xi,fc)
c cutoff function
implicit none
real*8 xi,fc
real*8 a
if (xi.ge.1.d0) then
fc = 1.d0
else if (xi.le.0.d0) then
fc = 0.d0
else
a = 1.d0-xi
a = a*a
a = a*a
a = 1.d0-a
fc = a*a
c fc = xi
endif
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine dfcut(xi,fc,dfc)
c cutoff function and its derivative
implicit none
real*8 xi,fc,dfc,a,a3,a4
if (xi.ge.1.d0) then
fc = 1.d0
dfc = 0.d0
else if (xi.le.0.d0) then
fc = 0.d0
dfc = 0.d0
else
a = 1.d0-xi
a3 = a*a*a
a4 = a*a3
fc = (1.d0-a4)**2
dfc = 8*(1.d0-a4)*a3
c fc = xi
c dfc = 1.d0
endif
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine dCfunc(rij2,rik2,rjk2,dCikj)
c Inputs: rij,rij2,rik2,rjk2
c Outputs: dCikj = derivative of Cikj w.r.t. rij
implicit none
real*8 rij2,rik2,rjk2,dCikj
real*8 rij4,a,b,denom
rij4 = rij2*rij2
a = rik2-rjk2
b = rik2+rjk2
denom = rij4 - a*a
denom = denom*denom
dCikj = -4*(-2*rij2*a*a + rij4*b + a*a*b)/denom
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
subroutine dCfunc2(rij2,rik2,rjk2,dCikj1,dCikj2)
c Inputs: rij,rij2,rik2,rjk2
c Outputs: dCikj1 = derivative of Cikj w.r.t. rik
c dCikj2 = derivative of Cikj w.r.t. rjk
implicit none
real*8 rij2,rik2,rjk2,dCikj1,dCikj2
real*8 rij4,rik4,rjk4,a,b,denom
rij4 = rij2*rij2
rik4 = rik2*rik2
rjk4 = rjk2*rjk2
a = rik2-rjk2
b = rik2+rjk2
denom = rij4 - a*a
denom = denom*denom
dCikj1 = 4*rij2*(rij4 + rik4 + 2*rik2*rjk2 - 3*rjk4 - 2*rij2*a)/
$ denom
dCikj2 = 4*rij2*(rij4 - 3*rik4 + 2*rik2*rjk2 + rjk4 + 2*rij2*a)/
$ denom
return
end
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc

View File

@ -1,608 +0,0 @@
c Extern "C" declaration has the form:
c
c void meam_force_(int *, int *, int *, double *, int *, int *, int *, double *,
c int *, int *, int *, int *, double *, double *,
c double *, double *, double *, double *, double *, double *,
c double *, double *, double *, double *, double *, double *,
c double *, double *, double *, double *, double *, double *, int *);
c
c Call from pair_meam.cpp has the form:
c
c meam_force_(&i,&nmax,&eflag_either,&eflag_global,&eflag_atom,&vflag_atom,
c &eng_vdwl,eatom,&ntype,type,fmap,&x[0][0],
c &numneigh[i],firstneigh[i],&numneigh_full[i],firstneigh_full[i],
c &scrfcn[offset],&dscrfcn[offset],&fcpair[offset],
c dgamma1,dgamma2,dgamma3,rho0,rho1,rho2,rho3,frhop,
c &arho1[0][0],&arho2[0][0],arho2b,&arho3[0][0],&arho3b[0][0],
c &t_ave[0][0],&tsq_ave[0][0],&f[0][0],&vatom[0][0],&errorflag);
c
subroutine meam_force(i, nmax,
$ eflag_either, eflag_global, eflag_atom, vflag_atom,
$ eng_vdwl, eatom, ntype, type, fmap, x,
$ numneigh, firstneigh, numneigh_full, firstneigh_full,
$ scrfcn, dscrfcn, fcpair,
$ dGamma1, dGamma2, dGamma3, rho0, rho1, rho2, rho3, fp,
$ Arho1, Arho2, Arho2b, Arho3, Arho3b, t_ave, tsq_ave, f,
$ vatom, errorflag)
use meam_data
implicit none
integer eflag_either, eflag_global, eflag_atom, vflag_atom
integer nmax, ntype, type, fmap
real*8 eng_vdwl, eatom, x
integer numneigh, firstneigh, numneigh_full, firstneigh_full
real*8 scrfcn, dscrfcn, fcpair
real*8 dGamma1, dGamma2, dGamma3
real*8 rho0, rho1, rho2, rho3, fp
real*8 Arho1, Arho2, Arho2b
real*8 Arho3, Arho3b
real*8 t_ave, tsq_ave, f, vatom
integer errorflag
dimension eatom(nmax)
dimension type(nmax), fmap(ntype)
dimension x(3,nmax)
dimension firstneigh(numneigh), firstneigh_full(numneigh_full)
dimension scrfcn(numneigh), dscrfcn(numneigh), fcpair(numneigh)
dimension dGamma1(nmax), dGamma2(nmax), dGamma3(nmax)
dimension rho0(nmax), rho1(nmax), rho2(nmax), rho3(nmax), fp(nmax)
dimension Arho1(3,nmax), Arho2(6,nmax), Arho2b(nmax)
dimension Arho3(10,nmax), Arho3b(3,nmax)
dimension t_ave(3,nmax), tsq_ave(3,nmax), f(3,nmax), vatom(6,nmax)
integer i,j,jn,k,kn,kk,m,n,p,q
integer nv2,nv3,elti,eltj,eltk,ind
real*8 xitmp,yitmp,zitmp,delij(3),delref(3),rij2,rij,rij3
real*8 delik(3),deljk(3),v(6),fi(3),fj(3)
real*8 Eu,astar,astarp,third,sixth
real*8 pp,phiforce,dUdrij,dUdsij,dUdrijm(3),force,forcem
real*8 B,r,recip,phi,phip,rhop,a
real*8 sij,fcij,dfcij,ds(3)
real*8 a0,a1,a1i,a1j,a2,a2i,a2j
real*8 a3i,a3j,a3i1,a3i2,a3j1,a3j2
real*8 G,dG,Gbar,dGbar,gam,shpi(3),shpj(3),Z,denom
real*8 ai,aj,ro0i,ro0j,invrei,invrej
real*8 b0,rhoa0j,drhoa0j,rhoa0i,drhoa0i
real*8 b1,rhoa1j,drhoa1j,rhoa1i,drhoa1i
real*8 b2,rhoa2j,drhoa2j,rhoa2i,drhoa2i
real*8 a3,a3a,b3,rhoa3j,drhoa3j,rhoa3i,drhoa3i
real*8 drho0dr1,drho0dr2,drho0ds1,drho0ds2
real*8 drho1dr1,drho1dr2,drho1ds1,drho1ds2
real*8 drho1drm1(3),drho1drm2(3)
real*8 drho2dr1,drho2dr2,drho2ds1,drho2ds2
real*8 drho2drm1(3),drho2drm2(3)
real*8 drho3dr1,drho3dr2,drho3ds1,drho3ds2
real*8 drho3drm1(3),drho3drm2(3)
real*8 dt1dr1,dt1dr2,dt1ds1,dt1ds2
real*8 dt2dr1,dt2dr2,dt2ds1,dt2ds2
real*8 dt3dr1,dt3dr2,dt3ds1,dt3ds2
real*8 drhodr1,drhodr2,drhods1,drhods2,drhodrm1(3),drhodrm2(3)
real*8 arg,arg1,arg2
real*8 arg1i1,arg1j1,arg1i2,arg1j2,arg2i2,arg2j2
real*8 arg1i3,arg1j3,arg2i3,arg2j3,arg3i3,arg3j3
real*8 dsij1,dsij2,force1,force2
real*8 t1i,t2i,t3i,t1j,t2j,t3j
errorflag = 0
third = 1.0/3.0
sixth = 1.0/6.0
c Compute forces atom i
elti = fmap(type(i))
if (elti.gt.0) then
xitmp = x(1,i)
yitmp = x(2,i)
zitmp = x(3,i)
c Treat each pair
do jn = 1,numneigh
j = firstneigh(jn)
eltj = fmap(type(j))
if (scrfcn(jn).ne.0.d0.and.eltj.gt.0) then
sij = scrfcn(jn)*fcpair(jn)
delij(1) = x(1,j) - xitmp
delij(2) = x(2,j) - yitmp
delij(3) = x(3,j) - zitmp
rij2 = delij(1)*delij(1) + delij(2)*delij(2)
$ + delij(3)*delij(3)
if (rij2.lt.cutforcesq) then
rij = sqrt(rij2)
r = rij
c Compute phi and phip
ind = eltind(elti,eltj)
pp = rij*rdrar + 1.0D0
kk = pp
kk = min(kk,nrar-1)
pp = pp - kk
pp = min(pp,1.0D0)
phi = ((phirar3(kk,ind)*pp + phirar2(kk,ind))*pp
$ + phirar1(kk,ind))*pp + phirar(kk,ind)
phip = (phirar6(kk,ind)*pp + phirar5(kk,ind))*pp
$ + phirar4(kk,ind)
recip = 1.0d0/r
if (eflag_either.ne.0) then
if (eflag_global.ne.0) eng_vdwl = eng_vdwl + phi*sij
if (eflag_atom.ne.0) then
eatom(i) = eatom(i) + 0.5*phi*sij
eatom(j) = eatom(j) + 0.5*phi*sij
endif
endif
c write(1,*) "force_meamf: phi: ",phi
c write(1,*) "force_meamf: phip: ",phip
c Compute pair densities and derivatives
invrei = 1.d0/re_meam(elti,elti)
ai = rij*invrei - 1.d0
ro0i = rho0_meam(elti)
rhoa0i = ro0i*fm_exp(-beta0_meam(elti)*ai)
drhoa0i = -beta0_meam(elti)*invrei*rhoa0i
rhoa1i = ro0i*fm_exp(-beta1_meam(elti)*ai)
drhoa1i = -beta1_meam(elti)*invrei*rhoa1i
rhoa2i = ro0i*fm_exp(-beta2_meam(elti)*ai)
drhoa2i = -beta2_meam(elti)*invrei*rhoa2i
rhoa3i = ro0i*fm_exp(-beta3_meam(elti)*ai)
drhoa3i = -beta3_meam(elti)*invrei*rhoa3i
if (elti.ne.eltj) then
invrej = 1.d0/re_meam(eltj,eltj)
aj = rij*invrej - 1.d0
ro0j = rho0_meam(eltj)
rhoa0j = ro0j*fm_exp(-beta0_meam(eltj)*aj)
drhoa0j = -beta0_meam(eltj)*invrej*rhoa0j
rhoa1j = ro0j*fm_exp(-beta1_meam(eltj)*aj)
drhoa1j = -beta1_meam(eltj)*invrej*rhoa1j
rhoa2j = ro0j*fm_exp(-beta2_meam(eltj)*aj)
drhoa2j = -beta2_meam(eltj)*invrej*rhoa2j
rhoa3j = ro0j*fm_exp(-beta3_meam(eltj)*aj)
drhoa3j = -beta3_meam(eltj)*invrej*rhoa3j
else
rhoa0j = rhoa0i
drhoa0j = drhoa0i
rhoa1j = rhoa1i
drhoa1j = drhoa1i
rhoa2j = rhoa2i
drhoa2j = drhoa2i
rhoa3j = rhoa3i
drhoa3j = drhoa3i
endif
if (ialloy.eq.1) then
rhoa1j = rhoa1j * t1_meam(eltj)
rhoa2j = rhoa2j * t2_meam(eltj)
rhoa3j = rhoa3j * t3_meam(eltj)
rhoa1i = rhoa1i * t1_meam(elti)
rhoa2i = rhoa2i * t2_meam(elti)
rhoa3i = rhoa3i * t3_meam(elti)
drhoa1j = drhoa1j * t1_meam(eltj)
drhoa2j = drhoa2j * t2_meam(eltj)
drhoa3j = drhoa3j * t3_meam(eltj)
drhoa1i = drhoa1i * t1_meam(elti)
drhoa2i = drhoa2i * t2_meam(elti)
drhoa3i = drhoa3i * t3_meam(elti)
endif
nv2 = 1
nv3 = 1
arg1i1 = 0.d0
arg1j1 = 0.d0
arg1i2 = 0.d0
arg1j2 = 0.d0
arg1i3 = 0.d0
arg1j3 = 0.d0
arg3i3 = 0.d0
arg3j3 = 0.d0
do n = 1,3
do p = n,3
do q = p,3
arg = delij(n)*delij(p)*delij(q)*v3D(nv3)
arg1i3 = arg1i3 + Arho3(nv3,i)*arg
arg1j3 = arg1j3 - Arho3(nv3,j)*arg
nv3 = nv3+1
enddo
arg = delij(n)*delij(p)*v2D(nv2)
arg1i2 = arg1i2 + Arho2(nv2,i)*arg
arg1j2 = arg1j2 + Arho2(nv2,j)*arg
nv2 = nv2+1
enddo
arg1i1 = arg1i1 + Arho1(n,i)*delij(n)
arg1j1 = arg1j1 - Arho1(n,j)*delij(n)
arg3i3 = arg3i3 + Arho3b(n,i)*delij(n)
arg3j3 = arg3j3 - Arho3b(n,j)*delij(n)
enddo
c rho0 terms
drho0dr1 = drhoa0j * sij
drho0dr2 = drhoa0i * sij
c rho1 terms
a1 = 2*sij/rij
drho1dr1 = a1*(drhoa1j-rhoa1j/rij)*arg1i1
drho1dr2 = a1*(drhoa1i-rhoa1i/rij)*arg1j1
a1 = 2.d0*sij/rij
do m = 1,3
drho1drm1(m) = a1*rhoa1j*Arho1(m,i)
drho1drm2(m) = -a1*rhoa1i*Arho1(m,j)
enddo
c rho2 terms
a2 = 2*sij/rij2
drho2dr1 = a2*(drhoa2j - 2*rhoa2j/rij)*arg1i2
$ - 2.d0/3.d0*Arho2b(i)*drhoa2j*sij
drho2dr2 = a2*(drhoa2i - 2*rhoa2i/rij)*arg1j2
$ - 2.d0/3.d0*Arho2b(j)*drhoa2i*sij
a2 = 4*sij/rij2
do m = 1,3
drho2drm1(m) = 0.d0
drho2drm2(m) = 0.d0
do n = 1,3
drho2drm1(m) = drho2drm1(m)
$ + Arho2(vind2D(m,n),i)*delij(n)
drho2drm2(m) = drho2drm2(m)
$ - Arho2(vind2D(m,n),j)*delij(n)
enddo
drho2drm1(m) = a2*rhoa2j*drho2drm1(m)
drho2drm2(m) = -a2*rhoa2i*drho2drm2(m)
enddo
c rho3 terms
rij3 = rij*rij2
a3 = 2*sij/rij3
a3a = 6.d0/5.d0*sij/rij
drho3dr1 = a3*(drhoa3j - 3*rhoa3j/rij)*arg1i3
$ - a3a*(drhoa3j - rhoa3j/rij)*arg3i3
drho3dr2 = a3*(drhoa3i - 3*rhoa3i/rij)*arg1j3
$ - a3a*(drhoa3i - rhoa3i/rij)*arg3j3
a3 = 6*sij/rij3
a3a = 6*sij/(5*rij)
do m = 1,3
drho3drm1(m) = 0.d0
drho3drm2(m) = 0.d0
nv2 = 1
do n = 1,3
do p = n,3
arg = delij(n)*delij(p)*v2D(nv2)
drho3drm1(m) = drho3drm1(m)
$ + Arho3(vind3D(m,n,p),i)*arg
drho3drm2(m) = drho3drm2(m)
$ + Arho3(vind3D(m,n,p),j)*arg
nv2 = nv2 + 1
enddo
enddo
drho3drm1(m) = (a3*drho3drm1(m) - a3a*Arho3b(m,i))
$ *rhoa3j
drho3drm2(m) = (-a3*drho3drm2(m) + a3a*Arho3b(m,j))
$ *rhoa3i
enddo
c Compute derivatives of weighting functions t wrt rij
t1i = t_ave(1,i)
t2i = t_ave(2,i)
t3i = t_ave(3,i)
t1j = t_ave(1,j)
t2j = t_ave(2,j)
t3j = t_ave(3,j)
if (ialloy.eq.1) then
a1i = 0.d0
a1j = 0.d0
a2i = 0.d0
a2j = 0.d0
a3i = 0.d0
a3j = 0.d0
if ( tsq_ave(1,i) .ne. 0.d0 ) then
a1i = drhoa0j*sij/tsq_ave(1,i)
endif
if ( tsq_ave(1,j) .ne. 0.d0 ) then
a1j = drhoa0i*sij/tsq_ave(1,j)
endif
if ( tsq_ave(2,i) .ne. 0.d0 ) then
a2i = drhoa0j*sij/tsq_ave(2,i)
endif
if ( tsq_ave(2,j) .ne. 0.d0 ) then
a2j = drhoa0i*sij/tsq_ave(2,j)
endif
if ( tsq_ave(3,i) .ne. 0.d0 ) then
a3i = drhoa0j*sij/tsq_ave(3,i)
endif
if ( tsq_ave(3,j) .ne. 0.d0 ) then
a3j = drhoa0i*sij/tsq_ave(3,j)
endif
dt1dr1 = a1i*(t1_meam(eltj)-t1i*t1_meam(eltj)**2)
dt1dr2 = a1j*(t1_meam(elti)-t1j*t1_meam(elti)**2)
dt2dr1 = a2i*(t2_meam(eltj)-t2i*t2_meam(eltj)**2)
dt2dr2 = a2j*(t2_meam(elti)-t2j*t2_meam(elti)**2)
dt3dr1 = a3i*(t3_meam(eltj)-t3i*t3_meam(eltj)**2)
dt3dr2 = a3j*(t3_meam(elti)-t3j*t3_meam(elti)**2)
else if (ialloy.eq.2) then
dt1dr1 = 0.d0
dt1dr2 = 0.d0
dt2dr1 = 0.d0
dt2dr2 = 0.d0
dt3dr1 = 0.d0
dt3dr2 = 0.d0
else
ai = 0.d0
if( rho0(i) .ne. 0.d0 ) then
ai = drhoa0j*sij/rho0(i)
end if
aj = 0.d0
if( rho0(j) .ne. 0.d0 ) then
aj = drhoa0i*sij/rho0(j)
end if
dt1dr1 = ai*(t1_meam(eltj)-t1i)
dt1dr2 = aj*(t1_meam(elti)-t1j)
dt2dr1 = ai*(t2_meam(eltj)-t2i)
dt2dr2 = aj*(t2_meam(elti)-t2j)
dt3dr1 = ai*(t3_meam(eltj)-t3i)
dt3dr2 = aj*(t3_meam(elti)-t3j)
endif
c Compute derivatives of total density wrt rij, sij and rij(3)
call get_shpfcn(shpi,lattce_meam(elti,elti))
call get_shpfcn(shpj,lattce_meam(eltj,eltj))
drhodr1 = dGamma1(i)*drho0dr1
$ + dGamma2(i)*
$ (dt1dr1*rho1(i)+t1i*drho1dr1
$ + dt2dr1*rho2(i)+t2i*drho2dr1
$ + dt3dr1*rho3(i)+t3i*drho3dr1)
$ - dGamma3(i)*
$ (shpi(1)*dt1dr1+shpi(2)*dt2dr1+shpi(3)*dt3dr1)
drhodr2 = dGamma1(j)*drho0dr2
$ + dGamma2(j)*
$ (dt1dr2*rho1(j)+t1j*drho1dr2
$ + dt2dr2*rho2(j)+t2j*drho2dr2
$ + dt3dr2*rho3(j)+t3j*drho3dr2)
$ - dGamma3(j)*
$ (shpj(1)*dt1dr2+shpj(2)*dt2dr2+shpj(3)*dt3dr2)
do m = 1,3
drhodrm1(m) = 0.d0
drhodrm2(m) = 0.d0
drhodrm1(m) = dGamma2(i)*
$ (t1i*drho1drm1(m)
$ + t2i*drho2drm1(m)
$ + t3i*drho3drm1(m))
drhodrm2(m) = dGamma2(j)*
$ (t1j*drho1drm2(m)
$ + t2j*drho2drm2(m)
$ + t3j*drho3drm2(m))
enddo
c Compute derivatives wrt sij, but only if necessary
if (dscrfcn(jn).ne.0.d0) then
drho0ds1 = rhoa0j
drho0ds2 = rhoa0i
a1 = 2.d0/rij
drho1ds1 = a1*rhoa1j*arg1i1
drho1ds2 = a1*rhoa1i*arg1j1
a2 = 2.d0/rij2
drho2ds1 = a2*rhoa2j*arg1i2
$ - 2.d0/3.d0*Arho2b(i)*rhoa2j
drho2ds2 = a2*rhoa2i*arg1j2
$ - 2.d0/3.d0*Arho2b(j)*rhoa2i
a3 = 2.d0/rij3
a3a = 6.d0/(5.d0*rij)
drho3ds1 = a3*rhoa3j*arg1i3 - a3a*rhoa3j*arg3i3
drho3ds2 = a3*rhoa3i*arg1j3 - a3a*rhoa3i*arg3j3
if (ialloy.eq.1) then
a1i = 0.d0
a1j = 0.d0
a2i = 0.d0
a2j = 0.d0
a3i = 0.d0
a3j = 0.d0
if ( tsq_ave(1,i) .ne. 0.d0 ) then
a1i = rhoa0j/tsq_ave(1,i)
endif
if ( tsq_ave(1,j) .ne. 0.d0 ) then
a1j = rhoa0i/tsq_ave(1,j)
endif
if ( tsq_ave(2,i) .ne. 0.d0 ) then
a2i = rhoa0j/tsq_ave(2,i)
endif
if ( tsq_ave(2,j) .ne. 0.d0 ) then
a2j = rhoa0i/tsq_ave(2,j)
endif
if ( tsq_ave(3,i) .ne. 0.d0 ) then
a3i = rhoa0j/tsq_ave(3,i)
endif
if ( tsq_ave(3,j) .ne. 0.d0 ) then
a3j = rhoa0i/tsq_ave(3,j)
endif
dt1ds1 = a1i*(t1_meam(eltj)-t1i*t1_meam(eltj)**2)
dt1ds2 = a1j*(t1_meam(elti)-t1j*t1_meam(elti)**2)
dt2ds1 = a2i*(t2_meam(eltj)-t2i*t2_meam(eltj)**2)
dt2ds2 = a2j*(t2_meam(elti)-t2j*t2_meam(elti)**2)
dt3ds1 = a3i*(t3_meam(eltj)-t3i*t3_meam(eltj)**2)
dt3ds2 = a3j*(t3_meam(elti)-t3j*t3_meam(elti)**2)
else if (ialloy.eq.2) then
dt1ds1 = 0.d0
dt1ds2 = 0.d0
dt2ds1 = 0.d0
dt2ds2 = 0.d0
dt3ds1 = 0.d0
dt3ds2 = 0.d0
else
ai = 0.d0
if( rho0(i) .ne. 0.d0 ) then
ai = rhoa0j/rho0(i)
end if
aj = 0.d0
if( rho0(j) .ne. 0.d0 ) then
aj = rhoa0i/rho0(j)
end if
dt1ds1 = ai*(t1_meam(eltj)-t1i)
dt1ds2 = aj*(t1_meam(elti)-t1j)
dt2ds1 = ai*(t2_meam(eltj)-t2i)
dt2ds2 = aj*(t2_meam(elti)-t2j)
dt3ds1 = ai*(t3_meam(eltj)-t3i)
dt3ds2 = aj*(t3_meam(elti)-t3j)
endif
drhods1 = dGamma1(i)*drho0ds1
$ + dGamma2(i)*
$ (dt1ds1*rho1(i)+t1i*drho1ds1
$ + dt2ds1*rho2(i)+t2i*drho2ds1
$ + dt3ds1*rho3(i)+t3i*drho3ds1)
$ - dGamma3(i)*
$ (shpi(1)*dt1ds1+shpi(2)*dt2ds1+shpi(3)*dt3ds1)
drhods2 = dGamma1(j)*drho0ds2
$ + dGamma2(j)*
$ (dt1ds2*rho1(j)+t1j*drho1ds2
$ + dt2ds2*rho2(j)+t2j*drho2ds2
$ + dt3ds2*rho3(j)+t3j*drho3ds2)
$ - dGamma3(j)*
$ (shpj(1)*dt1ds2+shpj(2)*dt2ds2+shpj(3)*dt3ds2)
endif
c Compute derivatives of energy wrt rij, sij and rij(3)
dUdrij = phip*sij
$ + fp(i)*drhodr1 + fp(j)*drhodr2
dUdsij = 0.d0
if (dscrfcn(jn).ne.0.d0) then
dUdsij = phi
$ + fp(i)*drhods1 + fp(j)*drhods2
endif
do m = 1,3
dUdrijm(m) = fp(i)*drhodrm1(m) + fp(j)*drhodrm2(m)
enddo
c Add the part of the force due to dUdrij and dUdsij
force = dUdrij*recip + dUdsij*dscrfcn(jn)
do m = 1,3
forcem = delij(m)*force + dUdrijm(m)
f(m,i) = f(m,i) + forcem
f(m,j) = f(m,j) - forcem
enddo
c Tabulate per-atom virial as symmetrized stress tensor
if (vflag_atom.ne.0) then
fi(1) = delij(1)*force + dUdrijm(1)
fi(2) = delij(2)*force + dUdrijm(2)
fi(3) = delij(3)*force + dUdrijm(3)
v(1) = -0.5 * (delij(1) * fi(1))
v(2) = -0.5 * (delij(2) * fi(2))
v(3) = -0.5 * (delij(3) * fi(3))
v(4) = -0.25 * (delij(1)*fi(2) + delij(2)*fi(1))
v(5) = -0.25 * (delij(1)*fi(3) + delij(3)*fi(1))
v(6) = -0.25 * (delij(2)*fi(3) + delij(3)*fi(2))
vatom(1,i) = vatom(1,i) + v(1)
vatom(2,i) = vatom(2,i) + v(2)
vatom(3,i) = vatom(3,i) + v(3)
vatom(4,i) = vatom(4,i) + v(4)
vatom(5,i) = vatom(5,i) + v(5)
vatom(6,i) = vatom(6,i) + v(6)
vatom(1,j) = vatom(1,j) + v(1)
vatom(2,j) = vatom(2,j) + v(2)
vatom(3,j) = vatom(3,j) + v(3)
vatom(4,j) = vatom(4,j) + v(4)
vatom(5,j) = vatom(5,j) + v(5)
vatom(6,j) = vatom(6,j) + v(6)
endif
c Now compute forces on other atoms k due to change in sij
if (sij.eq.0.d0.or.sij.eq.1.d0) goto 100
do kn = 1,numneigh_full
k = firstneigh_full(kn)
eltk = fmap(type(k))
if (k.ne.j.and.eltk.gt.0) then
call dsij(i,j,k,jn,nmax,numneigh,rij2,dsij1,dsij2,
$ ntype,type,fmap,x,scrfcn,fcpair)
if (dsij1.ne.0.d0.or.dsij2.ne.0.d0) then
force1 = dUdsij*dsij1
force2 = dUdsij*dsij2
do m = 1,3
delik(m) = x(m,k) - x(m,i)
deljk(m) = x(m,k) - x(m,j)
enddo
do m = 1,3
f(m,i) = f(m,i) + force1*delik(m)
f(m,j) = f(m,j) + force2*deljk(m)
f(m,k) = f(m,k) - force1*delik(m)
$ - force2*deljk(m)
enddo
c Tabulate per-atom virial as symmetrized stress tensor
if (vflag_atom.ne.0) then
fi(1) = force1*delik(1)
fi(2) = force1*delik(2)
fi(3) = force1*delik(3)
fj(1) = force2*deljk(1)
fj(2) = force2*deljk(2)
fj(3) = force2*deljk(3)
v(1) = -third * (delik(1)*fi(1) + deljk(1)*fj(1))
v(2) = -third * (delik(2)*fi(2) + deljk(2)*fj(2))
v(3) = -third * (delik(3)*fi(3) + deljk(3)*fj(3))
v(4) = -sixth * (delik(1)*fi(2) + deljk(1)*fj(2) +
$ delik(2)*fi(1) + deljk(2)*fj(1))
v(5) = -sixth * (delik(1)*fi(3) + deljk(1)*fj(3) +
$ delik(3)*fi(1) + deljk(3)*fj(1))
v(6) = -sixth * (delik(2)*fi(3) + deljk(2)*fj(3) +
$ delik(3)*fi(2) + deljk(3)*fj(2))
vatom(1,i) = vatom(1,i) + v(1)
vatom(2,i) = vatom(2,i) + v(2)
vatom(3,i) = vatom(3,i) + v(3)
vatom(4,i) = vatom(4,i) + v(4)
vatom(5,i) = vatom(5,i) + v(5)
vatom(6,i) = vatom(6,i) + v(6)
vatom(1,j) = vatom(1,j) + v(1)
vatom(2,j) = vatom(2,j) + v(2)
vatom(3,j) = vatom(3,j) + v(3)
vatom(4,j) = vatom(4,j) + v(4)
vatom(5,j) = vatom(5,j) + v(5)
vatom(6,j) = vatom(6,j) + v(6)
vatom(1,k) = vatom(1,k) + v(1)
vatom(2,k) = vatom(2,k) + v(2)
vatom(3,k) = vatom(3,k) + v(3)
vatom(4,k) = vatom(4,k) + v(4)
vatom(5,k) = vatom(5,k) + v(5)
vatom(6,k) = vatom(6,k) + v(6)
endif
endif
endif
c end of k loop
enddo
endif
100 continue
endif
c end of j loop
enddo
c else if elti=0, this is not a meam atom
endif
return
end

File diff suppressed because it is too large Load Diff

View File

@ -1,111 +0,0 @@
c
c declaration in pair_meam.h:
c
c void meam_setup_global(int *, int *, double *, int *, double *, double *,
c double *, double *, double *, double *, double *,
c double *, double *, double *, double *, double *,
c double *, double *, int *);
c
c call in pair_meam.cpp:
c
c meam_setup_global(&nelements,lat,z,ielement,atwt,alpha,b0,b1,b2,b3,
c alat,esub,asub,t0,t1,t2,t3,rozero,ibar);
c
c
subroutine meam_setup_global(nelt, lat, z, ielement, atwt, alpha,
$ b0, b1, b2, b3, alat, esub, asub,
$ t0, t1, t2, t3, rozero, ibar)
use meam_data
implicit none
integer nelt, lat, ielement, ibar
real*8 z, atwt, alpha, b0, b1, b2, b3
real*8 alat, esub, asub, t0, t1, t2, t3
real*8 rozero
dimension lat(nelt), ielement(nelt), ibar(nelt)
dimension z(nelt), atwt(nelt), alpha(nelt)
dimension b0(nelt), b1(nelt), b2(nelt), b3(nelt)
dimension alat(nelt), esub(nelt), asub(nelt)
dimension t0(nelt), t1(nelt), t2(nelt), t3(nelt), rozero(nelt)
integer i
real*8 tmplat(maxelt)
neltypes = nelt
do i = 1,nelt
if (lat(i).eq.0) then
lattce_meam(i,i) = 'fcc'
else if (lat(i).eq.1) then
lattce_meam(i,i) = 'bcc'
else if (lat(i).eq.2) then
lattce_meam(i,i) = 'hcp'
else if (lat(i).eq.3) then
lattce_meam(i,i) = 'dim'
else if (lat(i).eq.4) then
lattce_meam(i,i) = 'dia'
else
c unknown
endif
Z_meam(i) = z(i)
ielt_meam(i) = ielement(i)
alpha_meam(i,i) = alpha(i)
beta0_meam(i) = b0(i)
beta1_meam(i) = b1(i)
beta2_meam(i) = b2(i)
beta3_meam(i) = b3(i)
tmplat(i) = alat(i)
Ec_meam(i,i) = esub(i)
A_meam(i) = asub(i)
t0_meam(i) = t0(i)
t1_meam(i) = t1(i)
t2_meam(i) = t2(i)
t3_meam(i) = t3(i)
rho0_meam(i) = rozero(i)
ibar_meam(i) = ibar(i)
if (lattce_meam(i,i).eq.'fcc') then
re_meam(i,i) = tmplat(i)/sqrt(2.d0)
elseif (lattce_meam(i,i).eq.'bcc') then
re_meam(i,i) = tmplat(i)*sqrt(3.d0)/2.d0
elseif (lattce_meam(i,i).eq.'hcp') then
re_meam(i,i) = tmplat(i)
elseif (lattce_meam(i,i).eq.'dim') then
re_meam(i,i) = tmplat(i)
elseif (lattce_meam(i,i).eq.'dia') then
re_meam(i,i) = tmplat(i)*sqrt(3.d0)/4.d0
else
c error
endif
enddo
c Set some defaults
rc_meam = 4.0
delr_meam = 0.1
attrac_meam(:,:) = 0.0
repuls_meam(:,:) = 0.0
Cmax_meam(:,:,:) = 2.8
Cmin_meam(:,:,:) = 2.0
ebound_meam(:,:) = (2.8d0**2)/(4.d0*(2.8d0-1.d0))
delta_meam(:,:) = 0.0
nn2_meam(:,:) = 0
zbl_meam(:,:) = 1
gsmooth_factor = 99.0
augt1 = 1
ialloy = 0
mix_ref_t = 0
emb_lin_neg = 0
bkgd_dyn = 0
erose_form = 0
return
end

View File

@ -1,204 +0,0 @@
c
c do a sanity check on index parameters
subroutine meam_checkindex(num,lim,nidx,idx,ierr)
implicit none
integer i,num,lim,nidx,idx(3),ierr
ierr = 0
if (nidx.lt.num) then
ierr = 2
return
endif
do i=1,num
if ((idx(i).lt.1).or.(idx(i).gt.lim)) then
ierr = 3
return
endif
enddo
end
c
c Declaration in pair_meam.h:
c
c void meam_setup_param(int *, double *, int *, int *, int *);
c
c Call in pair_meam.cpp
c
c meam_setup_param(&which,&value,&nindex,index,&errorflag);
c
c
c
c The "which" argument corresponds to the index of the "keyword" array
c in pair_meam.cpp:
c
c 0 = Ec_meam
c 1 = alpha_meam
c 2 = rho0_meam
c 3 = delta_meam
c 4 = lattce_meam
c 5 = attrac_meam
c 6 = repuls_meam
c 7 = nn2_meam
c 8 = Cmin_meam
c 9 = Cmax_meam
c 10 = rc_meam
c 11 = delr_meam
c 12 = augt1
c 13 = gsmooth_factor
c 14 = re_meam
c 15 = ialloy
c 16 = mixture_ref_t
c 17 = erose_form
c 18 = zbl_meam
c 19 = emb_lin_neg
c 20 = bkgd_dyn
subroutine meam_setup_param(which, value, nindex,
$ index, errorflag)
use meam_data
implicit none
integer which, nindex, index(3), errorflag
real*8 value
integer i1, i2
errorflag = 0
c 0 = Ec_meam
if (which.eq.0) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
Ec_meam(index(1),index(2)) = value
c 1 = alpha_meam
else if (which.eq.1) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
alpha_meam(index(1),index(2)) = value
c 2 = rho0_meam
else if (which.eq.2) then
call meam_checkindex(1,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
rho0_meam(index(1)) = value
c 3 = delta_meam
else if (which.eq.3) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
delta_meam(index(1),index(2)) = value
c 4 = lattce_meam
else if (which.eq.4) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
if (value.eq.0) then
lattce_meam(index(1),index(2)) = "fcc"
else if (value.eq.1) then
lattce_meam(index(1),index(2)) = "bcc"
else if (value.eq.2) then
lattce_meam(index(1),index(2)) = "hcp"
else if (value.eq.3) then
lattce_meam(index(1),index(2)) = "dim"
else if (value.eq.4) then
lattce_meam(index(1),index(2)) = "dia"
else if (value.eq.5) then
lattce_meam(index(1),index(2)) = 'b1'
else if (value.eq.6) then
lattce_meam(index(1),index(2)) = 'c11'
else if (value.eq.7) then
lattce_meam(index(1),index(2)) = 'l12'
else if (value.eq.8) then
lattce_meam(index(1),index(2)) = 'b2'
endif
c 5 = attrac_meam
else if (which.eq.5) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
attrac_meam(index(1),index(2)) = value
c 6 = repuls_meam
else if (which.eq.6) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
repuls_meam(index(1),index(2)) = value
c 7 = nn2_meam
else if (which.eq.7) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
i1 = min(index(1),index(2))
i2 = max(index(1),index(2))
nn2_meam(i1,i2) = value
c 8 = Cmin_meam
else if (which.eq.8) then
call meam_checkindex(3,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
Cmin_meam(index(1),index(2),index(3)) = value
c 9 = Cmax_meam
else if (which.eq.9) then
call meam_checkindex(3,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
Cmax_meam(index(1),index(2),index(3)) = value
c 10 = rc_meam
else if (which.eq.10) then
rc_meam = value
c 11 = delr_meam
else if (which.eq.11) then
delr_meam = value
c 12 = augt1
else if (which.eq.12) then
augt1 = value
c 13 = gsmooth
else if (which.eq.13) then
gsmooth_factor = value
c 14 = re_meam
else if (which.eq.14) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
re_meam(index(1),index(2)) = value
c 15 = ialloy
else if (which.eq.15) then
ialloy = value
c 16 = mixture_ref_t
else if (which.eq.16) then
mix_ref_t = value
c 17 = erose_form
else if (which.eq.17) then
erose_form = value
c 18 = zbl_meam
else if (which.eq.18) then
call meam_checkindex(2,maxelt,nindex,index,errorflag)
if (errorflag.ne.0) return
i1 = min(index(1),index(2))
i2 = max(index(1),index(2))
zbl_meam(i1,i2) = value
c 19 = emb_lin_neg
else if (which.eq.19) then
emb_lin_neg = value
c 20 = bkgd_dyn
else if (which.eq.20) then
bkgd_dyn = value
else
errorflag = 1
endif
return
end

View File

@ -1 +0,0 @@
../Install.py

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.gfortran
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = g95
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.gfortran
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = gfortran
F90FLAGS = -O3 -fPIC -fno-second-underscore
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.ifort
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = ifort
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1,5 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
reax_SYSINC =
reax_SYSLIB =
reax_SYSPATH =

View File

@ -1,5 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
reax_SYSINC =
reax_SYSLIB = -lgfortran
reax_SYSPATH =

View File

@ -1,6 +0,0 @@
# Settings that the LAMMPS build will import when this package library is used
reax_SYSINC =
reax_SYSLIB = -lifcore
reax_SYSPATH =

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.empty
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = mpifort
F90FLAGS = -O3 -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.pgf90
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = pgf90
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.ifort
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = mpif90
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1 +0,0 @@
Makefile.gfortran

View File

@ -1,51 +0,0 @@
# *
# *_________________________________________________________________________*
# * Fortran Library for Reactive Force Field *
# * DESCRIPTION: SEE READ-ME *
# * FILE NAME: Makefile *
# * CONTRIBUTING AUTHORS: Hansohl Cho(MIT), Aidan Thompson(SNL) *
# * and Greg Wagner(SNL) *
# * CONTACT: hansohl@mit.edu, athompson@sandia.gov, gjwagne@sandia.gov *
# *_________________________________________________________________________*/
SHELL = /bin/sh
# which file will be copied to Makefile.lammps
EXTRAMAKE = Makefile.lammps.ifort
# ------ FILES ------
SRC = reax_connect.F reax_inout.F reax_lammps.F reax_poten.F reax_reac.F reax_charges.F
HEADERFILES = reax_defs.h *.blk
# ------ DEFINITIONS ------
LIB = libreax.a
OBJ = $(SRC:.F=.o)
# ------ SETTINGS ------
F90 = mpif90
F90FLAGS = -O -fPIC
ARCHIVE = ar
ARCHFLAG = -rc
USRLIB =
SYSLIB =
# ------ MAKE PROCEDURE ------
lib: $(OBJ)
$(ARCHIVE) $(ARFLAGS) $(LIB) $(OBJ)
@cp $(EXTRAMAKE) Makefile.lammps
# ------ COMPILE RULES ------
%.o:%.F $(HEADERFILES)
$(F90) $(F90FLAGS) -c $<
# ------ CLEAN ------
clean:
-rm *.o $(LIB)

View File

@ -1,78 +0,0 @@
ReaxFF library
Aidan Thompson, Sandia National Labs
athomps at sandia.gov
Jan 2008
This library is an implementation of the ReaxFF potential,
specifically designed to work with LAMMPS. It is derived from Adri van
Duin's original serial code, with intervening incarnations in CMDF and
GRASP.
-------------------------------------------------
This directory has source files to build a library that LAMMPS
links against when using the REAX package.
This library must be built with a F90 compiler, before LAMMPS is
built, so LAMMPS can link against it.
You can type "make lib-reax" from the src directory to see help on how
to build this library via make commands, or you can do the same thing
by typing "python Install.py" from within this directory, or you can
do it manually by following the instructions below.
Build the library using one of the provided Makefile.* files or create
your own, specific to your compiler and system. For example:
make -f Makefile.gfortran
When you are done building this library, two files should
exist in this directory:
libreax.a the library LAMMPS will link against
Makefile.lammps settings the LAMMPS Makefile will import
Makefile.lammps is created by the make command, by copying one of the
Makefile.lammps.* files. See the EXTRAMAKE setting at the top of the
Makefile.* files.
IMPORTANT: You must examine the final Makefile.lammps to insure it is
correct for your system, else the LAMMPS build will likely fail.
Makefile.lammps has settings for 3 variables:
user-reax_SYSINC = leave blank for this package
user-reax_SYSLIB = auxiliary F90 libs needed to link a F90 lib with
a C++ program (LAMMPS) via a C++ compiler
user-reax_SYSPATH = path(s) to where those libraries are
Because you have a F90 compiler on your system, you should have these
libraries. But you will have to figure out which ones are needed and
where they are. Examples of common configurations are in the
Makefile.lammps.* files.
-------------------------------------------------
Additional build notes:
The include file reax_defs.h is used by both the ReaxFF library source
files and the LAMMPS pair_reax.cpp source file (in package src/REAX).
It contains dimensions of statically-allocated arrays created by the
ReaxFF library. The size of these arrays must be set small enough to
avoid exceeding the available machine memory, and large enough to fit
the actual data generated by ReaxFF. If you change the values in
reax_defs.h, you must first rebuild the library and then rebuild
LAMMPS.
This library is called by functions in pair_reax.cpp. The C++ to
FORTRAN function calls in pair_reax.cpp assume that FORTRAN object
names are converted to C object names by appending an underscore
character. This is generally the case, but on machines that do not
conform to this convention, you will need to modify either the C++
code or your compiler settings. The name conversion is handled by the
preprocessor macro called FORTRAN in the file pair_reax_fortran.h,
which is included by pair_reax.cpp. Different definitions of this
macro can be obtained by adding a machine-specific macro definition to
the CCFLAGS variable in your your LAMMPS Makefile e.g. -D_IBM. See
pair_reax_fortran.h for more info.

View File

@ -1,116 +0,0 @@
#include "reax_defs.h"
implicit real*8 (a-h,o-z),integer(i-n)
parameter (nneighmax=NNEIGHMAXDEF)
parameter (nat=NATDEF) !Max number of atoms
parameter (nattot=NATTOTDEF) !Max number of global atoms
parameter (nsort=NSORTDEF) !Max number of atom types
parameter (mbond=MBONDDEF) !Max number of bonds connected to one atom
parameter (na1mx3=3*nat) !3*max number of atoms
parameter (navib=NAVIBDEF) !for 2nd derivatives
parameter (nbotym=NBOTYMDEF) !Max number of bond types
parameter (nvatym=NVATYMDEF) !Max number of valency angle types
parameter (ntotym=NTOTYMDEF) !Max number of torsion angle types
parameter (nhbtym=NHBTYMDEF) !Max number of hydrogen bond types
parameter (nodmtym=NODMTYMDEF) !Max number of off-diagonal Morse types
parameter (nboallmax=NBOALLMAXDEF) !Max number of all bonds
parameter (nbomax=NBOMAXDEF) !Max number of bonds
parameter (nhbmax=NHBMAXDEF) !Max number of hydrogen bonds
parameter (nvamax=NVAMAXDEF) !Max number of valency angles
parameter (nopmax=NOPMAXDEF) !Max number of out of plane angles
parameter (ntomax=NTOMAXDEF) !Max number of torsion angles
parameter (npamax=NPAMAXDEF) !Max number of general parameters in force field
parameter (nmolmax=NMOLMAXDEF) !Max number of molecules in system
parameter (nmolset=NMOLSETDEF) !Max number of molecules in training set
parameter (mrestra=MRESTRADEF) !Max number of restraints
parameter (mtreg=MTREGDEF) !Max number of temperature regimes
parameter (mtzone=MTZONEDEF) !Max number of temperature zones
parameter (mvreg=MVREGDEF) !Max number of volume regimes
parameter (mvzone=MVZONEDEF) !Max number of volume zones
parameter (mereg=MEREGDEF) !Max number of electric field regimes
parameter (mezone=MEZONEDEF) !Max number of electric field zones
character*1 qr,qrset,qresi2
character*2 qaset,qadd
character*3 qresi1
character*5 qlabel,qffty,qbgfaxes,qbgfsgn,qresi3
character*20 qkeyw
character*25 qfile
character*40 qffield,qformat,qstrana2
character*60 qremark,qremset,qmolset
character*200 qstrana1
common
$/cbka/ dhbdc(nhbmax,3,3),cp(nat,3),
$ cadd(nat,3),d2(3*navib,3*navib),
$ veladd(3,nat),
$ aold(3,nat),dic(3,nat),pvdw1(nsort,nsort),
$ pvdw2(nsort,nsort),angimp(nat,6),
$ yt(na1mx3),pt(na1mx3),gi(na1mx3),enmolset(nmolset),
$ ai(na1mx3),bi(na1mx3),yi(na1mx3),pn(na1mx3),tbo(nat),
$ chgbgf(nattot),
$ abo2(nat),bor4(nat),bos(nbomax),
$ eldef(nat),vradic(nat),
$ vmo2(nat),
$ ro(nbomax),dbondr(nbomax),
$ dbosidr(nbomax),thgo(nopmax),
$ elmol(nmolmax),
$ elaf(nsort),vpq(nsort),
$ rvdw(nsort),alf(nsort),eps(nsort),chat(nsort),
$ rcore(nsort,nsort),ecore(nsort,nsort),acore(nsort,nsort),
$ vlp2(nsort),
$ valp2(nsort),vincr(nsort),
$ vval3(nsort),
$ vuncor(nbotym),
$ vop(nsort),
$ sigqeq(nsort),
$ rrcha(mrestra),
$ rmstra3(mrestra),
$ rmcha(mrestra),
$ rtcha(mrestra),rvcha(mrestra),
$ v2bo(ntotym),v3bo(ntotym),
$ eel,fctor,elr,
$ presx2,presy2,presz2,
$ tset2,
$ enmol,formol,vvol,tpnrad,
$ delvib,
$ taut2,tincr,xmasmd,
$ gdicmax,parc1,parc2,sumelec,
$ xinh,fsnh,vqnh,snh,ham,errnh,sumhe,
$ swa,swb2,swc0,swc1,swc2,swc3,swc4,swc5,swc6,
$ swc7,plr,endpo2,ccpar,
$ c4,estrmin,endpo,accincr,
$ endpoold,xadd,yadd,zadd,addist,taddmol,
$ Hug_E0, Hug_P0, Hug_V0, xImpVcm, shock_vel,
$ shock_z_sep
common
$/cbka/
$ ioop(nopmax,9),ifreqset(nmolset),
$ ijk(nat,4),icgeopt(nmolset),
$ irap(50),irdo(50,2),
$ ityadd(nat),
$ nmoloo(nat),iradic(nat),idef(nsort),nasort(nsort),
$ ibgr1(nattot),ibgr2(nattot),idupc(6),
$ imolsta(nat),
$ ncent2(nbomax),irads,nrdd,nrddf,nbiolab,nuge,
$ nbon2,npar,nodmty,ngnh,irac,nincrop,
$ nboty,mdstep,
$ nreac,
$ nbonop,icelo2,
$ iaddfreq,iveladd,invt,
$ noop,ndtau,
$ nelc3,nfc,nsav2,nmmax,ibh2,
$ nmmaxold,nfcold,icellold,imodfile,
$ icelo2old,inmov1,inmov2,nchaold,naa,nadattempt,
$ nequi,iadj,
$ ntest,nmm,
$ nmolo5o,nradcount,nmollset,iflga,
$ iperiod,ibgfversion,iremark,iconne,
$ kx,ky,kz,iexco,iruid,ibity,nvlist,
$ ityrad,iredo,iexx,iexy,iexz,ncellopt,
$ ndata2,nprob,nit,i5758,ingeo,nmoloold,itemp,
$ icgeo,ishock_type,isymm,
$ qadd(nat),qlabel(nattot),qffty(nattot),qresi1(nattot),
$ qresi2(nattot),qresi3(nattot),
$ qremark(20),qformat(20),qr,qffield,
$ qstrana1,qstrana2,qmolset(nmolset)
***********************************************************************

View File

@ -1,4 +0,0 @@
common
$/cbkabo/ abo(nat)

View File

@ -1,3 +0,0 @@
common
$/cbkatomcoord/ id(nat,3),xmasat(nat),vel(3,nat),accel(3,nat)

View File

@ -1,3 +0,0 @@
common
$/cbkbo/ bo(nbomax)

View File

@ -1,5 +0,0 @@
common
$/cbkboncor/ dbosindc(nbomax,3,2*mbond+2),dbosidc(nbomax,3,2),
$ bo131(nsort),bo132(nsort),bo133(nsort),
$ ovc(nbotym),v13cor(nbotym)

View File

@ -1,3 +0,0 @@
common
$/cbkbopi/ bopi(nbomax)

View File

@ -1,3 +0,0 @@
common
$/cbkbopi2/ bopi2(nbomax)

View File

@ -1,4 +0,0 @@
common
$/cbkbosi/ bosi(nbomax)

View File

@ -1,5 +0,0 @@
common
$/cbkc/ c(nat,3),cglobal(nattot,3),itag(nat),
$chgglobal(nattot)

View File

@ -1,4 +0,0 @@
common
$/cbkch/ ch(nat)

View File

@ -1,5 +0,0 @@
common
$/cbkcha/ ech,syscha,chisys
$ vfieldx,vfieldy,vfieldz,nmcharge,ioldchg

View File

@ -1,4 +0,0 @@
common
$/cbkcharmol/ vmcha(nmolmax),
$ iat1mc(nmolmax),iat2mc(nmolmax)

View File

@ -1,3 +0,0 @@
common
$/cbkchb/ chi(nsort),eta(nsort),gam(nsort)

View File

@ -1,5 +0,0 @@
common
$/cbkconst/ dgrrdn,one,half,three,zero,caljou,rgasc,xjouca
$ convmd

View File

@ -1,7 +0,0 @@
common
$/cbkcovbon/ de2(nbotym),de3(nbotym),psi(nbotym),
$ psp(nbotym),
$ ltripstaball

View File

@ -1,7 +0,0 @@
integer Lvirial,Latomvirial
common
$/cbkd/ d(3,nat),estrain(nat)
common
$/cbkvirial/ atomvirial(6,nat),virial(6),Lvirial,Latomvirial

View File

@ -1,3 +0,0 @@
common
$/cbkdbodc/ dbodc(nbomax,3,2)

View File

@ -1,6 +0,0 @@
common
$/cbkdbopi2ndc/ dbopi2ndc(nbomax,3,2*mbond+2)

View File

@ -1,5 +0,0 @@
common
$/dbopidc/ dbopi2dc(nbomax,3,2),dbopidc(nbomax,3,2)

View File

@ -1,6 +0,0 @@
common
$/dbopindc/ dbopindc(nbomax,3,2*mbond+2)

View File

@ -1,5 +0,0 @@
common
$/cbkdcell/ dcell(3,nat,27)

View File

@ -1,5 +0,0 @@
common
$/cbkdhdc/ dhdc(nvamax,3,3)

View File

@ -1,4 +0,0 @@
common
$/cbkdistan/ axis(3),aaxh,baxh,caxh,iortho

View File

@ -1,5 +0,0 @@
common
$/cbkdrdc/ drdc(nbomax,3,2)

View File

@ -1,4 +0,0 @@
common
$/cbkefield/ efix,efiy,efiz,c1

Some files were not shown because too many files have changed in this diff Show More