more changes from Steve
This commit is contained in:
Steven J. Plimpton
@ -149,6 +149,7 @@ Package, Description, Doc page, Example, Library
|
||||
"USER-QTB"_#USER-QTB, quantum nuclear effects,"fix qtb"_fix_qtb.html "fix qbmsst"_fix_qbmsst.html, qtb, -
|
||||
"USER-QUIP"_#USER-QUIP, QUIP/libatoms interface,"pair_style quip"_pair_quip.html, USER/quip, ext
|
||||
"USER-REAXC"_#USER-REAXC, ReaxFF potential (C/C++) ,"pair_style reaxc"_pair_reaxc.html, reax, -
|
||||
"USER-SCAFACOS"_#USER-SCAFACOS, ScaFaCoS long-range Coulombics,"kspace_style scafacos"_kspace_style.html, USER/scafacos, ext
|
||||
"USER-SMD"_#USER-SMD, smoothed Mach dynamics,"SMD User Guide"_PDF/SMD_LAMMPS_userguide.pdf, USER/smd, ext
|
||||
"USER-SMTBQ"_#USER-SMTBQ, second moment tight binding QEq potential,"pair_style smtbq"_pair_smtbq.html, USER/smtbq, -
|
||||
"USER-SPH"_#USER-SPH, smoothed particle hydrodynamics,"SPH User Guide"_PDF/SPH_LAMMPS_userguide.pdf, USER/sph, -
|
||||
@ -2730,6 +2731,60 @@ examples/reax :ul
|
||||
|
||||
:line
|
||||
|
||||
RENE: check this section
|
||||
|
||||
USER-SCAFACOS package :link(USER-SCAFACOS),h4
|
||||
|
||||
[Contents:]
|
||||
|
||||
A KSpace style which wraps the "ScaFaCoS Coulomb solver
|
||||
library"_http://www.scafacos.de.
|
||||
|
||||
To use this package you must have the ScaFaCoS library available on
|
||||
your system. It is available for download at "http://scafacos.de" or
|
||||
can be cloned from the git-repository
|
||||
"git://github.com/scafacos/scafacos.git".
|
||||
|
||||
[Author:] Rene Halver (JSC) wrote the scafacos LAMMPS command.
|
||||
|
||||
ScaFaCoS was developed by a consortium of German research facilities
|
||||
with a BMBF (German Ministry of Science and Education) funded project
|
||||
in 2009-2012. Participants of the consortium were the Universities of
|
||||
Bonn, Chemnitz, Stuttgart, and Wuppertal as well as the
|
||||
Forschungszentrum Juelich.
|
||||
|
||||
[Install or un-install:]
|
||||
|
||||
Before building LAMMPS with this package, you must first download and
|
||||
build the ScaFaCoS library. You can do this manually if you prefer;
|
||||
follow the instructions in lib/scafacos/README. You can also do it in
|
||||
one step from the lammps/src dir, using a command like these, which
|
||||
simply invoke the lib/scafacos/Install.py script with the specified
|
||||
args:
|
||||
|
||||
make lib-scafacos # print help message
|
||||
make lib-scafacos args="-b" # download and build in lib/scafacos/scafacos-<version>
|
||||
make lib-voronoi args="-p $HOME/scafacos # use existing ScaFaCoS installation in $HOME/scafacos
|
||||
|
||||
You can then install/un-install the package and build LAMMPS in the
|
||||
usual manner:
|
||||
|
||||
make yes-user-scafacos
|
||||
make machine :pre
|
||||
|
||||
make no-user-scafacos
|
||||
make machine :pre
|
||||
|
||||
[Supporting info:]
|
||||
|
||||
src/USER-SCAFACOS: filenames -> commands
|
||||
src/USER-SCAFACOS/README
|
||||
"kspace_style scafacos"_kspace_style.html
|
||||
"kspace_modify"_kspace_modify.html
|
||||
examples/USER/scafacos :ul
|
||||
|
||||
:line
|
||||
|
||||
USER-SMD package :link(USER-SMD),h4
|
||||
|
||||
[Contents:]
|
||||
|
||||
@ -1,44 +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,Section_commands.html#comm)
|
||||
|
||||
:line
|
||||
|
||||
fix scafacos command :h3
|
||||
|
||||
[Syntax:]
|
||||
|
||||
fix ID group-ID scafacos solver tol_type tolerance :pre
|
||||
|
||||
ID, group-ID are documented in "fix"_fix.html command
|
||||
scafacos = style name of this fix command
|
||||
solver = NULL or name of solver for electro-statics computation :ul
|
||||
tol_type = NULL or one of (energy, energy_rel, field, field_rel, potential, potential_rel)
|
||||
tolerance = only of tol_type is given, value of tolerance
|
||||
|
||||
[Examples:]
|
||||
|
||||
fix 1 all scafacos fmm
|
||||
fix 1 all scafacos p3m field 0.001
|
||||
|
||||
[Description:]
|
||||
|
||||
This fix style is a wrapper for the Coulomb solver library ScaFaCoS,
|
||||
which provides a selection of different solvers for the compuation
|
||||
of electro-static interaction. If you download and build ScaFaCoS,
|
||||
it can be called as libray by LAMMPS via this fix, which allows the
|
||||
selection of a solver as well as setting the tolerance of a chosen
|
||||
parameter (if compatible with the chosen solver).
|
||||
|
||||
ScaFaCoS was developed by a consortium of German research facilities
|
||||
with a BMBF (German Ministry of Science and Education) funded project
|
||||
in 2009-2012. Participants of the consortium were the Universities of
|
||||
Bonn, Chemnitz, Stuttgart, and Wuppertal as well as the
|
||||
Forschungszentrum Juelich.
|
||||
|
||||
The library is available at "http://scafacos.de" or can be cloned
|
||||
from the git-repository "git://github.com/scafacos/scafacos.git".
|
||||
|
||||
:line
|
||||
@ -13,47 +13,51 @@ kspace_modify command :h3
|
||||
kspace_modify keyword value ... :pre
|
||||
|
||||
one or more keyword/value pairs may be listed :ulb,l
|
||||
keyword = {mesh} or {order} or {order/disp} or {mix/disp} or {overlap} or {minorder} or {force} or {gewald} or {gewald/disp} or {slab} or (nozforce} or {compute} or {cutoff/adjust} or {fftbench} or {collective} or {diff} or {kmax/ewald} or {force/disp/real} or {force/disp/kspace} or {splittol} or {disp/auto}:l
|
||||
{mesh} value = x y z
|
||||
x,y,z = grid size in each dimension for long-range Coulombics
|
||||
{mesh/disp} value = x y z
|
||||
x,y,z = grid size in each dimension for 1/r^6 dispersion
|
||||
{order} value = N
|
||||
N = extent of Gaussian for PPPM or MSM mapping of charge to grid
|
||||
{order/disp} value = N
|
||||
N = extent of Gaussian for PPPM mapping of dispersion term to grid
|
||||
{mix/disp} value = {pair} or {geom} or {none}
|
||||
{overlap} = {yes} or {no} = whether the grid stencil for PPPM is allowed to overlap into more than the nearest-neighbor processor
|
||||
{minorder} value = M
|
||||
M = min allowed extent of Gaussian when auto-adjusting to minimize grid communication
|
||||
keyword = {collective} or {compute} or {cutoff/adjust} or {diff} or {disp/auto} or {fftbench} or {force/disp/kspace} or {force/disp/real} or {force} or {gewald/disp} or {gewald} or {kmax/ewald} or {mesh} or {minorder} or {mix/disp} or {order/disp} or {order} or {overlap} or {scafacos} or {slab} or {splittol} :l
|
||||
{collective} value = {yes} or {no}
|
||||
{compute} value = {yes} or {no}
|
||||
{cutoff/adjust} value = {yes} or {no}
|
||||
{diff} value = {ad} or {ik} = 2 or 4 FFTs for PPPM in smoothed or non-smoothed mode
|
||||
{disp/auto} value = yes or no
|
||||
{fftbench} value = {yes} or {no}
|
||||
{force/disp/real} value = accuracy (force units)
|
||||
{force/disp/kspace} value = accuracy (force units)
|
||||
{force} value = accuracy (force units)
|
||||
{gewald} value = rinv (1/distance units)
|
||||
rinv = G-ewald parameter for Coulombics
|
||||
{gewald/disp} value = rinv (1/distance units)
|
||||
rinv = G-ewald parameter for dispersion
|
||||
{kmax/ewald} value = kx ky kz
|
||||
kx,ky,kz = number of Ewald sum kspace vectors in each dimension
|
||||
{mesh} value = x y z
|
||||
x,y,z = grid size in each dimension for long-range Coulombics
|
||||
{mesh/disp} value = x y z
|
||||
x,y,z = grid size in each dimension for 1/r^6 dispersion
|
||||
{minorder} value = M
|
||||
M = min allowed extent of Gaussian when auto-adjusting to minimize grid communication
|
||||
{mix/disp} value = {pair} or {geom} or {none}
|
||||
{order} value = N
|
||||
N = extent of Gaussian for PPPM or MSM mapping of charge to grid
|
||||
{order/disp} value = N
|
||||
N = extent of Gaussian for PPPM mapping of dispersion term to grid
|
||||
{overlap} = {yes} or {no} = whether the grid stencil for PPPM is allowed to overlap into more than the nearest-neighbor processor
|
||||
{pressure/scalar} value = {yes} or {no}
|
||||
{scafacos} values = option value1 value2 ...
|
||||
option = {tolerance}
|
||||
value = {energy} or {energy_rel} or {field} or {field_rel} or {potential} or {potential_rel}
|
||||
{slab} value = volfactor or {nozforce}
|
||||
volfactor = ratio of the total extended volume used in the
|
||||
2d approximation compared with the volume of the simulation domain
|
||||
{nozforce} turns off kspace forces in the z direction
|
||||
{compute} value = {yes} or {no}
|
||||
{cutoff/adjust} value = {yes} or {no}
|
||||
{pressure/scalar} value = {yes} or {no}
|
||||
{fftbench} value = {yes} or {no}
|
||||
{collective} value = {yes} or {no}
|
||||
{diff} value = {ad} or {ik} = 2 or 4 FFTs for PPPM in smoothed or non-smoothed mode
|
||||
{kmax/ewald} value = kx ky kz
|
||||
kx,ky,kz = number of Ewald sum kspace vectors in each dimension
|
||||
{force/disp/real} value = accuracy (force units)
|
||||
{force/disp/kspace} value = accuracy (force units)
|
||||
{splittol} value = tol
|
||||
tol = relative size of two eigenvalues (see discussion below)
|
||||
{disp/auto} value = yes or no :pre
|
||||
tol = relative size of two eigenvalues (see discussion below) :pre
|
||||
:ule
|
||||
|
||||
[Examples:]
|
||||
|
||||
kspace_modify mesh 24 24 30 order 6
|
||||
kspace_modify slab 3.0 :pre
|
||||
kspace_modify slab 3.0
|
||||
kspace_modify scafacos tolerance energy :pre
|
||||
|
||||
[Description:]
|
||||
|
||||
@ -61,6 +65,132 @@ Set parameters used by the kspace solvers defined by the
|
||||
"kspace_style"_kspace_style.html command. Not all parameters are
|
||||
relevant to all kspace styles.
|
||||
|
||||
:line
|
||||
|
||||
The {collective} keyword applies only to PPPM. It is set to {no} by
|
||||
default, except on IBM BlueGene machines. If this option is set to
|
||||
{yes}, LAMMPS will use MPI collective operations to remap data for
|
||||
3d-FFT operations instead of the default point-to-point communication.
|
||||
This is faster on IBM BlueGene machines, and may also be faster on
|
||||
other machines if they have an efficient implementation of MPI
|
||||
collective operations and adequate hardware.
|
||||
|
||||
:line
|
||||
|
||||
The {compute} keyword allows Kspace computations to be turned off,
|
||||
even though a "kspace_style"_kspace_style.html is defined. This is
|
||||
not useful for running a real simulation, but can be useful for
|
||||
debugging purposes or for computing only partial forces that do not
|
||||
include the Kspace contribution. You can also do this by simply not
|
||||
defining a "kspace_style"_kspace_style.html, but a Kspace-compatible
|
||||
"pair_style"_pair_style.html requires a kspace style to be defined.
|
||||
This keyword gives you that option.
|
||||
|
||||
:line
|
||||
|
||||
The {cutoff/adjust} keyword applies only to MSM. If this option is
|
||||
turned on, the Coulombic cutoff will be automatically adjusted at the
|
||||
beginning of the run to give the desired estimated error. Other
|
||||
cutoffs such as LJ will not be affected. If the grid is not set using
|
||||
the {mesh} command, this command will also attempt to use the optimal
|
||||
grid that minimizes cost using an estimate given by
|
||||
"(Hardy)"_#Hardy1. Note that this cost estimate is not exact, somewhat
|
||||
experimental, and still may not yield the optimal parameters.
|
||||
|
||||
:line
|
||||
|
||||
The {diff} keyword specifies the differentiation scheme used by the
|
||||
PPPM method to compute forces on particles given electrostatic
|
||||
potentials on the PPPM mesh. The {ik} approach is the default for
|
||||
PPPM and is the original formulation used in "(Hockney)"_#Hockney1. It
|
||||
performs differentiation in Kspace, and uses 3 FFTs to transfer each
|
||||
component of the computed fields back to real space for total of 4
|
||||
FFTs per timestep.
|
||||
|
||||
The analytic differentiation {ad} approach uses only 1 FFT to transfer
|
||||
information back to real space for a total of 2 FFTs per timestep. It
|
||||
then performs analytic differentiation on the single quantity to
|
||||
generate the 3 components of the electric field at each grid point.
|
||||
This is sometimes referred to as "smoothed" PPPM. This approach
|
||||
requires a somewhat larger PPPM mesh to achieve the same accuracy as
|
||||
the {ik} method. Currently, only the {ik} method (default) can be
|
||||
used for a triclinic simulation cell with PPPM. The {ad} method is
|
||||
always used for MSM.
|
||||
|
||||
NOTE: Currently, not all PPPM styles support the {ad} option. Support
|
||||
for those PPPM variants will be added later.
|
||||
|
||||
:line
|
||||
|
||||
The {disp/auto} option controls whether the pppm/disp is allowed to
|
||||
generate PPPM parameters automatically. If set to {no}, parameters have
|
||||
to be specified using the {gewald/disp}, {mesh/disp},
|
||||
{force/disp/real} or {force/disp/kspace} keywords, or
|
||||
the code will stop with an error message. When this option is set to
|
||||
{yes}, the error message will not appear and the simulation will start.
|
||||
For a typical application, using the automatic parameter generation
|
||||
will provide simulations that are either inaccurate or slow. Using this
|
||||
option is thus not recommended. For guidelines on how to obtain good
|
||||
parameters, see the "How-To"_Section_howto.html#howto_24 discussion.
|
||||
|
||||
:line
|
||||
|
||||
The {fftbench} keyword applies only to PPPM. It is off by default. If
|
||||
this option is turned on, LAMMPS will perform a short FFT benchmark
|
||||
computation and report its timings, and will thus finish a some seconds
|
||||
later than it would if this option were off.
|
||||
|
||||
:line
|
||||
|
||||
The {force/disp/real} and {force/disp/kspace} keywords set the force
|
||||
accuracy for the real and space computations for the dispersion part
|
||||
of pppm/disp. As shown in "(Isele-Holder)"_#Isele-Holder1, optimal
|
||||
performance and accuracy in the results is obtained when these values
|
||||
are different.
|
||||
|
||||
:line
|
||||
|
||||
The {force} keyword overrides the relative accuracy parameter set by
|
||||
the "kspace_style"_kspace_style.html command with an absolute
|
||||
accuracy. The accuracy determines the RMS error in per-atom forces
|
||||
calculated by the long-range solver and is thus specified in force
|
||||
units. A negative value for the accuracy setting means to use the
|
||||
relative accuracy parameter. The accuracy setting is used in
|
||||
conjunction with the pairwise cutoff to determine the number of
|
||||
K-space vectors for style {ewald}, the FFT grid size for style
|
||||
{pppm}, or the real space grid size for style {msm}.
|
||||
|
||||
:line
|
||||
|
||||
The {gewald} keyword sets the value of the Ewald or PPPM G-ewald
|
||||
parameter for charge as {rinv} in reciprocal distance units. Without
|
||||
this setting, LAMMPS chooses the parameter automatically as a function
|
||||
of cutoff, precision, grid spacing, etc. This means it can vary from
|
||||
one simulation to the next which may not be desirable for matching a
|
||||
KSpace solver to a pre-tabulated pairwise potential. This setting can
|
||||
also be useful if Ewald or PPPM fails to choose a good grid spacing
|
||||
and G-ewald parameter automatically. If the value is set to 0.0,
|
||||
LAMMPS will choose the G-ewald parameter automatically. MSM does not
|
||||
use the {gewald} parameter.
|
||||
|
||||
:line
|
||||
|
||||
The {gewald/disp} keyword sets the value of the Ewald or PPPM G-ewald
|
||||
parameter for dispersion as {rinv} in reciprocal distance units. It
|
||||
has the same meaning as the {gewald} setting for Coulombics.
|
||||
|
||||
:line
|
||||
|
||||
The {kmax/ewald} keyword sets the number of kspace vectors in each
|
||||
dimension for kspace style {ewald}. The three values must be positive
|
||||
integers, or else (0,0,0), which unsets the option. When this option
|
||||
is not set, the Ewald sum scheme chooses its own kspace vectors,
|
||||
consistent with the user-specified accuracy and pairwise cutoff. In
|
||||
any case, if kspace style {ewald} is invoked, the values used are
|
||||
printed to the screen and the log file at the start of the run.
|
||||
|
||||
:line
|
||||
|
||||
The {mesh} keyword sets the grid size for kspace style {pppm} or
|
||||
{msm}. In the case of PPPM, this is the FFT mesh, and each dimension
|
||||
must be factorizable into powers of 2, 3, and 5. In the case of MSM,
|
||||
@ -70,6 +200,8 @@ or MSM solver chooses its own grid size, consistent with the
|
||||
user-specified accuracy and pairwise cutoff. Values for x,y,z of
|
||||
0,0,0 unset the option.
|
||||
|
||||
:line
|
||||
|
||||
The {mesh/disp} keyword sets the grid size for kspace style
|
||||
{pppm/disp}. This is the FFT mesh for long-range dispersion and ach
|
||||
dimension must be factorizable into powers of 2, 3, and 5. When this
|
||||
@ -77,39 +209,7 @@ option is not set, the PPPM solver chooses its own grid size,
|
||||
consistent with the user-specified accuracy and pairwise cutoff.
|
||||
Values for x,y,z of 0,0,0 unset the option.
|
||||
|
||||
The {order} keyword determines how many grid spacings an atom's charge
|
||||
extends when it is mapped to the grid in kspace style {pppm} or {msm}.
|
||||
The default for this parameter is 5 for PPPM and 8 for MSM, which
|
||||
means each charge spans 5 or 8 grid cells in each dimension,
|
||||
respectively. For the LAMMPS implementation of MSM, the order can
|
||||
range from 4 to 10 and must be even. For PPPM, the minimum allowed
|
||||
setting is 2 and the maximum allowed setting is 7. The larger the
|
||||
value of this parameter, the smaller that LAMMPS will set the grid
|
||||
size, to achieve the requested accuracy. Conversely, the smaller the
|
||||
order value, the larger the grid size will be. Note that there is an
|
||||
inherent trade-off involved: a small grid will lower the cost of FFTs
|
||||
or MSM direct sum, but a larger order parameter will increase the cost
|
||||
of interpolating charge/fields to/from the grid.
|
||||
|
||||
The {order/disp} keyword determines how many grid spacings an atom's
|
||||
dispersion term extends when it is mapped to the grid in kspace style
|
||||
{pppm/disp}. It has the same meaning as the {order} setting for
|
||||
Coulombics.
|
||||
|
||||
The {overlap} keyword can be used in conjunction with the {minorder}
|
||||
keyword with the PPPM styles to adjust the amount of communication
|
||||
that occurs when values on the FFT grid are exchanged between
|
||||
processors. This communication is distinct from the communication
|
||||
inherent in the parallel FFTs themselves, and is required because
|
||||
processors interpolate charge and field values using grid point values
|
||||
owned by neighboring processors (i.e. ghost point communication). If
|
||||
the {overlap} keyword is set to {yes} then this communication is
|
||||
allowed to extend beyond nearest-neighbor processors, e.g. when using
|
||||
lots of processors on a small problem. If it is set to {no} then the
|
||||
communication will be limited to nearest-neighbor processors and the
|
||||
{order} setting will be reduced if necessary, as explained by the
|
||||
{minorder} keyword discussion. The {overlap} keyword is always set to
|
||||
{yes} in MSM.
|
||||
:line
|
||||
|
||||
The {minorder} keyword allows LAMMPS to reduce the {order} setting if
|
||||
necessary to keep the communication of ghost grid point limited to
|
||||
@ -126,6 +226,42 @@ error if the grid communication is non-nearest-neighbor and {overlap}
|
||||
is set to {no}. The {minorder} keyword is not currently supported in
|
||||
MSM.
|
||||
|
||||
:line
|
||||
|
||||
The {mix/disp} keyword selects the mixing rule for the dispersion
|
||||
coefficients. With {pair}, the dispersion coefficients of unlike
|
||||
types are computed as indicated with "pair_modify"_pair_modify.html.
|
||||
With {geom}, geometric mixing is enforced on the dispersion
|
||||
coefficients in the kspace coefficients. When using the arithmetic
|
||||
mixing rule, this will speed-up the simulations but introduces some
|
||||
error in the force computations, as shown in "(Wennberg)"_#Wennberg.
|
||||
With {none}, it is assumed that no mixing rule is
|
||||
applicable. Splitting of the dispersion coefficients will be performed
|
||||
as described in "(Isele-Holder)"_#Isele-Holder1.
|
||||
|
||||
This splitting can be influenced with the {splittol} keywords. Only
|
||||
the eigenvalues that are larger than tol compared to the largest
|
||||
eigenvalues are included. Using this keywords the original matrix of
|
||||
dispersion coefficients is approximated. This leads to faster
|
||||
computations, but the accuracy in the reciprocal space computations of
|
||||
the dispersion part is decreased.
|
||||
|
||||
:line
|
||||
|
||||
The {order} keyword determines how many grid spacings an atom's charge
|
||||
extends when it is mapped to the grid in kspace style {pppm} or {msm}.
|
||||
The default for this parameter is 5 for PPPM and 8 for MSM, which
|
||||
means each charge spans 5 or 8 grid cells in each dimension,
|
||||
respectively. For the LAMMPS implementation of MSM, the order can
|
||||
range from 4 to 10 and must be even. For PPPM, the minimum allowed
|
||||
setting is 2 and the maximum allowed setting is 7. The larger the
|
||||
value of this parameter, the smaller that LAMMPS will set the grid
|
||||
size, to achieve the requested accuracy. Conversely, the smaller the
|
||||
order value, the larger the grid size will be. Note that there is an
|
||||
inherent trade-off involved: a small grid will lower the cost of FFTs
|
||||
or MSM direct sum, but a larger order parameter will increase the cost
|
||||
of interpolating charge/fields to/from the grid.
|
||||
|
||||
The PPPM order parameter may be reset by LAMMPS when it sets up the
|
||||
FFT grid if the implied grid stencil extends beyond the grid cells
|
||||
owned by neighboring processors. Typically this will only occur when
|
||||
@ -134,30 +270,66 @@ be generated indicating the order parameter is being reduced to allow
|
||||
LAMMPS to run the problem. Automatic adjustment of the order parameter
|
||||
is not supported in MSM.
|
||||
|
||||
The {force} keyword overrides the relative accuracy parameter set by
|
||||
the "kspace_style"_kspace_style.html command with an absolute
|
||||
accuracy. The accuracy determines the RMS error in per-atom forces
|
||||
calculated by the long-range solver and is thus specified in force
|
||||
units. A negative value for the accuracy setting means to use the
|
||||
relative accuracy parameter. The accuracy setting is used in
|
||||
conjunction with the pairwise cutoff to determine the number of
|
||||
K-space vectors for style {ewald}, the FFT grid size for style
|
||||
{pppm}, or the real space grid size for style {msm}.
|
||||
:line
|
||||
|
||||
The {gewald} keyword sets the value of the Ewald or PPPM G-ewald
|
||||
parameter for charge as {rinv} in reciprocal distance units. Without
|
||||
this setting, LAMMPS chooses the parameter automatically as a function
|
||||
of cutoff, precision, grid spacing, etc. This means it can vary from
|
||||
one simulation to the next which may not be desirable for matching a
|
||||
KSpace solver to a pre-tabulated pairwise potential. This setting can
|
||||
also be useful if Ewald or PPPM fails to choose a good grid spacing
|
||||
and G-ewald parameter automatically. If the value is set to 0.0,
|
||||
LAMMPS will choose the G-ewald parameter automatically. MSM does not
|
||||
use the {gewald} parameter.
|
||||
The {order/disp} keyword determines how many grid spacings an atom's
|
||||
dispersion term extends when it is mapped to the grid in kspace style
|
||||
{pppm/disp}. It has the same meaning as the {order} setting for
|
||||
Coulombics.
|
||||
|
||||
The {gewald/disp} keyword sets the value of the Ewald or PPPM G-ewald
|
||||
parameter for dispersion as {rinv} in reciprocal distance units. It
|
||||
has the same meaning as the {gewald} setting for Coulombics.
|
||||
:line
|
||||
|
||||
The {overlap} keyword can be used in conjunction with the {minorder}
|
||||
keyword with the PPPM styles to adjust the amount of communication
|
||||
that occurs when values on the FFT grid are exchanged between
|
||||
processors. This communication is distinct from the communication
|
||||
inherent in the parallel FFTs themselves, and is required because
|
||||
processors interpolate charge and field values using grid point values
|
||||
owned by neighboring processors (i.e. ghost point communication). If
|
||||
the {overlap} keyword is set to {yes} then this communication is
|
||||
allowed to extend beyond nearest-neighbor processors, e.g. when using
|
||||
lots of processors on a small problem. If it is set to {no} then the
|
||||
communication will be limited to nearest-neighbor processors and the
|
||||
{order} setting will be reduced if necessary, as explained by the
|
||||
{minorder} keyword discussion. The {overlap} keyword is always set to
|
||||
{yes} in MSM.
|
||||
|
||||
:line
|
||||
|
||||
The {pressure/scalar} keyword applies only to MSM. If this option is
|
||||
turned on, only the scalar pressure (i.e. (Pxx + Pyy + Pzz)/3.0) will
|
||||
be computed, which can be used, for example, to run an isotropic barostat.
|
||||
Computing the full pressure tensor with MSM is expensive, and this option
|
||||
provides a faster alternative. The scalar pressure is computed using a
|
||||
relationship between the Coulombic energy and pressure "(Hummer)"_#Hummer
|
||||
instead of using the virial equation. This option cannot be used to access
|
||||
individual components of the pressure tensor, to compute per-atom virial,
|
||||
or with suffix kspace/pair styles of MSM, like OMP or GPU.
|
||||
|
||||
:line
|
||||
|
||||
RENE: check this section
|
||||
|
||||
The {scafacos} keyword is used for settings that are passed to the
|
||||
ScaFaCoS library when using "kspace_style scafacos"_kspace_style.html.
|
||||
|
||||
The {tolerance} option affects how the {accuracy} specified with
|
||||
the "kspace_style"_kspace_style.html command is interpreted by ScaFaCoS.
|
||||
The following values may be used:
|
||||
|
||||
energy = absolute accuracy in total Coulomic energy
|
||||
energy_rel = relative accuracy in total Coulomic energy
|
||||
field = absolute accuracy in electric field
|
||||
field_rel = relative accuracy in electric field
|
||||
potential = absolute accuracy in total Coulomic potential
|
||||
potential_rel = relative accuracy in total Coulomic potential
|
||||
|
||||
RENE: which one is closest to what LAMMPS used (see kspace_style
|
||||
command) - is it "field"? That's why I set it as the default.
|
||||
|
||||
RENE: How is potential different than energy for this ??
|
||||
|
||||
:line
|
||||
|
||||
The {slab} keyword allows an Ewald or PPPM solver to be used for a
|
||||
systems that are periodic in x,y but non-periodic in z - a
|
||||
@ -191,109 +363,11 @@ the "fix efield"_fix_efield.html command, it will not give the correct
|
||||
dielectric constant due to the Yeh/Berkowitz "(Yeh)"_#Yeh correction
|
||||
not being compatible with how "fix efield"_fix_efield.html works.
|
||||
|
||||
The {compute} keyword allows Kspace computations to be turned off,
|
||||
even though a "kspace_style"_kspace_style.html is defined. This is
|
||||
not useful for running a real simulation, but can be useful for
|
||||
debugging purposes or for computing only partial forces that do not
|
||||
include the Kspace contribution. You can also do this by simply not
|
||||
defining a "kspace_style"_kspace_style.html, but a Kspace-compatible
|
||||
"pair_style"_pair_style.html requires a kspace style to be defined.
|
||||
This keyword gives you that option.
|
||||
:line
|
||||
|
||||
The {cutoff/adjust} keyword applies only to MSM. If this option is
|
||||
turned on, the Coulombic cutoff will be automatically adjusted at the
|
||||
beginning of the run to give the desired estimated error. Other
|
||||
cutoffs such as LJ will not be affected. If the grid is not set using
|
||||
the {mesh} command, this command will also attempt to use the optimal
|
||||
grid that minimizes cost using an estimate given by
|
||||
"(Hardy)"_#Hardy1. Note that this cost estimate is not exact, somewhat
|
||||
experimental, and still may not yield the optimal parameters.
|
||||
The {splittol} keyword is disussed above with the {mix/disp} keyword.
|
||||
|
||||
The {pressure/scalar} keyword applies only to MSM. If this option is
|
||||
turned on, only the scalar pressure (i.e. (Pxx + Pyy + Pzz)/3.0) will
|
||||
be computed, which can be used, for example, to run an isotropic barostat.
|
||||
Computing the full pressure tensor with MSM is expensive, and this option
|
||||
provides a faster alternative. The scalar pressure is computed using a
|
||||
relationship between the Coulombic energy and pressure "(Hummer)"_#Hummer
|
||||
instead of using the virial equation. This option cannot be used to access
|
||||
individual components of the pressure tensor, to compute per-atom virial,
|
||||
or with suffix kspace/pair styles of MSM, like OMP or GPU.
|
||||
|
||||
The {fftbench} keyword applies only to PPPM. It is off by default. If
|
||||
this option is turned on, LAMMPS will perform a short FFT benchmark
|
||||
computation and report its timings, and will thus finish a some seconds
|
||||
later than it would if this option were off.
|
||||
|
||||
The {collective} keyword applies only to PPPM. It is set to {no} by
|
||||
default, except on IBM BlueGene machines. If this option is set to
|
||||
{yes}, LAMMPS will use MPI collective operations to remap data for
|
||||
3d-FFT operations instead of the default point-to-point communication.
|
||||
This is faster on IBM BlueGene machines, and may also be faster on
|
||||
other machines if they have an efficient implementation of MPI
|
||||
collective operations and adequate hardware.
|
||||
|
||||
The {diff} keyword specifies the differentiation scheme used by the
|
||||
PPPM method to compute forces on particles given electrostatic
|
||||
potentials on the PPPM mesh. The {ik} approach is the default for
|
||||
PPPM and is the original formulation used in "(Hockney)"_#Hockney1. It
|
||||
performs differentiation in Kspace, and uses 3 FFTs to transfer each
|
||||
component of the computed fields back to real space for total of 4
|
||||
FFTs per timestep.
|
||||
|
||||
The analytic differentiation {ad} approach uses only 1 FFT to transfer
|
||||
information back to real space for a total of 2 FFTs per timestep. It
|
||||
then performs analytic differentiation on the single quantity to
|
||||
generate the 3 components of the electric field at each grid point.
|
||||
This is sometimes referred to as "smoothed" PPPM. This approach
|
||||
requires a somewhat larger PPPM mesh to achieve the same accuracy as
|
||||
the {ik} method. Currently, only the {ik} method (default) can be
|
||||
used for a triclinic simulation cell with PPPM. The {ad} method is
|
||||
always used for MSM.
|
||||
|
||||
NOTE: Currently, not all PPPM styles support the {ad} option. Support
|
||||
for those PPPM variants will be added later.
|
||||
|
||||
The {kmax/ewald} keyword sets the number of kspace vectors in each
|
||||
dimension for kspace style {ewald}. The three values must be positive
|
||||
integers, or else (0,0,0), which unsets the option. When this option
|
||||
is not set, the Ewald sum scheme chooses its own kspace vectors,
|
||||
consistent with the user-specified accuracy and pairwise cutoff. In
|
||||
any case, if kspace style {ewald} is invoked, the values used are
|
||||
printed to the screen and the log file at the start of the run.
|
||||
|
||||
With the {mix/disp} keyword one can select the mixing rule for the
|
||||
dispersion coefficients. With {pair}, the dispersion coefficients of
|
||||
unlike types are computed as indicated with
|
||||
"pair_modify"_pair_modify.html. With {geom}, geometric mixing is
|
||||
enforced on the dispersion coefficients in the kspace
|
||||
coefficients. When using the arithmetic mixing rule, this will
|
||||
speed-up the simulations but introduces some error in the force
|
||||
computations, as shown in "(Wennberg)"_#Wennberg. With {none}, it is
|
||||
assumed that no mixing rule is applicable. Splitting of the dispersion
|
||||
coefficients will be performed as described in
|
||||
"(Isele-Holder)"_#Isele-Holder1. This splitting can be influenced with
|
||||
the {splittol} keywords. Only the eigenvalues that are larger than tol
|
||||
compared to the largest eigenvalues are included. Using this keywords
|
||||
the original matrix of dispersion coefficients is approximated. This
|
||||
leads to faster computations, but the accuracy in the reciprocal space
|
||||
computations of the dispersion part is decreased.
|
||||
|
||||
The {force/disp/real} and {force/disp/kspace} keywords set the force
|
||||
accuracy for the real and space computations for the dispersion part
|
||||
of pppm/disp. As shown in "(Isele-Holder)"_#Isele-Holder1, optimal
|
||||
performance and accuracy in the results is obtained when these values
|
||||
are different.
|
||||
|
||||
The {disp/auto} option controls whether the pppm/disp is allowed to
|
||||
generate PPPM parameters automatically. If set to {no}, parameters have
|
||||
to be specified using the {gewald/disp}, {mesh/disp},
|
||||
{force/disp/real} or {force/disp/kspace} keywords, or
|
||||
the code will stop with an error message. When this option is set to
|
||||
{yes}, the error message will not appear and the simulation will start.
|
||||
For a typical application, using the automatic parameter generation
|
||||
will provide simulations that are either inaccurate or slow. Using this
|
||||
option is thus not recommended. For guidelines on how to obtain good
|
||||
parameters, see the "How-To"_Section_howto.html#howto_24 discussion.
|
||||
:line
|
||||
|
||||
[Restrictions:] none
|
||||
|
||||
@ -306,10 +380,11 @@ parameters, see the "How-To"_Section_howto.html#howto_24 discussion.
|
||||
The option defaults are mesh = mesh/disp = 0 0 0, order = order/disp =
|
||||
5 (PPPM), order = 10 (MSM), minorder = 2, overlap = yes, force = -1.0,
|
||||
gewald = gewald/disp = 0.0, slab = 1.0, compute = yes, cutoff/adjust =
|
||||
yes (MSM), pressure/scalar = yes (MSM), fftbench = no (PPPM), diff = ik
|
||||
(PPPM), mix/disp = pair, force/disp/real = -1.0, force/disp/kspace = -1.0,
|
||||
split = 0, tol = 1.0e-6, and disp/auto = no. For pppm/intel, order =
|
||||
order/disp = 7.
|
||||
yes (MSM), pressure/scalar = yes (MSM), fftbench = no (PPPM), diff =
|
||||
ik (PPPM), mix/disp = pair, force/disp/real = -1.0, force/disp/kspace
|
||||
= -1.0, split = 0, tol = 1.0e-6, and disp/auto = no. For pppm/intel,
|
||||
order = order/disp = 7. For scafacos settings, scafacos tolerance
|
||||
field.
|
||||
|
||||
:line
|
||||
|
||||
|
||||
@ -12,7 +12,7 @@ kspace_style command :h3
|
||||
|
||||
kspace_style style value :pre
|
||||
|
||||
style = {none} or {ewald} or {ewald/disp} or {ewald/omp} or {pppm} or {pppm/cg} or {pppm/disp} or {pppm/tip4p} or {pppm/stagger} or {pppm/disp/tip4p} or {pppm/gpu} or {pppm/kk} or {pppm/omp} or {pppm/cg/omp} or {pppm/tip4p/omp} or {msm} or {msm/cg} or {msm/omp} or {msm/cg/omp} :ulb,l
|
||||
style = {none} or {ewald} or {ewald/disp} or {ewald/omp} or {pppm} or {pppm/cg} or {pppm/disp} or {pppm/tip4p} or {pppm/stagger} or {pppm/disp/tip4p} or {pppm/gpu} or {pppm/kk} or {pppm/omp} or {pppm/cg/omp} or {pppm/tip4p/omp} or {msm} or {msm/cg} or {msm/omp} or {msm/cg/omp} or {scafacos} :ulb,l
|
||||
{none} value = none
|
||||
{ewald} value = accuracy
|
||||
accuracy = desired relative error in forces
|
||||
@ -22,7 +22,7 @@ style = {none} or {ewald} or {ewald/disp} or {ewald/omp} or {pppm} or {pppm/cg}
|
||||
accuracy = desired relative error in forces
|
||||
{pppm} value = accuracy
|
||||
accuracy = desired relative error in forces
|
||||
{pppm/cg} value = accuracy (smallq)
|
||||
{pppm/cg} values = accuracy (smallq)
|
||||
accuracy = desired relative error in forces
|
||||
smallq = cutoff for charges to be considered (optional) (charge units)
|
||||
{pppm/disp} value = accuracy
|
||||
@ -56,7 +56,10 @@ style = {none} or {ewald} or {ewald/disp} or {ewald/omp} or {pppm} or {pppm/cg}
|
||||
accuracy = desired relative error in forces
|
||||
{msm/cg/omp} value = accuracy (smallq)
|
||||
accuracy = desired relative error in forces
|
||||
smallq = cutoff for charges to be considered (optional) (charge units) :pre
|
||||
smallq = cutoff for charges to be considered (optional) (charge units)
|
||||
{scafacos} values = method accuracy
|
||||
method = fmm or p3m or ... RENE
|
||||
accuracy = desired relative error in forces :pre
|
||||
:ule
|
||||
|
||||
[Examples:]
|
||||
@ -64,6 +67,7 @@ style = {none} or {ewald} or {ewald/disp} or {ewald/omp} or {pppm} or {pppm/cg}
|
||||
kspace_style pppm 1.0e-4
|
||||
kspace_style pppm/cg 1.0e-5 1.0e-6
|
||||
kspace style msm 1.0e-4
|
||||
kspace style scafacos fmm 1.0e-4
|
||||
kspace_style none :pre
|
||||
|
||||
[Description:]
|
||||
@ -210,6 +214,63 @@ pressure simulation with MSM will cause the code to run slower.
|
||||
|
||||
:line
|
||||
|
||||
RENE: edit this section
|
||||
|
||||
The {scafacos} style is a wrapper on the "ScaFaCoS Coulomb solver
|
||||
library"_http://www.scafacos.de which provides a variety of solver
|
||||
methods which can be used with LAMMPS. The paper by "(Who)"_#Who2012
|
||||
gives an overview of ScaFaCoS.
|
||||
|
||||
ScaFaCoS was developed by a consortium of German research facilities
|
||||
with a BMBF (German Ministry of Science and Education) funded project
|
||||
in 2009-2012. Participants of the consortium were the Universities of
|
||||
Bonn, Chemnitz, Stuttgart, and Wuppertal as well as the
|
||||
Forschungszentrum Juelich.
|
||||
|
||||
The library is available for download at "http://scafacos.de" or can
|
||||
be cloned from the git-repository
|
||||
"git://github.com/scafacos/scafacos.git".
|
||||
|
||||
In order to use this KSpace style, you must download and build the
|
||||
ScaFaCoS library, then build LAMMPS with the USER-SCAFACOS package
|
||||
installed package which links LAMMPS to the ScaFaCoS library.
|
||||
See details on "this page"_Section_packages.html#USER-SCAFACOS.
|
||||
|
||||
NOTE: Unlike other KSpace solvers in LAMMPS, ScaFaCoS computes all
|
||||
Coulombic interactions, both short- and long-range. Thus you should
|
||||
NOT use a Coulmbic pair style when using kspace_style scafacos. Thus
|
||||
the total Coulombic energy (short- and long-range) will be tallied for
|
||||
"thermodynamic output"_thermo_style.html command as part of the
|
||||
{elong} keyword; the {ecoul} keyword will be zero.
|
||||
|
||||
NOTE: See the restriction below about use of ScaFaCoS in LAMMPS with
|
||||
molecular charged systems or the TIP4P water model.
|
||||
|
||||
Unlike other KSpace solvers in LAMMPS, ScaFaCoS computes all
|
||||
Coulombic interactions, both short- and long-range. Thus you should
|
||||
NOT use a Coulmbic pair style when using kspace_style scafacos. Thus
|
||||
the total Coulombic energy (short- and long-range) will be tallied for
|
||||
"thermodynamic output"_thermo_style.html command as part of the
|
||||
{elong} keyword; the {ecoul} keyword will be zero.
|
||||
|
||||
The specified {method} determines which ScaFaCoS algorithm is used.
|
||||
The list of methods shown above are the ScaFaCoS methods currently
|
||||
available from LAMMPS.
|
||||
|
||||
RENE: explain this accuracy better compared to what
|
||||
it means for other LAMMPS methods?
|
||||
|
||||
The specified {accuracy} is similar to the accuracy setting for other
|
||||
LAMMPS KSpace styles, but is passed to ScaFaCoS, which can interpret
|
||||
it in different ways for different methods it supports. See the
|
||||
"kspace_modify scafacos accuracy"_kspace_modify.html command for
|
||||
details.
|
||||
|
||||
The "kspace_modify scafacos"_kspace_modify.html command also explains
|
||||
all the ScaFaCoS options currently exposed to LAMMPS.
|
||||
|
||||
:line
|
||||
|
||||
The specified {accuracy} determines the relative RMS error in per-atom
|
||||
forces calculated by the long-range solver. It is set as a
|
||||
dimensionless number, relative to the force that two unit point
|
||||
@ -321,12 +382,24 @@ dimensions. The only exception is if the slab option is set with
|
||||
"kspace_modify"_kspace_modify.html, in which case the xy dimensions
|
||||
must be periodic and the z dimension must be non-periodic.
|
||||
|
||||
The scafacos KSpace style will only be enabled if LAMMPS is built with
|
||||
the USER-SCAFACOS package. See the "Making
|
||||
LAMMPS"_Section_start.html#start_3 section for more info.
|
||||
|
||||
The use of ScaFaCos in LAMMPS does not yet support molecular charged
|
||||
systems where the short-range Coulombic interactions between atoms in
|
||||
the same bond/angle/dihedral are weighted by the
|
||||
"special_bonds"_special_bonds.html command. Likewise it does not
|
||||
support the "TIP4P water style" where a fictitious charge site is
|
||||
introduced in each water molecule.
|
||||
|
||||
[Related commands:]
|
||||
|
||||
"kspace_modify"_kspace_modify.html, "pair_style
|
||||
lj/cut/coul/long"_pair_lj.html, "pair_style
|
||||
lj/charmm/coul/long"_pair_charmm.html, "pair_style
|
||||
lj/long/coul/long"_pair_lj_long.html, "pair_style buck/coul/long"_pair_buck.html
|
||||
lj/long/coul/long"_pair_lj_long.html, "pair_style
|
||||
buck/coul/long"_pair_buck.html
|
||||
|
||||
[Default:]
|
||||
|
||||
@ -384,5 +457,11 @@ Evaluation of Forces for the Simulation of Biomolecules, University of
|
||||
Illinois at Urbana-Champaign, (2006).
|
||||
|
||||
:link(Hardy2009)
|
||||
[(Hardy2)] Hardy, Stone, Schulten, Parallel Computing 35 (2009)
|
||||
164-177.
|
||||
[(Hardy2)] Hardy, Stone, Schulten, Parallel Computing, 35, 164-177
|
||||
(2009).
|
||||
|
||||
RENE:
|
||||
|
||||
:link(Who2012)
|
||||
[(Who)] Who, Author2, Author3, J of Long Range Solvers, 35, 164-177
|
||||
(2012).
|
||||
|
||||
10
examples/USER/scafacos/README
Normal file
10
examples/USER/scafacos/README
Normal file
@ -0,0 +1,10 @@
|
||||
RENE: this dir needs a few short-running examples,
|
||||
for both 1 and 4 procs, of using different Scafacos methods.
|
||||
with sample log file outputs - so users can verify
|
||||
they are running it correctly
|
||||
|
||||
For each input script, there should be files like
|
||||
|
||||
in.scafacos.fmm
|
||||
log.20Jul18.scafacos.fmm.g++.1
|
||||
log.20Jul18.scafacos.fmm.g++.4
|
||||
@ -29,7 +29,7 @@ make lib-scafacos args="-p $HOME/scafacos" # use existing Scafacos installation
|
||||
|
||||
# settings
|
||||
|
||||
#version = "voro++-0.4.6"
|
||||
version = "scafacos-1.0"
|
||||
#url = "http://math.lbl.gov/voro++/download/dir/%s.tar.gz" % version
|
||||
|
||||
# print error message or help
|
||||
@ -89,7 +89,6 @@ def geturl(url,fname):
|
||||
|
||||
args = sys.argv[1:]
|
||||
nargs = len(args)
|
||||
#if nargs == 0: error()
|
||||
|
||||
homepath = "."
|
||||
|
||||
@ -114,8 +113,7 @@ while iarg < nargs:
|
||||
else: error()
|
||||
|
||||
homepath = fullpath(homepath)
|
||||
#homedir = "%s/%s" % (homepath,version)
|
||||
homedir = homepath
|
||||
homedir = "%s/%s" % (homepath,version)
|
||||
|
||||
if (pathflag):
|
||||
if not os.path.isdir(scafacospath): error("Scafacos path does not exist")
|
||||
@ -124,16 +122,13 @@ if (pathflag):
|
||||
if (buildflag and pathflag):
|
||||
error("Cannot use -b and -p flag at the same time")
|
||||
|
||||
#if (not buildflag and not pathflag):
|
||||
# error("Have to use either -b or -p flag")
|
||||
|
||||
# download and unpack Scafacos tarball
|
||||
|
||||
if buildflag:
|
||||
print("Downloading Scafacos ...")
|
||||
geturl(url,"%s/%s.tar.gz" % (homepath,version))
|
||||
|
||||
print("Unpacking Voro++ tarball ...")
|
||||
print("Unpacking Scafacos tarball ...")
|
||||
if os.path.exists("%s/%s" % (homepath,version)):
|
||||
cmd = 'rm -rf "%s/%s"' % (homepath,version)
|
||||
subprocess.check_output(cmd,stderr=subprocess.STDOUT,shell=True)
|
||||
@ -162,7 +157,7 @@ if linkflag:
|
||||
os.remove("includelink")
|
||||
if os.path.isfile("liblink") or os.path.islink("liblink"):
|
||||
os.remove("liblink")
|
||||
cmd = 'ln -s "scafacos/include" includelink'
|
||||
cmd = 'ln -s "%s/include" includelink' % homedir
|
||||
subprocess.check_output(cmd,stderr=subprocess.STDOUT,shell=True)
|
||||
cmd = 'ln -s "scafacos/lib" liblink'
|
||||
cmd = 'ln -s "%s/lib" liblink' % homedir
|
||||
subprocess.check_output(cmd,stderr=subprocess.STDOUT,shell=True)
|
||||
|
||||
5
src/USER-SCAFACOS/README
Normal file
5
src/USER-SCAFACOS/README
Normal file
@ -0,0 +1,5 @@
|
||||
RENE: need a short README file - see other examples in the
|
||||
USER package dirs within src
|
||||
|
||||
give short background on ScaFaCos, authors, contact info
|
||||
for LAMMPS users to have if needed, citation
|
||||
@ -1,456 +0,0 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Rene Halver (JSC)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <cstdio>
|
||||
#include <cstring>
|
||||
#include <cstdlib>
|
||||
#include "fix_scafacos.h"
|
||||
#include "atom.h"
|
||||
#include "comm.h"
|
||||
#include "update.h"
|
||||
#include "neighbor.h"
|
||||
#include "domain.h"
|
||||
#include "force.h"
|
||||
#include "pair.h"
|
||||
#include "modify.h"
|
||||
#include "compute.h"
|
||||
#include "memory.h"
|
||||
#include "error.h"
|
||||
|
||||
// ScaFaCoS library
|
||||
#include <sstream>
|
||||
#include <string>
|
||||
#include "fcs.h"
|
||||
#include "universe.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
using namespace FixConst;
|
||||
|
||||
#define INVOKED_PERATOM 8
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
FixScafacos::FixScafacos(LAMMPS *lmp, int narg, char **arg) :
|
||||
Fix(lmp, narg, arg)
|
||||
{
|
||||
// form of fix string:
|
||||
// fix <#> <group> scafacos <method> [ tolerance <type> <value> ]
|
||||
|
||||
|
||||
if (narg < 4) error->all(FLERR,"Illegal fix scafacos command");
|
||||
|
||||
// read method from fix input
|
||||
method = arg[3];
|
||||
|
||||
int arg_index = 4;
|
||||
|
||||
tolerance_set = false;
|
||||
|
||||
while(arg_index < narg)
|
||||
{
|
||||
// check if tolerance option is set
|
||||
if (strcmp(arg[arg_index],"tolerance") == 0)
|
||||
{
|
||||
tolerance_set = true;
|
||||
++arg_index;
|
||||
if (strcmp(arg[arg_index],"energy") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_ENERGY;
|
||||
else if (strcmp(arg[arg_index],"energy_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_ENERGY_REL;
|
||||
else if (strcmp(arg[arg_index],"field") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD;
|
||||
else if (strcmp(arg[arg_index],"field_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD_REL;
|
||||
else if (strcmp(arg[arg_index],"potential") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_POTENTIAL;
|
||||
else if (strcmp(arg[arg_index],"potential_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_POTENTIAL_REL;
|
||||
else
|
||||
error->all(FLERR,"Illegal fix scafacos command");
|
||||
++arg_index;
|
||||
tolerance_value = atof(arg[arg_index]);
|
||||
++arg_index;
|
||||
}
|
||||
else
|
||||
error->all(FLERR,"Illegal fix scafacos command");
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
FixScafacos::~FixScafacos()
|
||||
{
|
||||
memory->destroy(pot);
|
||||
memory->destroy(field);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
int FixScafacos::setmask()
|
||||
{
|
||||
int mask = 0;
|
||||
//mask |= INITIAL_INTEGRATE;
|
||||
//mask |= FINAL_INTEGRATE;
|
||||
mask |= PRE_REVERSE;
|
||||
mask |= POST_FORCE;
|
||||
mask |= MIN_POST_FORCE;
|
||||
mask |= THERMO_ENERGY;
|
||||
return mask;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::init()
|
||||
{
|
||||
// error checks
|
||||
|
||||
if (domain->dimension == 2)
|
||||
error->all(FLERR,"Fix scafacos requires 3d problem");
|
||||
|
||||
rank = comm->me;
|
||||
|
||||
int nlocal = 0;
|
||||
int proc_grid[3] = {comm->procgrid[0], comm->procgrid[1], comm->procgrid[2]};
|
||||
int myloc[3] = {comm->myloc[0], comm->myloc[1], comm->myloc[2]};
|
||||
|
||||
// call ScaFaCoS init
|
||||
// TODO: check if universe->uworld is a good idea for computations
|
||||
result = fcs_init(&fcs,method.c_str(),universe->uworld);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
setup_handle();
|
||||
|
||||
nlocal = atom -> nlocal;
|
||||
// get particle data
|
||||
x = &atom->x[0][0];
|
||||
q = atom->q;
|
||||
|
||||
if (tolerance_set)
|
||||
{
|
||||
result = fcs_set_tolerance(fcs,tolerance_type,tolerance_value);
|
||||
if (!check_result(result, rank)) return;
|
||||
}
|
||||
|
||||
// print out parameters within handle
|
||||
if (rank == 0) fcs_print_parameters(fcs);
|
||||
|
||||
// call the tuning routine (needs to be redone, if critical system
|
||||
// parameters should change)
|
||||
result = fcs_tune(fcs,nlocal,x,q);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
// allocate arrays larger in order to avoid fast
|
||||
// reallocation due to fluctuations
|
||||
local_array_size = (int)(1.25 * (double)nlocal);
|
||||
|
||||
// allocate new result arrays
|
||||
memory->create(pot,local_array_size,"scafacos:potential");
|
||||
memory->create(field,local_array_size*3,"scafacos:field");
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::init_list(int id, NeighList *ptr)
|
||||
{
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::setup(int vflag)
|
||||
{
|
||||
post_force(vflag);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::min_setup(int vflag)
|
||||
{
|
||||
post_force(vflag);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::setup_pre_reverse(int eflag, int vflag)
|
||||
{
|
||||
pre_reverse(eflag,vflag);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
integrate electronic degrees of freedom
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::initial_integrate(int vflag) {}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
store eflag, so can use it in post_force to tally per-atom energies
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::pre_reverse(int eflag, int vflag)
|
||||
{
|
||||
int eflag_caller = eflag;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::post_force(int vflag)
|
||||
{
|
||||
|
||||
int nlocal;
|
||||
|
||||
nlocal = atom->nlocal;
|
||||
x = &atom->x[0][0];
|
||||
q = atom->q;
|
||||
|
||||
// check if box has changed since the last call of fcs_tune,
|
||||
// if it has, call fcs_tune
|
||||
if (box_has_changed())
|
||||
{
|
||||
setup_handle();
|
||||
|
||||
// print out parameters within handle TODO: should be done in C style
|
||||
/*
|
||||
if (rank == 0)
|
||||
{
|
||||
std::cout << " updated ScaFaCoS handle: " << std::endl;
|
||||
fcs_print_parameters(fcs);
|
||||
}
|
||||
*/
|
||||
|
||||
// call the tuning routine (needs to be redone, if critical
|
||||
// system parameters should change)
|
||||
result = fcs_tune(fcs,nlocal,x,q);
|
||||
if (!check_result(result, rank)) return;
|
||||
}
|
||||
|
||||
// check if arrays for potentials and field are still large enough
|
||||
// for the number of particles on process
|
||||
if (nlocal > local_array_size)
|
||||
{
|
||||
// allocate arrays larger in order to avoid fast
|
||||
// reallocation due to fluctuations
|
||||
local_array_size = (int)(1.25 * (double)nlocal);
|
||||
|
||||
// destroy old result arrays
|
||||
memory->destroy(pot);
|
||||
memory->destroy(field);
|
||||
|
||||
// allocate result arrays
|
||||
memory->create(pot,local_array_size,"scafacos:potential");
|
||||
memory->create(field,3*local_array_size,"scafacos:field");
|
||||
|
||||
}
|
||||
|
||||
// set result vectors to zero
|
||||
for ( int i = 0; i < nlocal; ++i)
|
||||
{
|
||||
pot[i] = 0.0;
|
||||
field[3*i] = field[3*i+1] = field[3*i+2] = 0.0;
|
||||
}
|
||||
|
||||
// compute Coulomb
|
||||
fcs_run(fcs, nlocal, x, q, field, pot);
|
||||
if(!check_result(result,rank)) return;
|
||||
|
||||
double E_coul_loc = 0.0;
|
||||
|
||||
for (int i = 0; i < atom->nlocal; ++i)
|
||||
{
|
||||
//std::cout << atom->f[i][0] << " " << field[3*i] << " " << q[i] << std::endl;
|
||||
atom->f[i][0] += field[3*i] * q[i];
|
||||
atom->f[i][1] += field[3*i+1] * q[i];
|
||||
atom->f[i][2] += field[3*i+2] * q[i];
|
||||
E_coul_loc += 0.5 * q[i] * pot[i];
|
||||
}
|
||||
|
||||
force->pair->eng_coul += E_coul_loc;
|
||||
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::min_post_force(int vflag)
|
||||
{
|
||||
post_force(vflag);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
integrate electronic degrees of freedom
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::final_integrate() {}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixScafacos::reset_dt()
|
||||
{
|
||||
//dtv = update->dt;
|
||||
//dtf = 0.5 * update->dt * force->ftm2v;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
DFTB energy from LATTE
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
double FixScafacos::compute_scalar()
|
||||
{
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
memory usage of local arrays
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
double FixScafacos::memory_usage()
|
||||
{
|
||||
double bytes = 0.0;
|
||||
bytes += local_array_size * sizeof(double);
|
||||
bytes += local_array_size * 3 * sizeof(double);
|
||||
return bytes;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
setup of ScaFaCoS handle with common parameters
|
||||
------------------------------------------------------------------------- */
|
||||
void FixScafacos::setup_handle()
|
||||
{
|
||||
// store periodicity
|
||||
periodicity[0] = domain->xperiodic;
|
||||
periodicity[1] = domain->yperiodic;
|
||||
periodicity[2] = domain->zperiodic;
|
||||
|
||||
// store offset of the system
|
||||
offset[0] = domain->boundary[0][0];
|
||||
offset[1] = domain->boundary[1][0];
|
||||
offset[2] = domain->boundary[2][0];
|
||||
|
||||
// calculate box vectors
|
||||
box_x[0] = domain->prd[0];
|
||||
box_x[1] = box_x[2] = 0.0;
|
||||
|
||||
box_y[1] = domain->prd[1];
|
||||
box_y[0] = box_y[2] = 0.0;
|
||||
|
||||
box_z[2] = domain->prd[2];
|
||||
box_z[1] = box_z[0] = 0.0;
|
||||
|
||||
total_particles = atom->natoms;
|
||||
|
||||
// TODO: for now disable short-range calculations within LAMMPS
|
||||
near_field_flag = 0;
|
||||
|
||||
// enter all parameters required to ScaFaCoS handle
|
||||
result = fcs_set_box_a(fcs, box_x);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
result = fcs_set_box_b(fcs, box_y);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
result = fcs_set_box_c(fcs, box_z);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
result = fcs_set_box_origin(fcs, offset);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
result = fcs_set_periodicity(fcs, periodicity);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
result = fcs_set_near_field_flag(fcs, near_field_flag);
|
||||
if (!check_result(result, rank)) return;
|
||||
|
||||
result = fcs_set_total_particles(fcs, atom->natoms);
|
||||
if (!check_result(result, rank)) return;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
check if box parameters changed, requiring a new call to fcs_tune
|
||||
------------------------------------------------------------------------- */
|
||||
bool FixScafacos::box_has_changed()
|
||||
{
|
||||
bool changed = false;
|
||||
|
||||
double n_periodicity[3];
|
||||
double n_offset[3];
|
||||
double n_box_x[3];
|
||||
double n_box_y[3];
|
||||
double n_box_z[3];
|
||||
|
||||
int n_total_particles;
|
||||
|
||||
// store periodicity
|
||||
n_periodicity[0] = domain->xperiodic;
|
||||
n_periodicity[1] = domain->yperiodic;
|
||||
n_periodicity[2] = domain->zperiodic;
|
||||
|
||||
// store offset of the system
|
||||
n_offset[0] = domain->boundary[0][0];
|
||||
n_offset[1] = domain->boundary[1][0];
|
||||
n_offset[2] = domain->boundary[2][0];
|
||||
|
||||
// calculate box vectors
|
||||
n_box_x[0] = domain->prd[0];
|
||||
n_box_x[1] = n_box_x[2] = 0.0;
|
||||
|
||||
n_box_y[1] = domain->prd[1];
|
||||
n_box_y[0] = n_box_y[2] = 0.0;
|
||||
|
||||
n_box_z[2] = domain->prd[2];
|
||||
n_box_z[1] = n_box_z[0] = 0.0;
|
||||
|
||||
n_total_particles = atom->natoms;
|
||||
|
||||
changed = changed ||
|
||||
( n_periodicity[0] != periodicity[0] ) ||
|
||||
( n_periodicity[1] != periodicity[1] ) ||
|
||||
( n_periodicity[2] != periodicity[2] ) ||
|
||||
( n_offset[0] != offset[0] ) ||
|
||||
( n_offset[1] != offset[1] ) ||
|
||||
( n_offset[2] != offset[2] ) ||
|
||||
( n_box_x[0] != box_x[0] ) ||
|
||||
( n_box_x[1] != box_x[1] ) ||
|
||||
( n_box_x[2] != box_x[2] ) ||
|
||||
( n_box_y[0] != box_y[0] ) ||
|
||||
( n_box_y[1] != box_y[1] ) ||
|
||||
( n_box_y[2] != box_y[2] ) ||
|
||||
( n_box_z[0] != box_z[0] ) ||
|
||||
( n_box_z[1] != box_z[1] ) ||
|
||||
( n_box_z[2] != box_z[2] ) ||
|
||||
( n_total_particles != total_particles );
|
||||
return changed;
|
||||
}
|
||||
|
||||
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
check of ScaFaCoS result
|
||||
------------------------------------------------------------------------- */
|
||||
bool FixScafacos::check_result(FCSResult result, int comm_rank)
|
||||
{
|
||||
if (result)
|
||||
{
|
||||
printf("ScaFaCoS Error: Caught error on task %d.\n", comm_rank);
|
||||
std::string err_msg;
|
||||
std::stringstream ss;
|
||||
|
||||
ss << fcs_result_get_function(result) << "\n"
|
||||
<< fcs_result_get_message(result) << "\n";
|
||||
err_msg = ss.str();
|
||||
|
||||
error -> all(FLERR, err_msg.c_str());
|
||||
fcs_result_destroy(result);
|
||||
}
|
||||
return true;
|
||||
}
|
||||
|
||||
@ -1,164 +0,0 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef FIX_CLASS
|
||||
|
||||
FixStyle(scafacos,FixScafacos)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_FIX_SCAFACOS_H
|
||||
#define LMP_FIX_SCAFACOS_H
|
||||
|
||||
#include "fix.h"
|
||||
#include "fcs.h"
|
||||
#include <string>
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class FixScafacos : public Fix {
|
||||
public:
|
||||
FixScafacos(class LAMMPS *, int, char **);
|
||||
virtual ~FixScafacos();
|
||||
int setmask();
|
||||
void init();
|
||||
void init_list(int, NeighList*);
|
||||
void setup(int);
|
||||
void min_setup(int);
|
||||
void setup_pre_reverse(int, int);
|
||||
void initial_integrate(int);
|
||||
void pre_reverse(int, int);
|
||||
void post_force(int);
|
||||
void min_post_force(int);
|
||||
void final_integrate();
|
||||
void reset_dt();
|
||||
double compute_scalar();
|
||||
double memory_usage();
|
||||
|
||||
protected:
|
||||
std::string method;
|
||||
|
||||
// MPI rank
|
||||
int rank;
|
||||
|
||||
// source arrays for positions and charges
|
||||
double *x, *q;
|
||||
// result arrays for potentials and field
|
||||
double *pot, *field;
|
||||
|
||||
// box vectors for each dimension
|
||||
fcs_float box_x[3], box_y[3], box_z[3];
|
||||
// offset of the box from the origin
|
||||
fcs_float offset[3];
|
||||
|
||||
// periodicity of the system
|
||||
fcs_int periodicity[3];
|
||||
|
||||
// ScaFaCoS handle
|
||||
FCS fcs;
|
||||
|
||||
// ScaFaCoS result variable
|
||||
FCSResult result;
|
||||
|
||||
// function to check results
|
||||
bool check_result(FCSResult, int);
|
||||
|
||||
// function to set up handle with common parameters
|
||||
void setup_handle();
|
||||
|
||||
// function to check if the box parameters changed, so that a new tuning step is required
|
||||
bool box_has_changed();
|
||||
|
||||
// store total number of particles (to check if tune needs to be called again)
|
||||
fcs_int total_particles;
|
||||
|
||||
// store number of local particles (to be able to readjust the size of result arrays, when needed)
|
||||
int local_array_size;
|
||||
|
||||
// should the near field calculations be computed by LAMMPS?
|
||||
fcs_int near_field_flag;
|
||||
|
||||
// type of accuracy chosen (if customized)
|
||||
fcs_int tolerance_type;
|
||||
|
||||
// value of tolerance
|
||||
fcs_float tolerance_value;
|
||||
|
||||
// is tolerance set?
|
||||
bool tolerance_set;
|
||||
|
||||
// check if fmm is chosen (ghost particles, since the implementation needs at least 1 particle on each process!)
|
||||
bool fmm_chosen;
|
||||
|
||||
// FMM: fmm particle array size
|
||||
int fmm_array_size;
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Must use units metal with fix latte command
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Fix latte currently runs only in serial
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: LAMMPS is linked against incompatible LATTE library
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Illegal ... command
|
||||
|
||||
Self-explanatory. Check the input script syntax and compare to the
|
||||
documentation for the command. You can use -echo screen as a
|
||||
command-line option when running LAMMPS to see the offending line.
|
||||
|
||||
E: Fix latte does not yet support a LAMMPS calculation of a Coulomb potential
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Could not find fix latte compute ID
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Fix latte compute ID does not compute pe/atom
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Fix latte requires 3d problem
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Fix latte cannot compute Coulomb potential
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Fix latte requires 3d simulation
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
W: Fix latte should come after all other integration fixes
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
E: Internal LATTE problem
|
||||
|
||||
UNDOCUMENTED
|
||||
|
||||
*/
|
||||
@ -38,40 +38,16 @@ using namespace LAMMPS_NS;
|
||||
|
||||
Scafacos::Scafacos(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
|
||||
{
|
||||
if (narg < 2) error->all(FLERR,"Illegal scafacos command");
|
||||
if (narg != 2) error->all(FLERR,"Illegal scafacos command");
|
||||
|
||||
int n = strlen(arg[0]) + 1;
|
||||
method = new char[n];
|
||||
strcpy(method,arg[0]);
|
||||
tolerance = force->numeric(FLERR,arg[1]);
|
||||
|
||||
// additional args
|
||||
// optional ScaFaCoS library setting defaults
|
||||
|
||||
int tflag = 0;
|
||||
|
||||
int iarg = 1;
|
||||
while (iarg < narg) {
|
||||
if (strcmp(arg[iarg],"tolerance") == 0) {
|
||||
if (iarg+3 > narg) error->all(FLERR,"Illegal scafacos command");
|
||||
tflag = 1;
|
||||
if (strcmp(arg[iarg+1],"energy") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_ENERGY;
|
||||
else if (strcmp(arg[iarg+1],"energy_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_ENERGY_REL;
|
||||
else if (strcmp(arg[iarg+1],"field") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD;
|
||||
else if (strcmp(arg[iarg+1],"field_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD_REL;
|
||||
else if (strcmp(arg[iarg+1],"potential") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_POTENTIAL;
|
||||
else if (strcmp(arg[iarg+1],"potential_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_POTENTIAL_REL;
|
||||
else error->all(FLERR,"Illegal scafacos command");
|
||||
tolerance = force->numeric(FLERR,arg[iarg+2]);
|
||||
iarg += 3;
|
||||
} else error->all(FLERR,"Illegal scafacos command");
|
||||
}
|
||||
|
||||
if (!tflag) error->all(FLERR,"Must set tolerance in scafacos command");
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD;
|
||||
|
||||
// initializations
|
||||
|
||||
@ -92,7 +68,7 @@ Scafacos::~Scafacos()
|
||||
memory->destroy(epot);
|
||||
memory->destroy(efield);
|
||||
|
||||
// NOTE: any clean-up call to ScaFaCoS needed?
|
||||
// RENE: any clean-up/shut-down call to ScaFaCoS needed?
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
@ -100,7 +76,6 @@ Scafacos::~Scafacos()
|
||||
void Scafacos::init()
|
||||
{
|
||||
// error checks
|
||||
// NOTE: allow triclinic at some point
|
||||
|
||||
if (domain->dimension == 2)
|
||||
error->all(FLERR,"Cannot use ScaFaCoS with 2d simulation");
|
||||
@ -108,11 +83,13 @@ void Scafacos::init()
|
||||
if (domain->triclinic)
|
||||
error->all(FLERR,"Cannot use ScaFaCoS with triclinic domain yet");
|
||||
|
||||
if (atom->natoms > INT_MAX && sizeof(fcs_int) != 8)
|
||||
error->all(FLERR,"Scafacos atom count exceeds 2B");
|
||||
|
||||
// one-time initialization of ScaFaCoS
|
||||
|
||||
scale = 1.0;
|
||||
qqrd2e = force->qqrd2e;
|
||||
//qsum_qsq();
|
||||
|
||||
if (!initialized) {
|
||||
result = fcs_init(&fcs,method,world);
|
||||
@ -147,13 +124,14 @@ void Scafacos::compute(int eflag, int vflag)
|
||||
double *q = atom->q;
|
||||
int nlocal = atom->nlocal;
|
||||
|
||||
// RENE: why is scale needed?
|
||||
|
||||
const double qscale = qqrd2e * scale;
|
||||
|
||||
// if simluation box has changed, call fcs_tune()
|
||||
|
||||
if (eflag || vflag) ev_setup(eflag,vflag);
|
||||
else
|
||||
eflag_atom = 0;
|
||||
else eflag_atom = 0;
|
||||
|
||||
// if simulation box has changed, call fcs_tune()
|
||||
|
||||
if (box_has_changed()) {
|
||||
setup_handle();
|
||||
@ -173,7 +151,7 @@ void Scafacos::compute(int eflag, int vflag)
|
||||
}
|
||||
|
||||
// initialize epot & efield
|
||||
// NOTE: is this necessary?
|
||||
// RENE: is this necessary? or does Scafacos just set them
|
||||
|
||||
for (int i = 0; i < nlocal; i++) {
|
||||
epot[i] = 0.0;
|
||||
@ -201,17 +179,47 @@ void Scafacos::compute(int eflag, int vflag)
|
||||
myeng += 0.5 * qone * epot[i];
|
||||
}
|
||||
|
||||
if (eflag_atom)
|
||||
{
|
||||
if (eflag_atom) {
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
{
|
||||
eatom[i] = qscale * epot[i];
|
||||
}
|
||||
}
|
||||
|
||||
MPI_Allreduce(&myeng,&energy,1,MPI_DOUBLE,MPI_SUM,world);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
int Scafacos::modify_param(int narg, char **arg)
|
||||
{
|
||||
// RENE: add any Scafacos options here you want to expose to LAMMPS
|
||||
// syntax: kspace_modify scafacos keyword value1 value2 ...
|
||||
// keyword = tolerance
|
||||
// value1 = energy, energy_rel, etc
|
||||
// everyone of these should have a default, so user doesn't need to set
|
||||
|
||||
if (strcmp(arg[0],"scafacos") != 0) return 0;
|
||||
|
||||
if (strcmp(arg[1],"tolerance") == 0) {
|
||||
if (narg < 2) error->all(FLERR,"Illegal kspace_modify command");
|
||||
if (strcmp(arg[2],"energy") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_ENERGY;
|
||||
else if (strcmp(arg[2],"energy_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_ENERGY_REL;
|
||||
else if (strcmp(arg[2],"field") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD;
|
||||
else if (strcmp(arg[2],"field_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_FIELD_REL;
|
||||
else if (strcmp(arg[2],"potential") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_POTENTIAL;
|
||||
else if (strcmp(arg[2],"potential_rel") == 0)
|
||||
tolerance_type = FCS_TOLERANCE_TYPE_POTENTIAL_REL;
|
||||
else error->all(FLERR,"Illegal kspace_modify command");
|
||||
return 3;
|
||||
}
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
memory usage of local arrays
|
||||
------------------------------------------------------------------------- */
|
||||
@ -231,15 +239,14 @@ double Scafacos::memory_usage()
|
||||
void Scafacos::setup_handle()
|
||||
{
|
||||
// store simulation box params
|
||||
// NOTE: this assumes orthogonal box
|
||||
// NOTE: total particles may not be a 4-byte int
|
||||
// NOTE: what does SCFCS mean by offset?
|
||||
// it's an integer flag in LAMMPS
|
||||
|
||||
old_periodicity[0] = domain->xperiodic;
|
||||
old_periodicity[1] = domain->yperiodic;
|
||||
old_periodicity[2] = domain->zperiodic;
|
||||
|
||||
// RENE: what does SCFCS mean by offset?
|
||||
// it's an integer flag in LAMMPS, but being stored in a float?
|
||||
|
||||
old_offset[0] = domain->boundary[0][0];
|
||||
old_offset[1] = domain->boundary[1][0];
|
||||
old_offset[2] = domain->boundary[2][0];
|
||||
@ -273,7 +280,9 @@ void Scafacos::setup_handle()
|
||||
result = fcs_set_total_particles(fcs,old_natoms);
|
||||
check_result(result);
|
||||
|
||||
// NOTE: for now disable short-range calculations within LAMMPS
|
||||
// RENE: disable short-range calculations within LAMMPS
|
||||
// not sure what this is doing
|
||||
// is this the correct thing to do for now?
|
||||
|
||||
int near_field_flag = 0;
|
||||
result = fcs_set_near_field_flag(fcs,near_field_flag);
|
||||
@ -306,7 +315,6 @@ bool Scafacos::box_has_changed()
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
check ScaFaCoS result for error condition
|
||||
NOTE: will all procs have same error?
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void Scafacos::check_result(FCSResult result)
|
||||
@ -320,6 +328,9 @@ void Scafacos::check_result(FCSResult result)
|
||||
std::string err_msg = ss.str();
|
||||
const char *str = err_msg.c_str();
|
||||
|
||||
// RENE: will all procs have same error?
|
||||
// if so, then should call error->all(FLERR,str)
|
||||
|
||||
error->one(FLERR,str);
|
||||
}
|
||||
|
||||
|
||||
@ -32,6 +32,7 @@ class Scafacos : public KSpace {
|
||||
void init();
|
||||
void setup();
|
||||
void compute(int, int);
|
||||
int modify_param(int, char **);
|
||||
double memory_usage();
|
||||
|
||||
private:
|
||||
|
||||
@ -580,7 +580,11 @@ void KSpace::modify_params(int narg, char **arg)
|
||||
else if (strcmp(arg[iarg+1],"no") == 0) auto_disp_flag = 0;
|
||||
else error->all(FLERR,"Illegal kspace_modify command");
|
||||
iarg += 2;
|
||||
} else error->all(FLERR,"Illegal kspace_modify command");
|
||||
} else {
|
||||
int n = modify_param(narg-iarg,&arg[iarg]);
|
||||
if (n == 0) error->all(FLERR,"Illegal kspace_modify command");
|
||||
iarg += n;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@ -126,6 +126,8 @@ class KSpace : protected Pointers {
|
||||
virtual int timing(int, double &, double &) {return 0;}
|
||||
virtual int timing_1d(int, double &) {return 0;}
|
||||
virtual int timing_3d(int, double &) {return 0;}
|
||||
|
||||
virtual int modify_param(int, char **) {return 0;}
|
||||
virtual double memory_usage() {return 0.0;}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
|
||||
Reference in New Issue
Block a user