Merge branch 'master' of https://github.com/lammps/lammps into kk_pair_fix_comm

This commit is contained in:
Stan Moore
2020-12-22 09:19:18 -07:00
66 changed files with 1645 additions and 699 deletions

View File

@ -84,13 +84,15 @@ information is available, then also a heuristic based on that bond length
is computed. It is used as communication cutoff, if there is no pair is computed. It is used as communication cutoff, if there is no pair
style present and no *comm_modify cutoff* command used. Otherwise a style present and no *comm_modify cutoff* command used. Otherwise a
warning is printed, if this bond based estimate is larger than the warning is printed, if this bond based estimate is larger than the
communication cutoff used. A communication cutoff used.
The *cutoff/multi* option is equivalent to *cutoff*\ , but applies to The *cutoff/multi* option is equivalent to *cutoff*\ , but applies to
communication mode *multi* instead. Since in this case the communication communication mode *multi* instead. Since in this case the communication
cutoffs are determined per atom type, a type specifier is needed and cutoffs are determined per atom type, a type specifier is needed and
cutoff for one or multiple types can be extended. Also ranges of types cutoff for one or multiple types can be extended. Also ranges of types
using the usual asterisk notation can be given. using the usual asterisk notation can be given. For granular pair styles,
the default cutoff is set to the sum of the current maximum atomic radii
for each type.
These are simulation scenarios in which it may be useful or even These are simulation scenarios in which it may be useful or even
necessary to set a ghost cutoff > neighbor cutoff: necessary to set a ghost cutoff > neighbor cutoff:

View File

@ -49,7 +49,9 @@ sometimes be faster. Either style should give the same answers.
The *multi* style is a modified binning algorithm that is useful for The *multi* style is a modified binning algorithm that is useful for
systems with a wide range of cutoff distances, e.g. due to different systems with a wide range of cutoff distances, e.g. due to different
size particles. For the *bin* style, the bin size is set to 1/2 of size particles. For granular pair styles, cutoffs are set to the
sum of the maximum atomic radii for each atom type.
For the *bin* style, the bin size is set to 1/2 of
the largest cutoff distance between any pair of atom types and a the largest cutoff distance between any pair of atom types and a
single set of bins is defined to search over for all atom types. This single set of bins is defined to search over for all atom types. This
can be inefficient if one pair of types has a very long cutoff, but can be inefficient if one pair of types has a very long cutoff, but
@ -57,8 +59,10 @@ other type pairs have a much shorter cutoff. For style *multi* the
bin size is set to 1/2 of the shortest cutoff distance and multiple bin size is set to 1/2 of the shortest cutoff distance and multiple
sets of bins are defined to search over for different atom types. sets of bins are defined to search over for different atom types.
This imposes some extra setup overhead, but the searches themselves This imposes some extra setup overhead, but the searches themselves
may be much faster for the short-cutoff cases. See the :doc:`comm_modify mode multi <comm_modify>` command for a communication option may be much faster for the short-cutoff cases.
that may also be beneficial for simulations of this kind. See the :doc:`comm_modify mode multi <comm_modify>` command for a
communication option that may also be beneficial for simulations of
this kind.
The :doc:`neigh_modify <neigh_modify>` command has additional options The :doc:`neigh_modify <neigh_modify>` command has additional options
that control how often neighbor lists are built and which pairs are that control how often neighbor lists are built and which pairs are

View File

@ -77,7 +77,7 @@ class lammps(object):
modpath = dirname(abspath(getsourcefile(lambda:0))) modpath = dirname(abspath(getsourcefile(lambda:0)))
# for windows installers the shared library is in a different folder # for windows installers the shared library is in a different folder
winpath = abspath(os.path.join(modpath,'..','bin')) winpath = abspath(os.path.join(modpath,'..','..','bin'))
self.lib = None self.lib = None
self.lmp = None self.lmp = None

View File

@ -16,10 +16,11 @@
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
#include "pair_lj_cubic_gpu.h" #include "pair_lj_cubic_gpu.h"
#include <cmath> #include <cmath>
#include <cstdio> #include <cstdio>
#include <cstring> #include <cstring>
#include "atom.h" #include "atom.h"
#include "atom_vec.h" #include "atom_vec.h"
#include "comm.h" #include "comm.h"
@ -36,6 +37,8 @@
#include "gpu_extra.h" #include "gpu_extra.h"
#include "suffix.h" #include "suffix.h"
#include "pair_lj_cubic_const.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace PairLJCubicConstants; using namespace PairLJCubicConstants;

View File

@ -21,6 +21,7 @@
#include "atom_masks.h" #include "atom_masks.h"
#include "kokkos.h" #include "kokkos.h"
#include "math_const.h" #include "math_const.h"
#include "math_special.h"
#include "memory_kokkos.h" #include "memory_kokkos.h"
#include "neigh_list.h" #include "neigh_list.h"
#include "neigh_request.h" #include "neigh_request.h"
@ -32,6 +33,7 @@
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace MathConst; using namespace MathConst;
using namespace MathSpecial;
using namespace std; using namespace std;
#ifdef DBL_EPSILON #ifdef DBL_EPSILON

View File

@ -30,6 +30,8 @@
#include "atom_masks.h" #include "atom_masks.h"
#include "kokkos.h" #include "kokkos.h"
#include "pair_zbl_const.h"
// From J.F. Zeigler, J. P. Biersack and U. Littmark, // From J.F. Zeigler, J. P. Biersack and U. Littmark,
// "The Stopping and Range of Ions in Matter" volume 1, Pergamon, 1985. // "The Stopping and Range of Ions in Matter" volume 1, Pergamon, 1985.

View File

@ -38,7 +38,7 @@ class PairComb : public Pair {
virtual double yasu_char(double *, int &); virtual double yasu_char(double *, int &);
double enegtot; double enegtot;
static const int NPARAMS_PER_LINE = 49; static constexpr int NPARAMS_PER_LINE = 49;
protected: protected:
struct Param { struct Param {

View File

@ -37,7 +37,7 @@ class PairComb3 : public Pair {
virtual double combqeq(double *, int &); virtual double combqeq(double *, int &);
double enegtot; double enegtot;
static const int NPARAMS_PER_LINE = 74; static constexpr int NPARAMS_PER_LINE = 74;
protected: protected:
// general potential parameters // general potential parameters

View File

@ -34,7 +34,7 @@ class PairGW : public Pair {
void init_style(); void init_style();
double init_one(int, int); double init_one(int, int);
static const int NPARAMS_PER_LINE = 17; static constexpr int NPARAMS_PER_LINE = 17;
protected: protected:
struct Param { struct Param {

View File

@ -29,7 +29,7 @@ class PairGWZBL : public PairGW {
PairGWZBL(class LAMMPS *); PairGWZBL(class LAMMPS *);
~PairGWZBL() {} ~PairGWZBL() {}
static const int NPARAMS_PER_LINE = 21; static constexpr int NPARAMS_PER_LINE = 21;
private: private:
double global_a_0; // Bohr radius for Coulomb repulsion double global_a_0; // Bohr radius for Coulomb repulsion

View File

@ -34,7 +34,7 @@ class PairNb3bHarmonic : public Pair {
double init_one(int, int); double init_one(int, int);
void init_style(); void init_style();
static const int NPARAMS_PER_LINE = 6; static constexpr int NPARAMS_PER_LINE = 6;
protected: protected:
struct Param { struct Param {

View File

@ -34,7 +34,7 @@ class PairSW : public Pair {
virtual double init_one(int, int); virtual double init_one(int, int);
virtual void init_style(); virtual void init_style();
static const int NPARAMS_PER_LINE = 14; static constexpr int NPARAMS_PER_LINE = 14;
struct Param { struct Param {
double epsilon,sigma; double epsilon,sigma;

View File

@ -34,7 +34,7 @@ class PairTersoff : public Pair {
virtual void init_style(); virtual void init_style();
double init_one(int, int); double init_one(int, int);
static const int NPARAMS_PER_LINE = 17; static constexpr int NPARAMS_PER_LINE = 17;
protected: protected:

View File

@ -30,7 +30,7 @@ class PairTersoffMOD : public PairTersoff {
PairTersoffMOD(class LAMMPS *); PairTersoffMOD(class LAMMPS *);
~PairTersoffMOD() {} ~PairTersoffMOD() {}
static const int NPARAMS_PER_LINE = 20; static constexpr int NPARAMS_PER_LINE = 20;
protected: protected:
virtual void read_file(char *); virtual void read_file(char *);

View File

@ -29,7 +29,7 @@ class PairTersoffMODC : public PairTersoffMOD {
PairTersoffMODC(class LAMMPS *lmp) : PairTersoffMOD(lmp) {}; PairTersoffMODC(class LAMMPS *lmp) : PairTersoffMOD(lmp) {};
~PairTersoffMODC() {} ~PairTersoffMODC() {}
static const int NPARAMS_PER_LINE = 21; static constexpr int NPARAMS_PER_LINE = 21;
protected: protected:
void read_file(char *); void read_file(char *);

View File

@ -29,7 +29,7 @@ class PairTersoffZBL : public PairTersoff {
PairTersoffZBL(class LAMMPS *); PairTersoffZBL(class LAMMPS *);
~PairTersoffZBL() {} ~PairTersoffZBL() {}
static const int NPARAMS_PER_LINE = 21; static constexpr int NPARAMS_PER_LINE = 21;
protected: protected:
double global_a_0; // Bohr radius for Coulomb repulsion double global_a_0; // Bohr radius for Coulomb repulsion

View File

@ -34,7 +34,7 @@ class PairVashishta : public Pair {
double init_one(int, int); double init_one(int, int);
void init_style(); void init_style();
static const int NPARAMS_PER_LINE = 17; static constexpr int NPARAMS_PER_LINE = 17;
struct Param { struct Param {
double bigb,gamma,r0,bigc,costheta; double bigb,gamma,r0,bigc,costheta;

View File

@ -33,21 +33,21 @@ namespace LAMMPS_NS {
TORQUE = 1<<8 TORQUE = 1<<8
}; };
static const double TOLERANCE = 1.0e-6; static constexpr double TOLERANCE = 1.0e-6;
static const double EPSILON = 1.0e-7; static constexpr double EPSILON = 1.0e-7;
static const double BIG = 1.0e20; static constexpr double BIG = 1.0e20;
// moment of inertia prefactor for sphere // moment of inertia prefactor for sphere
static const double SINERTIA = 0.4; static constexpr double SINERTIA = 0.4;
// moment of inertia prefactor for ellipsoid // moment of inertia prefactor for ellipsoid
static const double EINERTIA = 0.2; static constexpr double EINERTIA = 0.2;
// moment of inertia prefactor for line segment // moment of inertia prefactor for line segment
static const double LINERTIA = 1.0/12.0; static constexpr double LINERTIA = 1.0/12.0;
static const int MAXLINE = 1024; static constexpr int MAXLINE = 1024;
static const int CHUNK = 1024; static constexpr int CHUNK = 1024;
static const int DELTA_BODY = 10000; static constexpr int DELTA_BODY = 10000;
static const int ATTRIBUTE_PERBODY = 20; static constexpr int ATTRIBUTE_PERBODY = 20;
} }
} }

View File

@ -18,6 +18,7 @@
#include "sna.h" #include "sna.h"
#include <cmath> #include <cmath>
#include "math_const.h" #include "math_const.h"
#include "math_special.h"
#include "memory.h" #include "memory.h"
#include "error.h" #include "error.h"
#include "comm.h" #include "comm.h"
@ -25,6 +26,7 @@
using namespace std; using namespace std;
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace MathConst; using namespace MathConst;
using namespace MathSpecial;
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -1363,196 +1365,6 @@ void SNA::destroy_twojmax_arrays()
} }
/* ----------------------------------------------------------------------
factorial n, wrapper for precomputed table
------------------------------------------------------------------------- */
double SNA::factorial(int n)
{
if (n < 0 || n > nmaxfactorial) {
char str[128];
sprintf(str, "Invalid argument to factorial %d", n);
error->all(FLERR, str);
}
return nfac_table[n];
}
/* ----------------------------------------------------------------------
factorial n table, size SNA::nmaxfactorial+1
------------------------------------------------------------------------- */
const double SNA::nfac_table[] = {
1,
1,
2,
6,
24,
120,
720,
5040,
40320,
362880,
3628800,
39916800,
479001600,
6227020800,
87178291200,
1307674368000,
20922789888000,
355687428096000,
6.402373705728e+15,
1.21645100408832e+17,
2.43290200817664e+18,
5.10909421717094e+19,
1.12400072777761e+21,
2.5852016738885e+22,
6.20448401733239e+23,
1.5511210043331e+25,
4.03291461126606e+26,
1.08888694504184e+28,
3.04888344611714e+29,
8.8417619937397e+30,
2.65252859812191e+32,
8.22283865417792e+33,
2.63130836933694e+35,
8.68331761881189e+36,
2.95232799039604e+38,
1.03331479663861e+40,
3.71993326789901e+41,
1.37637530912263e+43,
5.23022617466601e+44,
2.03978820811974e+46,
8.15915283247898e+47,
3.34525266131638e+49,
1.40500611775288e+51,
6.04152630633738e+52,
2.65827157478845e+54,
1.1962222086548e+56,
5.50262215981209e+57,
2.58623241511168e+59,
1.24139155925361e+61,
6.08281864034268e+62,
3.04140932017134e+64,
1.55111875328738e+66,
8.06581751709439e+67,
4.27488328406003e+69,
2.30843697339241e+71,
1.26964033536583e+73,
7.10998587804863e+74,
4.05269195048772e+76,
2.35056133128288e+78,
1.3868311854569e+80,
8.32098711274139e+81,
5.07580213877225e+83,
3.14699732603879e+85,
1.98260831540444e+87,
1.26886932185884e+89,
8.24765059208247e+90,
5.44344939077443e+92,
3.64711109181887e+94,
2.48003554243683e+96,
1.71122452428141e+98,
1.19785716699699e+100,
8.50478588567862e+101,
6.12344583768861e+103,
4.47011546151268e+105,
3.30788544151939e+107,
2.48091408113954e+109,
1.88549470166605e+111,
1.45183092028286e+113,
1.13242811782063e+115,
8.94618213078297e+116,
7.15694570462638e+118,
5.79712602074737e+120,
4.75364333701284e+122,
3.94552396972066e+124,
3.31424013456535e+126,
2.81710411438055e+128,
2.42270953836727e+130,
2.10775729837953e+132,
1.85482642257398e+134,
1.65079551609085e+136,
1.48571596448176e+138,
1.3520015276784e+140,
1.24384140546413e+142,
1.15677250708164e+144,
1.08736615665674e+146,
1.03299784882391e+148,
9.91677934870949e+149,
9.61927596824821e+151,
9.42689044888324e+153,
9.33262154439441e+155,
9.33262154439441e+157,
9.42594775983835e+159,
9.61446671503512e+161,
9.90290071648618e+163,
1.02990167451456e+166,
1.08139675824029e+168,
1.14628056373471e+170,
1.22652020319614e+172,
1.32464181945183e+174,
1.44385958320249e+176,
1.58824554152274e+178,
1.76295255109024e+180,
1.97450685722107e+182,
2.23119274865981e+184,
2.54355973347219e+186,
2.92509369349301e+188,
3.3931086844519e+190,
3.96993716080872e+192,
4.68452584975429e+194,
5.5745857612076e+196,
6.68950291344912e+198,
8.09429852527344e+200,
9.8750442008336e+202,
1.21463043670253e+205,
1.50614174151114e+207,
1.88267717688893e+209,
2.37217324288005e+211,
3.01266001845766e+213,
3.8562048236258e+215,
4.97450422247729e+217,
6.46685548922047e+219,
8.47158069087882e+221,
1.118248651196e+224,
1.48727070609069e+226,
1.99294274616152e+228,
2.69047270731805e+230,
3.65904288195255e+232,
5.01288874827499e+234,
6.91778647261949e+236,
9.61572319694109e+238,
1.34620124757175e+241,
1.89814375907617e+243,
2.69536413788816e+245,
3.85437071718007e+247,
5.5502938327393e+249,
8.04792605747199e+251,
1.17499720439091e+254,
1.72724589045464e+256,
2.55632391787286e+258,
3.80892263763057e+260,
5.71338395644585e+262,
8.62720977423323e+264,
1.31133588568345e+267,
2.00634390509568e+269,
3.08976961384735e+271,
4.78914290146339e+273,
7.47106292628289e+275,
1.17295687942641e+278,
1.85327186949373e+280,
2.94670227249504e+282,
4.71472363599206e+284,
7.59070505394721e+286,
1.22969421873945e+289,
2.0044015765453e+291,
3.28721858553429e+293,
5.42391066613159e+295,
9.00369170577843e+297,
1.503616514865e+300, // nmaxfactorial = 167
};
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
the function delta given by VMK Eq. 8.2(1) the function delta given by VMK Eq. 8.2(1)
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */

View File

@ -100,10 +100,6 @@ private:
double** dulist_r, ** dulist_i; double** dulist_r, ** dulist_i;
int elem_duarray; // element of j in derivative int elem_duarray; // element of j in derivative
static const int nmaxfactorial = 167;
static const double nfac_table[];
double factorial(int);
void create_twojmax_arrays(); void create_twojmax_arrays();
void destroy_twojmax_arrays(); void destroy_twojmax_arrays();
void init_clebsch_gordan(); void init_clebsch_gordan();

View File

@ -195,6 +195,8 @@ void DihedralQuadratic::compute(int eflag, int vflag)
siinv = 1.0/si; siinv = 1.0/si;
double dphi = phi-phi0[type]; double dphi = phi-phi0[type];
if (dphi > MY_PI) dphi -= 2*MY_PI;
else if (dphi < -MY_PI) dphi += 2*MY_PI;
p = k[type]*dphi; p = k[type]*dphi;
pd = - 2.0 * p * siinv; pd = - 2.0 * p * siinv;
p = p * dphi; p = p * dphi;

View File

@ -46,7 +46,7 @@ using namespace LAMMPS_NS;
// max number of interaction per atom for f(Z) environment potential // max number of interaction per atom for f(Z) environment potential
#define leadDimInteractionList 64 static constexpr int leadDimInteractionList = 64;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

View File

@ -32,13 +32,11 @@
#include "error.h" #include "error.h"
#include "citeme.h" #include "citeme.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
#define MAXLINE 1024 #define MAXLINE 1024
#define DELTA 4 #define DELTA 4
static const char cite_pair_edip[] = static const char cite_pair_edip[] =
"@article{cjiang2012\n" "@article{cjiang2012\n"
" author = {Jian, Chao and Morgan, Dane, and Szlufarska, Izabella},\n" " author = {Jian, Chao and Morgan, Dane, and Szlufarska, Izabella},\n"
@ -56,7 +54,9 @@ static const char cite_pair_edip[] =
" year = {2010},\n" " year = {2010},\n"
"}\n\n"; "}\n\n";
// max number of interaction per atom for f(Z) environment potential
static constexpr int leadDimInteractionList = 64;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
@ -94,7 +94,6 @@ PairEDIPMulti::~PairEDIPMulti()
memory->destroy(cutsq); memory->destroy(cutsq);
delete [] map; delete [] map;
//XXX deallocateGrids();
deallocatePreLoops(); deallocatePreLoops();
} }
} }

View File

@ -62,10 +62,6 @@ class PairEDIPMulti : public Pair {
int maxparam; // max # of parameter sets int maxparam; // max # of parameter sets
Param *params; // parameter set for an I-J-K interaction Param *params; // parameter set for an I-J-K interaction
// max number of interaction per atom for f(Z) environment potential
static const int leadDimInteractionList = 64;
void allocate(); void allocate();
void allocatePreLoops(void); void allocatePreLoops(void);
void deallocatePreLoops(void); void deallocatePreLoops(void);

View File

@ -43,7 +43,7 @@ class PairTersoffTable : public Pair {
void init_style(); void init_style();
double init_one(int, int); double init_one(int, int);
static const int NPARAMS_PER_LINE = 17; static constexpr int NPARAMS_PER_LINE = 17;
protected: protected:
struct Param { struct Param {

View File

@ -23,11 +23,13 @@
#include "neighbor.h" #include "neighbor.h"
#include "force.h" #include "force.h"
#include "update.h" #include "update.h"
#include "math_const.h"
#include "error.h" #include "error.h"
#include "suffix.h" #include "suffix.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace MathConst;
#define TOLERANCE 0.05 #define TOLERANCE 0.05
#define SMALL 0.001 #define SMALL 0.001
@ -218,6 +220,8 @@ void DihedralQuadraticOMP::eval(int nfrom, int nto, ThrData * const thr)
siinv = 1.0/si; siinv = 1.0/si;
double dphi = phi-phi0[type]; double dphi = phi-phi0[type];
if (dphi > MY_PI) dphi -= 2*MY_PI;
else if (dphi < -MY_PI) dphi += 2*MY_PI;
p = k[type]*dphi; p = k[type]*dphi;
pd = - 2.0 * p * siinv; pd = - 2.0 * p * siinv;
p = p * dphi; p = p * dphi;

View File

@ -0,0 +1,126 @@
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
#include "omp_compat.h"
#include "npair_half_size_multi_newtoff_omp.h"
#include "npair_omp.h"
#include "neigh_list.h"
#include "atom.h"
#include "atom_vec.h"
#include "my_page.h"
#include "error.h"
using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
NPairHalfSizeMultiNewtoffOmp::NPairHalfSizeMultiNewtoffOmp(LAMMPS *lmp) :
NPair(lmp) {}
/* ----------------------------------------------------------------------
size particles
binned neighbor list construction with partial Newton's 3rd law
each owned atom i checks own bin and other bins in stencil
multi-type stencil is itype dependent and is distance checked
pair stored once if i,j are both owned and i < j
pair stored by me if j is ghost (also stored by proc owning j)
------------------------------------------------------------------------- */
void NPairHalfSizeMultiNewtoffOmp::build(NeighList *list)
{
const int nlocal = (includegroup) ? atom->nfirst : atom->nlocal;
const int history = list->history;
const int mask_history = 3 << SBBITS;
NPAIR_OMP_INIT;
#if defined(_OPENMP)
#pragma omp parallel LMP_DEFAULT_NONE LMP_SHARED(list)
#endif
NPAIR_OMP_SETUP(nlocal);
int i,j,k,m,n,itype,jtype,ibin,ns;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
double radi,radsum,cutdistsq;
int *neighptr,*s;
double *cutsq,*distsq;
double **x = atom->x;
double *radius = atom->radius;
int *type = atom->type;
int *mask = atom->mask;
tagint *molecule = atom->molecule;
int *ilist = list->ilist;
int *numneigh = list->numneigh;
int **firstneigh = list->firstneigh;
// each thread has its own page allocator
MyPage<int> &ipage = list->ipage[tid];
ipage.reset();
for (i = ifrom; i < ito; i++) {
n = 0;
neighptr = ipage.vget();
itype = type[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
radi = radius[i];
// loop over all atoms in other bins in stencil including self
// only store pair if i < j
// skip if i,j neighbor cutoff is less than bin distance
// stores own/own pairs only once
// stores own/ghost pairs on both procs
ibin = atom2bin[i];
s = stencil_multi[itype];
distsq = distsq_multi[itype];
cutsq = cutneighsq[itype];
ns = nstencil_multi[itype];
for (k = 0; k < ns; k++) {
for (j = binhead[ibin+s[k]]; j >= 0; j = bins[j]) {
if (j <= i) continue;
jtype = type[j];
if (cutsq[jtype] < distsq[k]) continue;
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
}
ilist[i] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
ipage.vgot(n);
if (ipage.status())
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NPAIR_OMP_CLOSE;
list->inum = nlocal;
}

View File

@ -0,0 +1,44 @@
/* -*- 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 NPAIR_CLASS
NPairStyle(half/size/multi/newtoff/omp,
NPairHalfSizeMultiNewtoffOmp,
NP_HALF | NP_SIZE | NP_MULTI | NP_NEWTOFF | NP_OMP |
NP_ORTHO | NP_TRI)
#else
#ifndef LMP_NPAIR_HALF_SIZE_MULTI_NEWTOFF_OMP_H
#define LMP_NPAIR_HALF_SIZE_MULTI_NEWTOFF_OMP_H
#include "npair.h"
namespace LAMMPS_NS {
class NPairHalfSizeMultiNewtoffOmp : public NPair {
public:
NPairHalfSizeMultiNewtoffOmp(class LAMMPS *);
~NPairHalfSizeMultiNewtoffOmp() {}
void build(class NeighList *);
};
}
#endif
#endif
/* ERROR/WARNING messages:
*/

View File

@ -0,0 +1,151 @@
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
#include "omp_compat.h"
#include "npair_half_size_multi_newton_omp.h"
#include "npair_omp.h"
#include "neigh_list.h"
#include "atom.h"
#include "atom_vec.h"
#include "my_page.h"
#include "error.h"
using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
NPairHalfSizeMultiNewtonOmp::NPairHalfSizeMultiNewtonOmp(LAMMPS *lmp) :
NPair(lmp) {}
/* ----------------------------------------------------------------------
size particles
binned neighbor list construction with full Newton's 3rd law
each owned atom i checks its own bin and other bins in Newton stencil
multi-type stencil is itype dependent and is distance checked
every pair stored exactly once by some processor
------------------------------------------------------------------------- */
void NPairHalfSizeMultiNewtonOmp::build(NeighList *list)
{
const int nlocal = (includegroup) ? atom->nfirst : atom->nlocal;
const int history = list->history;
const int mask_history = 3 << SBBITS;
NPAIR_OMP_INIT;
#if defined(_OPENMP)
#pragma omp parallel LMP_DEFAULT_NONE LMP_SHARED(list)
#endif
NPAIR_OMP_SETUP(nlocal);
int i,j,k,m,n,itype,jtype,ibin,ns;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
double radi,radsum,cutdistsq;
int *neighptr,*s;
double *cutsq,*distsq;
double **x = atom->x;
double *radius = atom->radius;
int *type = atom->type;
int *mask = atom->mask;
tagint *molecule = atom->molecule;
int *ilist = list->ilist;
int *numneigh = list->numneigh;
int **firstneigh = list->firstneigh;
// each thread has its own page allocator
MyPage<int> &ipage = list->ipage[tid];
ipage.reset();
for (i = ifrom; i < ito; i++) {
n = 0;
neighptr = ipage.vget();
itype = type[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
radi = radius[i];
// loop over rest of atoms in i's bin, ghosts are at end of linked list
// if j is owned atom, store it, since j is beyond i in linked list
// if j is ghost, only store if j coords are "above and to the right" of i
for (j = bins[i]; j >= 0; j = bins[j]) {
if (j >= nlocal) {
if (x[j][2] < ztmp) continue;
if (x[j][2] == ztmp) {
if (x[j][1] < ytmp) continue;
if (x[j][1] == ytmp && x[j][0] < xtmp) continue;
}
}
jtype = type[j];
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
// loop over all atoms in other bins in stencil, store every pair
// skip if i,j neighbor cutoff is less than bin distance
ibin = atom2bin[i];
s = stencil_multi[itype];
distsq = distsq_multi[itype];
cutsq = cutneighsq[itype];
ns = nstencil_multi[itype];
for (k = 0; k < ns; k++) {
for (j = binhead[ibin+s[k]]; j >= 0; j = bins[j]) {
jtype = type[j];
if (cutsq[jtype] < distsq[k]) continue;
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
}
ilist[i] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
ipage.vgot(n);
if (ipage.status())
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NPAIR_OMP_CLOSE;
list->inum = nlocal;
}

View File

@ -0,0 +1,43 @@
/* -*- 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 NPAIR_CLASS
NPairStyle(half/size/multi/newton/omp,
NPairHalfSizeMultiNewtonOmp,
NP_HALF | NP_SIZE | NP_MULTI | NP_NEWTON | NP_OMP | NP_ORTHO)
#else
#ifndef LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_OMP_H
#define LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_OMP_H
#include "npair.h"
namespace LAMMPS_NS {
class NPairHalfSizeMultiNewtonOmp : public NPair {
public:
NPairHalfSizeMultiNewtonOmp(class LAMMPS *);
~NPairHalfSizeMultiNewtonOmp() {}
void build(class NeighList *);
};
}
#endif
#endif
/* ERROR/WARNING messages:
*/

View File

@ -0,0 +1,134 @@
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
#include "omp_compat.h"
#include "npair_half_size_multi_newton_tri_omp.h"
#include "npair_omp.h"
#include "neigh_list.h"
#include "atom.h"
#include "atom_vec.h"
#include "my_page.h"
#include "error.h"
using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
NPairHalfSizeMultiNewtonTriOmp::NPairHalfSizeMultiNewtonTriOmp(LAMMPS *lmp) :
NPair(lmp) {}
/* ----------------------------------------------------------------------
size particles
binned neighbor list construction with Newton's 3rd law for triclinic
each owned atom i checks its own bin and other bins in triclinic stencil
multi-type stencil is itype dependent and is distance checked
every pair stored exactly once by some processor
------------------------------------------------------------------------- */
void NPairHalfSizeMultiNewtonTriOmp::build(NeighList *list)
{
const int nlocal = (includegroup) ? atom->nfirst : atom->nlocal;
const int history = list->history;
const int mask_history = 3 << SBBITS;
NPAIR_OMP_INIT;
#if defined(_OPENMP)
#pragma omp parallel LMP_DEFAULT_NONE LMP_SHARED(list)
#endif
NPAIR_OMP_SETUP(nlocal);
int i,j,k,m,n,itype,jtype,ibin,ns;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
double radi,radsum,cutdistsq;
int *neighptr,*s;
double *cutsq,*distsq;
double **x = atom->x;
double *radius = atom->radius;
int *type = atom->type;
int *mask = atom->mask;
tagint *molecule = atom->molecule;
int *ilist = list->ilist;
int *numneigh = list->numneigh;
int **firstneigh = list->firstneigh;
// each thread has its own page allocator
MyPage<int> &ipage = list->ipage[tid];
ipage.reset();
for (i = ifrom; i < ito; i++) {
n = 0;
neighptr = ipage.vget();
itype = type[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
radi = radius[i];
// loop over all atoms in bins, including self, in stencil
// skip if i,j neighbor cutoff is less than bin distance
// bins below self are excluded from stencil
// pairs for atoms j "below" i are excluded
// below = lower z or (equal z and lower y) or (equal zy and lower x)
// (equal zyx and j <= i)
// latter excludes self-self interaction but allows superposed atoms
ibin = atom2bin[i];
s = stencil_multi[itype];
distsq = distsq_multi[itype];
cutsq = cutneighsq[itype];
ns = nstencil_multi[itype];
for (k = 0; k < ns; k++) {
for (j = binhead[ibin+s[k]]; j >= 0; j = bins[j]) {
jtype = type[j];
if (cutsq[jtype] < distsq[k]) continue;
if (x[j][2] < ztmp) continue;
if (x[j][2] == ztmp) {
if (x[j][1] < ytmp) continue;
if (x[j][1] == ytmp) {
if (x[j][0] < xtmp) continue;
if (x[j][0] == xtmp && j <= i) continue;
}
}
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
}
ilist[i] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
ipage.vgot(n);
if (ipage.status())
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NPAIR_OMP_CLOSE;
list->inum = nlocal;
}

View File

@ -0,0 +1,43 @@
/* -*- 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 NPAIR_CLASS
NPairStyle(half/size/multi/newton/tri/omp,
NPairHalfSizeMultiNewtonTriOmp,
NP_HALF | NP_SIZE | NP_MULTI | NP_NEWTON | NP_TRI | NP_OMP)
#else
#ifndef LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_TRI_OMP_H
#define LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_TRI_OMP_H
#include "npair.h"
namespace LAMMPS_NS {
class NPairHalfSizeMultiNewtonTriOmp : public NPair {
public:
NPairHalfSizeMultiNewtonTriOmp(class LAMMPS *);
~NPairHalfSizeMultiNewtonTriOmp() {}
void build(class NeighList *);
};
}
#endif
#endif
/* ERROR/WARNING messages:
*/

View File

@ -28,7 +28,7 @@ using namespace LAMMPS_NS;
// max number of interaction per atom for f(Z) environment potential // max number of interaction per atom for f(Z) environment potential
#define leadDimInteractionList 64 static constexpr int leadDimInteractionList = 64;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

View File

@ -23,6 +23,8 @@
#include <cmath> #include <cmath>
#include "omp_compat.h" #include "omp_compat.h"
#include "pair_lj_cubic_const.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace PairLJCubicConstants; using namespace PairLJCubicConstants;

View File

@ -24,7 +24,6 @@
#include "omp_compat.h" #include "omp_compat.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace PairZBLConstants;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

View File

@ -72,17 +72,6 @@ using namespace MathSpecial;
#define PGDELTA 1 #define PGDELTA 1
#define MAXNEIGH 24 #define MAXNEIGH 24
/* ------------------------------------------------------------------------------------
Calculates the factorial of an integer n via recursion.
------------------------------------------------------------------------------------ */
static double factorial(int n)
{
if (n <= 1) return 1.0;
else return static_cast<double>(n)*factorial(n-1);
}
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
PairSMTBQ::PairSMTBQ(LAMMPS *lmp) : Pair(lmp) PairSMTBQ::PairSMTBQ(LAMMPS *lmp) : Pair(lmp)

View File

@ -23,6 +23,7 @@
#include "error.h" #include "error.h"
#include "force.h" #include "force.h"
#include "math_const.h" #include "math_const.h"
#include "math_special.h"
#include "memory.h" #include "memory.h"
#include "modify.h" #include "modify.h"
#include "neigh_list.h" #include "neigh_list.h"
@ -36,6 +37,8 @@
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace MathConst; using namespace MathConst;
using namespace MathSpecial;
#ifdef DBL_EPSILON #ifdef DBL_EPSILON
#define MY_EPSILON (10.0*DBL_EPSILON) #define MY_EPSILON (10.0*DBL_EPSILON)
@ -708,191 +711,3 @@ void ComputeOrientOrderAtom::init_clebsch_gordan()
} }
} }
} }
/* ----------------------------------------------------------------------
factorial n, wrapper for precomputed table
------------------------------------------------------------------------- */
double ComputeOrientOrderAtom::factorial(int n)
{
if (n < 0 || n > nmaxfactorial)
error->all(FLERR,fmt::format("Invalid argument to factorial {}", n));
return nfac_table[n];
}
/* ----------------------------------------------------------------------
factorial n table, size SNA::nmaxfactorial+1
------------------------------------------------------------------------- */
const double ComputeOrientOrderAtom::nfac_table[] = {
1,
1,
2,
6,
24,
120,
720,
5040,
40320,
362880,
3628800,
39916800,
479001600,
6227020800,
87178291200,
1307674368000,
20922789888000,
355687428096000,
6.402373705728e+15,
1.21645100408832e+17,
2.43290200817664e+18,
5.10909421717094e+19,
1.12400072777761e+21,
2.5852016738885e+22,
6.20448401733239e+23,
1.5511210043331e+25,
4.03291461126606e+26,
1.08888694504184e+28,
3.04888344611714e+29,
8.8417619937397e+30,
2.65252859812191e+32,
8.22283865417792e+33,
2.63130836933694e+35,
8.68331761881189e+36,
2.95232799039604e+38,
1.03331479663861e+40,
3.71993326789901e+41,
1.37637530912263e+43,
5.23022617466601e+44,
2.03978820811974e+46,
8.15915283247898e+47,
3.34525266131638e+49,
1.40500611775288e+51,
6.04152630633738e+52,
2.65827157478845e+54,
1.1962222086548e+56,
5.50262215981209e+57,
2.58623241511168e+59,
1.24139155925361e+61,
6.08281864034268e+62,
3.04140932017134e+64,
1.55111875328738e+66,
8.06581751709439e+67,
4.27488328406003e+69,
2.30843697339241e+71,
1.26964033536583e+73,
7.10998587804863e+74,
4.05269195048772e+76,
2.35056133128288e+78,
1.3868311854569e+80,
8.32098711274139e+81,
5.07580213877225e+83,
3.14699732603879e+85,
1.98260831540444e+87,
1.26886932185884e+89,
8.24765059208247e+90,
5.44344939077443e+92,
3.64711109181887e+94,
2.48003554243683e+96,
1.71122452428141e+98,
1.19785716699699e+100,
8.50478588567862e+101,
6.12344583768861e+103,
4.47011546151268e+105,
3.30788544151939e+107,
2.48091408113954e+109,
1.88549470166605e+111,
1.45183092028286e+113,
1.13242811782063e+115,
8.94618213078297e+116,
7.15694570462638e+118,
5.79712602074737e+120,
4.75364333701284e+122,
3.94552396972066e+124,
3.31424013456535e+126,
2.81710411438055e+128,
2.42270953836727e+130,
2.10775729837953e+132,
1.85482642257398e+134,
1.65079551609085e+136,
1.48571596448176e+138,
1.3520015276784e+140,
1.24384140546413e+142,
1.15677250708164e+144,
1.08736615665674e+146,
1.03299784882391e+148,
9.91677934870949e+149,
9.61927596824821e+151,
9.42689044888324e+153,
9.33262154439441e+155,
9.33262154439441e+157,
9.42594775983835e+159,
9.61446671503512e+161,
9.90290071648618e+163,
1.02990167451456e+166,
1.08139675824029e+168,
1.14628056373471e+170,
1.22652020319614e+172,
1.32464181945183e+174,
1.44385958320249e+176,
1.58824554152274e+178,
1.76295255109024e+180,
1.97450685722107e+182,
2.23119274865981e+184,
2.54355973347219e+186,
2.92509369349301e+188,
3.3931086844519e+190,
3.96993716080872e+192,
4.68452584975429e+194,
5.5745857612076e+196,
6.68950291344912e+198,
8.09429852527344e+200,
9.8750442008336e+202,
1.21463043670253e+205,
1.50614174151114e+207,
1.88267717688893e+209,
2.37217324288005e+211,
3.01266001845766e+213,
3.8562048236258e+215,
4.97450422247729e+217,
6.46685548922047e+219,
8.47158069087882e+221,
1.118248651196e+224,
1.48727070609069e+226,
1.99294274616152e+228,
2.69047270731805e+230,
3.65904288195255e+232,
5.01288874827499e+234,
6.91778647261949e+236,
9.61572319694109e+238,
1.34620124757175e+241,
1.89814375907617e+243,
2.69536413788816e+245,
3.85437071718007e+247,
5.5502938327393e+249,
8.04792605747199e+251,
1.17499720439091e+254,
1.72724589045464e+256,
2.55632391787286e+258,
3.80892263763057e+260,
5.71338395644585e+262,
8.62720977423323e+264,
1.31133588568345e+267,
2.00634390509568e+269,
3.08976961384735e+271,
4.78914290146339e+273,
7.47106292628289e+275,
1.17295687942641e+278,
1.85327186949373e+280,
2.94670227249504e+282,
4.71472363599206e+284,
7.59070505394721e+286,
1.22969421873945e+289,
2.0044015765453e+291,
3.28721858553429e+293,
5.42391066613159e+295,
9.00369170577843e+297,
1.503616514865e+300, // nmaxfactorial = 167
};

View File

@ -56,9 +56,6 @@ class ComputeOrientOrderAtom : public Compute {
double polar_prefactor(int, int, double); double polar_prefactor(int, int, double);
double associated_legendre(int, int, double); double associated_legendre(int, int, double);
static const int nmaxfactorial = 167;
static const double nfac_table[];
double factorial(int);
virtual void init_clebsch_gordan(); virtual void init_clebsch_gordan();
double *cglist; // Clebsch-Gordan coeffs double *cglist; // Clebsch-Gordan coeffs
int idxcg_max; int idxcg_max;

View File

@ -251,31 +251,32 @@ class Fix : protected Pointers {
}; };
namespace FixConst { namespace FixConst {
static const int INITIAL_INTEGRATE = 1<<0; enum {
static const int POST_INTEGRATE = 1<<1; INITIAL_INTEGRATE = 1<<0,
static const int PRE_EXCHANGE = 1<<2; POST_INTEGRATE = 1<<1,
static const int PRE_NEIGHBOR = 1<<3; PRE_EXCHANGE = 1<<2,
static const int POST_NEIGHBOR = 1<<4; PRE_NEIGHBOR = 1<<3,
static const int PRE_FORCE = 1<<5; POST_NEIGHBOR = 1<<4,
static const int PRE_REVERSE = 1<<6; PRE_FORCE = 1<<5,
static const int POST_FORCE = 1<<7; PRE_REVERSE = 1<<6,
static const int FINAL_INTEGRATE = 1<<8; POST_FORCE = 1<<7,
static const int END_OF_STEP = 1<<9; FINAL_INTEGRATE = 1<<8,
static const int POST_RUN = 1<<10; END_OF_STEP = 1<<9,
static const int THERMO_ENERGY = 1<<11; POST_RUN = 1<<10,
static const int INITIAL_INTEGRATE_RESPA = 1<<12; THERMO_ENERGY = 1<<11,
static const int POST_INTEGRATE_RESPA = 1<<13; INITIAL_INTEGRATE_RESPA = 1<<12,
static const int PRE_FORCE_RESPA = 1<<14; POST_INTEGRATE_RESPA = 1<<13,
static const int POST_FORCE_RESPA = 1<<15; PRE_FORCE_RESPA = 1<<14,
static const int FINAL_INTEGRATE_RESPA = 1<<16; POST_FORCE_RESPA = 1<<15,
static const int MIN_PRE_EXCHANGE = 1<<17; FINAL_INTEGRATE_RESPA = 1<<16,
static const int MIN_PRE_NEIGHBOR = 1<<18; MIN_PRE_EXCHANGE = 1<<17,
static const int MIN_POST_NEIGHBOR = 1<<19; MIN_PRE_NEIGHBOR = 1<<18,
static const int MIN_PRE_FORCE = 1<<20; MIN_POST_NEIGHBOR = 1<<19,
static const int MIN_PRE_REVERSE = 1<<21; MIN_PRE_FORCE = 1<<20,
static const int MIN_POST_FORCE = 1<<22; MIN_PRE_REVERSE = 1<<21,
static const int MIN_ENERGY = 1<<23; MIN_POST_FORCE = 1<<22,
static const int FIX_CONST_LAST = 1<<24; MIN_ENERGY = 1<<23
};
} }
} }

View File

@ -305,16 +305,17 @@ double FixTempCSLD::compute_scalar()
void FixTempCSLD::write_restart(FILE *fp) void FixTempCSLD::write_restart(FILE *fp)
{ {
int nsize = (98+2+3)*comm->nprocs+2; // pRNG state per proc + nprocs + energy const int PRNGSIZE = 98+2+3;
int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + energy
double *list = nullptr; double *list = nullptr;
if (comm->me == 0) { if (comm->me == 0) {
list = new double[nsize]; list = new double[nsize];
list[0] = energy; list[0] = energy;
list[1] = comm->nprocs; list[1] = comm->nprocs;
} }
double state[103]; double state[PRNGSIZE];
random->get_state(state); random->get_state(state);
MPI_Gather(state,103,MPI_DOUBLE,list+2,103*comm->nprocs,MPI_DOUBLE,0,world); MPI_Gather(state,PRNGSIZE,MPI_DOUBLE,list+2,PRNGSIZE,MPI_DOUBLE,0,world);
if (comm->me == 0) { if (comm->me == 0) {
int size = nsize * sizeof(double); int size = nsize * sizeof(double);

View File

@ -338,16 +338,17 @@ double FixTempCSVR::compute_scalar()
void FixTempCSVR::write_restart(FILE *fp) void FixTempCSVR::write_restart(FILE *fp)
{ {
int nsize = (98+2+3)*comm->nprocs+2; // pRNG state per proc + nprocs + energy const int PRNGSIZE = 98+2+3;
int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + energy
double *list = nullptr; double *list = nullptr;
if (comm->me == 0) { if (comm->me == 0) {
list = new double[nsize]; list = new double[nsize];
list[0] = energy; list[0] = energy;
list[1] = comm->nprocs; list[1] = comm->nprocs;
} }
double state[103]; double state[PRNGSIZE];
random->get_state(state); random->get_state(state);
MPI_Gather(state,103,MPI_DOUBLE,list+2,103*comm->nprocs,MPI_DOUBLE,0,world); MPI_Gather(state,PRNGSIZE,MPI_DOUBLE,list+2,PRNGSIZE,MPI_DOUBLE,0,world);
if (comm->me == 0) { if (comm->me == 0) {
int size = nsize * sizeof(double); int size = nsize * sizeof(double);

View File

@ -17,18 +17,18 @@
namespace LAMMPS_NS { namespace LAMMPS_NS {
namespace MathConst { namespace MathConst {
static const double THIRD = 1.0/3.0; static constexpr double THIRD = 1.0/3.0;
static const double TWOTHIRDS = 2.0/3.0; static constexpr double TWOTHIRDS = 2.0/3.0;
static const double MY_PI = 3.14159265358979323846; // pi static constexpr double MY_PI = 3.14159265358979323846; // pi
static const double MY_2PI = 6.28318530717958647692; // 2pi static constexpr double MY_2PI = 6.28318530717958647692; // 2pi
static const double MY_3PI = 9.42477796076937971538; // 3pi static constexpr double MY_3PI = 9.42477796076937971538; // 3pi
static const double MY_4PI = 12.56637061435917295384; // 4pi static constexpr double MY_4PI = 12.56637061435917295384; // 4pi
static const double MY_PI2 = 1.57079632679489661923; // pi/2 static constexpr double MY_PI2 = 1.57079632679489661923; // pi/2
static const double MY_PI4 = 0.78539816339744830962; // pi/4 static constexpr double MY_PI4 = 0.78539816339744830962; // pi/4
static const double MY_PIS = 1.77245385090551602729; // sqrt(pi) static constexpr double MY_PIS = 1.77245385090551602729; // sqrt(pi)
static const double MY_ISPI4 = 1.12837916709551257389; // 1/sqrt(pi/4) static constexpr double MY_ISPI4 = 1.12837916709551257389; // 1/sqrt(pi/4)
static const double MY_SQRT2 = 1.41421356237309504880; // sqrt(2) static constexpr double MY_SQRT2 = 1.41421356237309504880; // sqrt(2)
static const double MY_CBRT2 = 1.25992104989487316476; // 2*(1/3) static constexpr double MY_CBRT2 = 1.25992104989487316476; // 2*(1/3)
} }
} }

View File

@ -1,9 +1,203 @@
#include "math_special.h" #include "math_special.h"
#include <cmath> #include <cmath>
#include <cstdint> // IWYU pragma: keep #include <cstdint> // IWYU pragma: keep
#include <limits>
#include "error.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
static constexpr int nmaxfactorial = 167;
/* ----------------------------------------------------------------------
factorial n table, size nmaxfactorial+1
------------------------------------------------------------------------- */
static const double nfac_table[] = {
1,
1,
2,
6,
24,
120,
720,
5040,
40320,
362880,
3628800,
39916800,
479001600,
6227020800,
87178291200,
1307674368000,
20922789888000,
355687428096000,
6.402373705728e+15,
1.21645100408832e+17,
2.43290200817664e+18,
5.10909421717094e+19,
1.12400072777761e+21,
2.5852016738885e+22,
6.20448401733239e+23,
1.5511210043331e+25,
4.03291461126606e+26,
1.08888694504184e+28,
3.04888344611714e+29,
8.8417619937397e+30,
2.65252859812191e+32,
8.22283865417792e+33,
2.63130836933694e+35,
8.68331761881189e+36,
2.95232799039604e+38,
1.03331479663861e+40,
3.71993326789901e+41,
1.37637530912263e+43,
5.23022617466601e+44,
2.03978820811974e+46,
8.15915283247898e+47,
3.34525266131638e+49,
1.40500611775288e+51,
6.04152630633738e+52,
2.65827157478845e+54,
1.1962222086548e+56,
5.50262215981209e+57,
2.58623241511168e+59,
1.24139155925361e+61,
6.08281864034268e+62,
3.04140932017134e+64,
1.55111875328738e+66,
8.06581751709439e+67,
4.27488328406003e+69,
2.30843697339241e+71,
1.26964033536583e+73,
7.10998587804863e+74,
4.05269195048772e+76,
2.35056133128288e+78,
1.3868311854569e+80,
8.32098711274139e+81,
5.07580213877225e+83,
3.14699732603879e+85,
1.98260831540444e+87,
1.26886932185884e+89,
8.24765059208247e+90,
5.44344939077443e+92,
3.64711109181887e+94,
2.48003554243683e+96,
1.71122452428141e+98,
1.19785716699699e+100,
8.50478588567862e+101,
6.12344583768861e+103,
4.47011546151268e+105,
3.30788544151939e+107,
2.48091408113954e+109,
1.88549470166605e+111,
1.45183092028286e+113,
1.13242811782063e+115,
8.94618213078297e+116,
7.15694570462638e+118,
5.79712602074737e+120,
4.75364333701284e+122,
3.94552396972066e+124,
3.31424013456535e+126,
2.81710411438055e+128,
2.42270953836727e+130,
2.10775729837953e+132,
1.85482642257398e+134,
1.65079551609085e+136,
1.48571596448176e+138,
1.3520015276784e+140,
1.24384140546413e+142,
1.15677250708164e+144,
1.08736615665674e+146,
1.03299784882391e+148,
9.91677934870949e+149,
9.61927596824821e+151,
9.42689044888324e+153,
9.33262154439441e+155,
9.33262154439441e+157,
9.42594775983835e+159,
9.61446671503512e+161,
9.90290071648618e+163,
1.02990167451456e+166,
1.08139675824029e+168,
1.14628056373471e+170,
1.22652020319614e+172,
1.32464181945183e+174,
1.44385958320249e+176,
1.58824554152274e+178,
1.76295255109024e+180,
1.97450685722107e+182,
2.23119274865981e+184,
2.54355973347219e+186,
2.92509369349301e+188,
3.3931086844519e+190,
3.96993716080872e+192,
4.68452584975429e+194,
5.5745857612076e+196,
6.68950291344912e+198,
8.09429852527344e+200,
9.8750442008336e+202,
1.21463043670253e+205,
1.50614174151114e+207,
1.88267717688893e+209,
2.37217324288005e+211,
3.01266001845766e+213,
3.8562048236258e+215,
4.97450422247729e+217,
6.46685548922047e+219,
8.47158069087882e+221,
1.118248651196e+224,
1.48727070609069e+226,
1.99294274616152e+228,
2.69047270731805e+230,
3.65904288195255e+232,
5.01288874827499e+234,
6.91778647261949e+236,
9.61572319694109e+238,
1.34620124757175e+241,
1.89814375907617e+243,
2.69536413788816e+245,
3.85437071718007e+247,
5.5502938327393e+249,
8.04792605747199e+251,
1.17499720439091e+254,
1.72724589045464e+256,
2.55632391787286e+258,
3.80892263763057e+260,
5.71338395644585e+262,
8.62720977423323e+264,
1.31133588568345e+267,
2.00634390509568e+269,
3.08976961384735e+271,
4.78914290146339e+273,
7.47106292628289e+275,
1.17295687942641e+278,
1.85327186949373e+280,
2.94670227249504e+282,
4.71472363599206e+284,
7.59070505394721e+286,
1.22969421873945e+289,
2.0044015765453e+291,
3.28721858553429e+293,
5.42391066613159e+295,
9.00369170577843e+297,
1.503616514865e+300, // nmaxfactorial = 167
};
/* ----------------------------------------------------------------------
factorial n vial lookup from precomputed table
------------------------------------------------------------------------- */
double MathSpecial::factorial(const int n)
{
if (n < 0 || n > nmaxfactorial)
return std::numeric_limits<double>::quiet_NaN();
return nfac_table[n];
}
/* Library libcerf: /* Library libcerf:
* Compute complex error functions, based on a new implementation of * Compute complex error functions, based on a new implementation of
* Faddeeva's w_of_z. Also provide Dawson and Voigt functions. * Faddeeva's w_of_z. Also provide Dawson and Voigt functions.

View File

@ -20,6 +20,10 @@ namespace LAMMPS_NS {
namespace MathSpecial { namespace MathSpecial {
// tabulated factorial function
extern double factorial(const int);
// support function for scaled error function complement // support function for scaled error function complement
extern double erfcx_y100(const double y100); extern double erfcx_y100(const double y100);

View File

@ -235,49 +235,55 @@ class Neighbor : protected Pointers {
}; };
namespace NeighConst { namespace NeighConst {
static const int NB_INTEL = 1<<0; enum {
static const int NB_KOKKOS_DEVICE = 1<<1; NB_INTEL = 1<<0,
static const int NB_KOKKOS_HOST = 1<<2; NB_KOKKOS_DEVICE = 1<<1,
static const int NB_SSA = 1<<3; NB_KOKKOS_HOST = 1<<2,
NB_SSA = 1<<3
};
static const int NS_BIN = 1<<0; enum {
static const int NS_MULTI = 1<<1; NS_BIN = 1<<0,
static const int NS_HALF = 1<<2; NS_MULTI = 1<<1,
static const int NS_FULL = 1<<3; NS_HALF = 1<<2,
static const int NS_2D = 1<<4; NS_FULL = 1<<3,
static const int NS_3D = 1<<5; NS_2D = 1<<4,
static const int NS_NEWTON = 1<<6; NS_3D = 1<<5,
static const int NS_NEWTOFF = 1<<7; NS_NEWTON = 1<<6,
static const int NS_ORTHO = 1<<8; NS_NEWTOFF = 1<<7,
static const int NS_TRI = 1<<9; NS_ORTHO = 1<<8,
static const int NS_GHOST = 1<<10; NS_TRI = 1<<9,
static const int NS_SSA = 1<<11; NS_GHOST = 1<<10,
NS_SSA = 1<<11
};
static const int NP_NSQ = 1<<0; enum {
static const int NP_BIN = 1<<1; NP_NSQ = 1<<0,
static const int NP_MULTI = 1<<2; NP_BIN = 1<<1,
static const int NP_HALF = 1<<3; NP_MULTI = 1<<2,
static const int NP_FULL = 1<<4; NP_HALF = 1<<3,
static const int NP_ORTHO = 1<<5; NP_FULL = 1<<4,
static const int NP_TRI = 1<<6; NP_ORTHO = 1<<5,
static const int NP_ATOMONLY = 1<<7; NP_TRI = 1<<6,
static const int NP_MOLONLY = 1<<8; NP_ATOMONLY = 1<<7,
static const int NP_NEWTON = 1<<9; NP_MOLONLY = 1<<8,
static const int NP_NEWTOFF = 1<<10; NP_NEWTON = 1<<9,
static const int NP_GHOST = 1<<11; NP_NEWTOFF = 1<<10,
static const int NP_SIZE = 1<<12; NP_GHOST = 1<<11,
static const int NP_ONESIDE = 1<<13; NP_SIZE = 1<<12,
static const int NP_RESPA = 1<<14; NP_ONESIDE = 1<<13,
static const int NP_BOND = 1<<15; NP_RESPA = 1<<14,
static const int NP_OMP = 1<<16; NP_BOND = 1<<15,
static const int NP_INTEL = 1<<17; NP_OMP = 1<<16,
static const int NP_KOKKOS_DEVICE = 1<<18; NP_INTEL = 1<<17,
static const int NP_KOKKOS_HOST = 1<<19; NP_KOKKOS_DEVICE = 1<<18,
static const int NP_SSA = 1<<20; NP_KOKKOS_HOST = 1<<19,
static const int NP_COPY = 1<<21; NP_SSA = 1<<20,
static const int NP_SKIP = 1<<22; NP_COPY = 1<<21,
static const int NP_HALF_FULL = 1<<23; NP_SKIP = 1<<22,
static const int NP_OFF2ON = 1<<24; NP_HALF_FULL = 1<<23,
NP_OFF2ON = 1<<24
};
} }
} }

View File

@ -0,0 +1,117 @@
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
#include "npair_half_size_multi_newtoff.h"
#include "neigh_list.h"
#include "atom.h"
#include "atom_vec.h"
#include "my_page.h"
#include "error.h"
using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
NPairHalfSizeMultiNewtoff::NPairHalfSizeMultiNewtoff(LAMMPS *lmp) : NPair(lmp) {}
/* ----------------------------------------------------------------------
size particles
binned neighbor list construction with partial Newton's 3rd law
each owned atom i checks own bin and other bins in stencil
multi-type stencil is itype dependent and is distance checked
pair stored once if i,j are both owned and i < j
pair stored by me if j is ghost (also stored by proc owning j)
------------------------------------------------------------------------- */
void NPairHalfSizeMultiNewtoff::build(NeighList *list)
{
int i,j,k,m,n,itype,jtype,ibin,ns;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
double radi,radsum,cutdistsq;
int *neighptr,*s;
double *cutsq,*distsq;
double **x = atom->x;
double *radius = atom->radius;
int *type = atom->type;
int *mask = atom->mask;
tagint *molecule = atom->molecule;
int nlocal = atom->nlocal;
if (includegroup) nlocal = atom->nfirst;
int history = list->history;
int *ilist = list->ilist;
int *numneigh = list->numneigh;
int **firstneigh = list->firstneigh;
MyPage<int> *ipage = list->ipage;
int mask_history = 3 << SBBITS;
int inum = 0;
ipage->reset();
for (i = 0; i < nlocal; i++) {
n = 0;
neighptr = ipage->vget();
itype = type[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
radi = radius[i];
// loop over all atoms in other bins in stencil including self
// only store pair if i < j
// skip if i,j neighbor cutoff is less than bin distance
// stores own/own pairs only once
// stores own/ghost pairs on both procs
ibin = atom2bin[i];
s = stencil_multi[itype];
distsq = distsq_multi[itype];
cutsq = cutneighsq[itype];
ns = nstencil_multi[itype];
for (k = 0; k < ns; k++) {
for (j = binhead[ibin+s[k]]; j >= 0; j = bins[j]) {
if (j <= i) continue;
jtype = type[j];
if (cutsq[jtype] < distsq[k]) continue;
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
}
ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
ipage->vgot(n);
if (ipage->status())
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
list->inum = inum;
}

View File

@ -0,0 +1,46 @@
/* -*- 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 NPAIR_CLASS
NPairStyle(half/size/multi/newtoff,
NPairHalfSizeMultiNewtoff,
NP_HALF | NP_SIZE | NP_MULTI | NP_NEWTOFF | NP_ORTHO | NP_TRI)
#else
#ifndef LMP_NPAIR_HALF_SIZE_MULTI_NEWTOFF_H
#define LMP_NPAIR_HALF_SIZE_MULTI_NEWTOFF_H
#include "npair.h"
namespace LAMMPS_NS {
class NPairHalfSizeMultiNewtoff : public NPair {
public:
NPairHalfSizeMultiNewtoff(class LAMMPS *);
~NPairHalfSizeMultiNewtoff() {}
void build(class NeighList *);
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Neighbor list overflow, boost neigh_modify one
UNDOCUMENTED
*/

View File

@ -0,0 +1,143 @@
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
#include "npair_half_size_multi_newton.h"
#include "neigh_list.h"
#include "atom.h"
#include "atom_vec.h"
#include "my_page.h"
#include "error.h"
using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
NPairHalfSizeMultiNewton::NPairHalfSizeMultiNewton(LAMMPS *lmp) : NPair(lmp) {}
/* ----------------------------------------------------------------------
size particles
binned neighbor list construction with full Newton's 3rd law
each owned atom i checks its own bin and other bins in Newton stencil
multi-type stencil is itype dependent and is distance checked
every pair stored exactly once by some processor
------------------------------------------------------------------------- */
void NPairHalfSizeMultiNewton::build(NeighList *list)
{
int i,j,k,m,n,itype,jtype,ibin,ns;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
double radi,radsum,cutdistsq;
int *neighptr,*s;
double *cutsq,*distsq;
double **x = atom->x;
double *radius = atom->radius;
int *type = atom->type;
int *mask = atom->mask;
tagint *molecule = atom->molecule;
int nlocal = atom->nlocal;
if (includegroup) nlocal = atom->nfirst;
int history = list->history;
int *ilist = list->ilist;
int *numneigh = list->numneigh;
int **firstneigh = list->firstneigh;
MyPage<int> *ipage = list->ipage;
int mask_history = 3 << SBBITS;
int inum = 0;
ipage->reset();
for (i = 0; i < nlocal; i++) {
n = 0;
neighptr = ipage->vget();
itype = type[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
radi = radius[i];
// loop over rest of atoms in i's bin, ghosts are at end of linked list
// if j is owned atom, store it, since j is beyond i in linked list
// if j is ghost, only store if j coords are "above and to the right" of i
for (j = bins[i]; j >= 0; j = bins[j]) {
if (j >= nlocal) {
if (x[j][2] < ztmp) continue;
if (x[j][2] == ztmp) {
if (x[j][1] < ytmp) continue;
if (x[j][1] == ytmp && x[j][0] < xtmp) continue;
}
}
jtype = type[j];
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
// loop over all atoms in other bins in stencil, store every pair
// skip if i,j neighbor cutoff is less than bin distance
ibin = atom2bin[i];
s = stencil_multi[itype];
distsq = distsq_multi[itype];
cutsq = cutneighsq[itype];
ns = nstencil_multi[itype];
for (k = 0; k < ns; k++) {
for (j = binhead[ibin+s[k]]; j >= 0; j = bins[j]) {
jtype = type[j];
if (cutsq[jtype] < distsq[k]) continue;
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
}
ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
ipage->vgot(n);
if (ipage->status())
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
list->inum = inum;
}

View File

@ -0,0 +1,46 @@
/* -*- 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 NPAIR_CLASS
NPairStyle(half/size/multi/newton,
NPairHalfSizeMultiNewton,
NP_HALF | NP_SIZE | NP_MULTI | NP_NEWTON | NP_ORTHO)
#else
#ifndef LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_H
#define LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_H
#include "npair.h"
namespace LAMMPS_NS {
class NPairHalfSizeMultiNewton : public NPair {
public:
NPairHalfSizeMultiNewton(class LAMMPS *);
~NPairHalfSizeMultiNewton() {}
void build(class NeighList *);
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Neighbor list overflow, boost neigh_modify one
UNDOCUMENTED
*/

View File

@ -0,0 +1,125 @@
/* ----------------------------------------------------------------------
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.
------------------------------------------------------------------------- */
#include "npair_half_size_multi_newton_tri.h"
#include "neigh_list.h"
#include "atom.h"
#include "atom_vec.h"
#include "my_page.h"
#include "error.h"
using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */
NPairHalfSizeMultiNewtonTri::NPairHalfSizeMultiNewtonTri(LAMMPS *lmp) : NPair(lmp) {}
/* ----------------------------------------------------------------------
binned neighbor list construction with Newton's 3rd law for triclinic
each owned atom i checks its own bin and other bins in triclinic stencil
multi-type stencil is itype dependent and is distance checked
every pair stored exactly once by some processor
------------------------------------------------------------------------- */
void NPairHalfSizeMultiNewtonTri::build(NeighList *list)
{
int i,j,k,m,n,itype,jtype,ibin,ns;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
double radi,radsum,cutdistsq;
int *neighptr,*s;
double *cutsq,*distsq;
double **x = atom->x;
double *radius = atom->radius;
int *type = atom->type;
int *mask = atom->mask;
tagint *molecule = atom->molecule;
int nlocal = atom->nlocal;
if (includegroup) nlocal = atom->nfirst;
int history = list->history;
int *ilist = list->ilist;
int *numneigh = list->numneigh;
int **firstneigh = list->firstneigh;
MyPage<int> *ipage = list->ipage;
int mask_history = 3 << SBBITS;
int inum = 0;
ipage->reset();
for (i = 0; i < nlocal; i++) {
n = 0;
neighptr = ipage->vget();
itype = type[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
radi = radius[i];
// loop over all atoms in bins, including self, in stencil
// skip if i,j neighbor cutoff is less than bin distance
// bins below self are excluded from stencil
// pairs for atoms j "below" i are excluded
// below = lower z or (equal z and lower y) or (equal zy and lower x)
// (equal zyx and j <= i)
// latter excludes self-self interaction but allows superposed atoms
ibin = atom2bin[i];
s = stencil_multi[itype];
distsq = distsq_multi[itype];
cutsq = cutneighsq[itype];
ns = nstencil_multi[itype];
for (k = 0; k < ns; k++) {
for (j = binhead[ibin+s[k]]; j >= 0; j = bins[j]) {
jtype = type[j];
if (cutsq[jtype] < distsq[k]) continue;
if (x[j][2] < ztmp) continue;
if (x[j][2] == ztmp) {
if (x[j][1] < ytmp) continue;
if (x[j][1] == ytmp) {
if (x[j][0] < xtmp) continue;
if (x[j][0] == xtmp && j <= i) continue;
}
}
if (exclude && exclusion(i,j,itype,jtype,mask,molecule)) continue;
delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j];
cutdistsq = (radsum+skin) * (radsum+skin);
if (rsq <= cutdistsq) {
if (history && rsq < radsum*radsum)
neighptr[n++] = j ^ mask_history;
else
neighptr[n++] = j;
}
}
}
ilist[inum++] = i;
firstneigh[i] = neighptr;
numneigh[i] = n;
ipage->vgot(n);
if (ipage->status())
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
list->inum = inum;
}

View File

@ -0,0 +1,46 @@
/* -*- 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 NPAIR_CLASS
NPairStyle(half/size/multi/newton/tri,
NPairHalfSizeMultiNewtonTri,
NP_HALF | NP_SIZE | NP_MULTI | NP_NEWTON | NP_TRI)
#else
#ifndef LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_TRI_H
#define LMP_NPAIR_HALF_SIZE_MULTI_NEWTON_TRI_H
#include "npair.h"
namespace LAMMPS_NS {
class NPairHalfSizeMultiNewtonTri : public NPair {
public:
NPairHalfSizeMultiNewtonTri(class LAMMPS *);
~NPairHalfSizeMultiNewtonTri() {}
void build(class NeighList *);
};
}
#endif
#endif
/* ERROR/WARNING messages:
E: Neighbor list overflow, boost neigh_modify one
UNDOCUMENTED
*/

View File

@ -36,7 +36,7 @@ class PairCoulStreitz : public Pair {
double memory_usage(); double memory_usage();
virtual void *extract(const char *, int &); virtual void *extract(const char *, int &);
static const int NPARAMS_PER_LINE = 6; static constexpr int NPARAMS_PER_LINE = 6;
protected: protected:
struct Param { struct Param {

View File

@ -26,6 +26,7 @@
#include "memory.h" #include "memory.h"
#include "error.h" #include "error.h"
#include "pair_lj_cubic_const.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace PairLJCubicConstants; using namespace PairLJCubicConstants;

View File

@ -46,20 +46,6 @@ class PairLJCubic : public Pair {
void allocate(); void allocate();
}; };
namespace PairLJCubicConstants {
// LJ quantities scaled by epsilon and rmin = sigma*2^1/6
static const double RT6TWO = 1.1224621; // 2^1/6
static const double SS = 1.1086834; // inflection point (13/7)^1/6
static const double PHIS = -0.7869823; // energy at s
static const double DPHIDS = 2.6899009; // gradient at s
static const double A3 = 27.93357; // cubic coefficient
static const double SM = 1.5475375; // cubic cutoff = s*67/48
}
} }
#endif #endif

48
src/pair_lj_cubic_const.h Normal file
View File

@ -0,0 +1,48 @@
/* -*- 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.
------------------------------------------------------------------------- */
#ifndef LMP_PAIR_LJ_CUBIC_CONST_H
#define LMP_PAIR_LJ_CUBIC_CONST_H
namespace LAMMPS_NS {
namespace PairLJCubicConstants {
// LJ quantities scaled by epsilon and rmin = sigma*2^1/6
static constexpr double RT6TWO = 1.1224620483093730; // 2^1/6
static constexpr double SS = 1.1086834179687215; // inflection point (13/7)^1/6
static constexpr double PHIS = -0.7869822485207097; // energy at s
static constexpr double DPHIDS = 2.6899008972047196; // gradient at s
static constexpr double A3 = 27.9335700460986445; // cubic coefficient
static constexpr double SM = 1.5475372709146737; // cubic cutoff = s*67/48}
}
}
#endif
// python script to compute the constants
//
// sixth = 1.0/6.0
// rmin = pow(2.0,sixth)
// rs = pow(26.0/7.0,sixth)
// ss = rs/rmin
// pow6 = pow(1.0/rs,6.0)
// phis = 4.0*pow6*(pow6-1.0)
// dpds = -24.0*pow6*(2.0*pow6-1.0)/rs*rmin
// a3 = 8.0*pow(dpds,3)/(9.0*phis*phis)
// sm = 67.0/48.0*ss
// print("static constexpr double RT6TWO = %19.16f; // 2^1/6" % rmin)
// print("static constexpr double SS = %19.16f; // inflection point (13/7)^1/6" % ss)
// print("static constexpr double PHIS = %19.16f; // energy at s" % phis)
// print("static constexpr double DPHIDS = %19.16f; // gradient at s" % dpds)
// print("static constexpr double A3 = %19.16f; // cubic coefficient" % a3)
// print("static constexpr double SM = %19.16f; // cubic cutoff = s*67/48" % sm)

View File

@ -16,15 +16,18 @@
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
#include "pair_zbl.h" #include "pair_zbl.h"
#include <cmath> #include <cmath>
#include "atom.h" #include "atom.h"
#include "comm.h" #include "comm.h"
#include "error.h"
#include "force.h" #include "force.h"
#include "memory.h"
#include "neighbor.h" #include "neighbor.h"
#include "neigh_list.h" #include "neigh_list.h"
#include "memory.h"
#include "error.h"
#include "pair_zbl_const.h"
// From J.F. Zeigler, J. P. Biersack and U. Littmark, // From J.F. Zeigler, J. P. Biersack and U. Littmark,
// "The Stopping and Range of Ions in Matter" volume 1, Pergamon, 1985. // "The Stopping and Range of Ions in Matter" volume 1, Pergamon, 1985.

View File

@ -54,23 +54,6 @@ class PairZBL : public Pair {
double d2zbldr2(double, int, int); double d2zbldr2(double, int, int);
void set_coeff(int, int, double, double); void set_coeff(int, int, double, double);
}; };
namespace PairZBLConstants {
// ZBL constants
static const double pzbl = 0.23;
static const double a0 = 0.46850;
static const double c1 = 0.02817;
static const double c2 = 0.28022;
static const double c3 = 0.50986;
static const double c4 = 0.18175;
static const double d1 = 0.20162;
static const double d2 = 0.40290;
static const double d3 = 0.94229;
static const double d4 = 3.19980;
}
} }
#endif #endif

34
src/pair_zbl_const.h Normal file
View File

@ -0,0 +1,34 @@
/* -*- 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.
------------------------------------------------------------------------- */
#ifndef LMP_PAIR_ZBL_CONST_H
#define LMP_PAIR_ZBL_CONST_H
namespace LAMMPS_NS {
namespace PairZBLConstants {
// ZBL constants
static constexpr double pzbl = 0.23;
static constexpr double a0 = 0.46850;
static constexpr double c1 = 0.02817;
static constexpr double c2 = 0.28022;
static constexpr double c3 = 0.50986;
static constexpr double c4 = 0.18175;
static constexpr double d1 = 0.20162;
static constexpr double d2 = 0.40290;
static constexpr double d3 = 0.94229;
static constexpr double d4 = 3.19980;
}
}
#endif

View File

@ -27,7 +27,7 @@ namespace LAMMPS_NS
class TextFileReader { class TextFileReader {
std::string filename; std::string filename;
std::string filetype; std::string filetype;
static const int MAXLINE = 1024; static constexpr int MAXLINE = 1024;
char line[MAXLINE]; char line[MAXLINE];
FILE *fp; FILE *fp;

View File

@ -173,7 +173,7 @@ TEST(MPI, multi_partition)
MPI_Comm_size(MPI_COMM_WORLD, &nprocs); MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
MPI_Comm_rank(MPI_COMM_WORLD, &me); MPI_Comm_rank(MPI_COMM_WORLD, &me);
const char *args[] = {"LAMMPS_test", "-log", "none", "-partition", "2x2", const char *args[] = {"LAMMPS_test", "-log", "none", "-partition", "4x1",
"-echo", "screen", "-nocite", "-in", "none"}; "-echo", "screen", "-nocite", "-in", "none"};
char **argv = (char **)args; char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *); int argc = sizeof(args) / sizeof(char *);
@ -183,19 +183,15 @@ TEST(MPI, multi_partition)
lammps_command(lmp, "atom_style atomic"); lammps_command(lmp, "atom_style atomic");
lammps_command(lmp, "region box block 0 2 0 2 0 2"); lammps_command(lmp, "region box block 0 2 0 2 0 2");
lammps_command(lmp, "create_box 1 box"); lammps_command(lmp, "create_box 1 box");
lammps_command(lmp, "variable partition universe 1 2"); lammps_command(lmp, "variable partition universe 1 2 3 4");
EXPECT_EQ(lammps_extract_setting(lmp, "universe_size"), nprocs); EXPECT_EQ(lammps_extract_setting(lmp, "universe_size"), nprocs);
EXPECT_EQ(lammps_extract_setting(lmp, "universe_rank"), me); EXPECT_EQ(lammps_extract_setting(lmp, "universe_rank"), me);
EXPECT_EQ(lammps_extract_setting(lmp, "world_size"), nprocs / 2); EXPECT_EQ(lammps_extract_setting(lmp, "world_size"), 1);
EXPECT_EQ(lammps_extract_setting(lmp, "world_rank"), me % 2); EXPECT_EQ(lammps_extract_setting(lmp, "world_rank"), 0);
char *part_id = (char *)lammps_extract_variable(lmp, "partition", nullptr); char *part_id = (char *)lammps_extract_variable(lmp, "partition", nullptr);
if (me < 2) { ASSERT_THAT(part_id, StrEq(std::to_string(me+1)));
ASSERT_THAT(part_id, StrEq("1"));
} else {
ASSERT_THAT(part_id, StrEq("2"));
}
lammps_close(lmp); lammps_close(lmp);
}; };

View File

@ -1,7 +1,7 @@
--- ---
lammps_version: 24 Aug 2020 lammps_version: 24 Aug 2020
date_generated: Tue Sep 15 09:44:41 202 date_generated: Tue Sep 15 09:44:41 202
epsilon: 2e-14 epsilon: 4e-14
prerequisites: ! | prerequisites: ! |
atom full atom full
fix oneway fix oneway

View File

@ -1,7 +1,7 @@
--- ---
lammps_version: 24 Aug 2020 lammps_version: 30 Nov 2020
date_generated: Tue Sep 15 09:44:15 202 date_generated: Fri Dec 18 22:30:42 202
epsilon: 2e-13 epsilon: 1e-13
prerequisites: ! | prerequisites: ! |
atom full atom full
pair lj/cubic pair lj/cubic
@ -19,72 +19,72 @@ pair_coeff: ! |
5 5 0.015 3.1 5 5 0.015 3.1
extract: ! "" extract: ! ""
natoms: 29 natoms: 29
init_vdwl: 166.005266308562 init_vdwl: 166.005266338502
init_coul: 0 init_coul: 0
init_stress: ! |2- init_stress: ! |2-
4.5860676757908186e+02 4.8091912919212928e+02 1.0767204080701006e+03 -2.1005546139122362e+02 -2.9491286717936713e+00 1.6145675857120941e+02 4.5860676760613336e+02 4.8091912923656065e+02 1.0767204080957010e+03 -2.1005546139899724e+02 -2.9491286649299440e+00 1.6145675855773089e+02
init_forces: ! |2 init_forces: ! |2
1 9.1849370411551270e+00 7.6268937957720553e+01 6.1726872441625311e+01 1 9.1849370400954342e+00 7.6268937960282486e+01 6.1726872440953251e+01
2 2.2858712118514426e+01 1.8809274242266209e+01 -2.6905829837199740e+01 2 2.2858712118557477e+01 1.8809274242339992e+01 -2.6905829837209875e+01
3 -3.2016987482543328e+01 -9.4135849525427091e+01 -3.4799279593035926e+01 3 -3.2016987482586380e+01 -9.4135849525500873e+01 -3.4799279593025787e+01
4 -5.5341015869901478e-01 1.5206999898436971e-01 -3.9418368928369890e-01 4 -5.5341015917141478e-01 1.5206999930354809e-01 -3.9418368927471259e-01
5 -1.8042057425348118e-01 -3.0459951056385326e-01 8.7068483241007189e-01 5 -1.8042057490785210e-01 -3.0459951013385317e-01 8.7068483295449139e-01
6 -2.0038994438822397e+02 2.3344446299945159e+02 2.8487343926572851e+02 6 -2.0038994438983289e+02 2.3344446299845976e+02 2.8487343926562755e+02
7 8.0909912172413883e+00 -7.8410849891085633e+01 -4.3214084684451740e+02 7 8.0909912152126484e+00 -7.8410849888970290e+01 -4.3214084684693944e+02
8 4.7943581255133857e+01 -2.1287511456246008e+01 1.4094503445180061e+02 8 4.7943581253518161e+01 -2.1287511458762772e+01 1.4094503445043824e+02
9 1.1447552368270737e+01 1.2328709806786962e+01 5.0656476982000299e+01 9 1.1447552368334570e+01 1.2328709806695397e+01 5.0656476981938759e+01
10 1.3071496571967870e+02 -1.4589264560693914e+02 -4.4748155922123622e+01 10 1.3071496571895125e+02 -1.4589264560750766e+02 -4.4748155921681018e+01
11 -1.6551880116149281e-01 -4.1534332040572380e-01 -6.8284765241715795e-01 11 -1.6551880130741767e-01 -4.1534332035668403e-01 -6.8284765271893910e-01
12 1.7721533626133388e+00 6.3456329073685158e-01 -8.2372301448028962e-01 12 1.7721533633976141e+00 6.3456328808910523e-01 -8.2372301296999062e-01
13 5.6789360334118277e-01 -2.2634410312439054e-01 -9.7536738055328392e-03 13 5.6789360517170440e-01 -2.2634410322310763e-01 -9.7536744156455080e-03
14 -2.4337021468262635e-01 4.6659433642728905e-02 -6.1110664501270184e-01 14 -2.4337021462140446e-01 4.6659433737428486e-02 -6.1110664517858926e-01
15 -2.1936997101927893e-02 5.9238263972968364e-01 2.1493099548264527e-01 15 -2.1936996885396808e-02 5.9238264018028652e-01 2.1493099542010694e-01
16 1.1121534968449923e+02 -7.8056927924992834e+01 -2.9249212971206231e+02 16 1.1121534968641315e+02 -7.8056927926703892e+01 -2.9249212971001145e+02
17 -1.1020604609843586e+02 7.6481296254913858e+01 2.9430701446263464e+02 17 -1.1020604609471918e+02 7.6481296251709367e+01 2.9430701446464468e+02
18 -1.6570656719723909e-02 -2.7996966177077785e-02 2.6456326954440619e-02 18 -1.6570656172995385e-02 -2.7996961634495443e-02 2.6456324093123054e-02
19 7.4243353115058947e-04 6.3524893127716046e-04 1.8675586277048476e-04 19 7.4243314752471426e-04 6.3524860303507949e-04 1.8675576627109685e-04
20 -7.4243353115058947e-04 -6.3524893127716046e-04 -1.8675586277048476e-04 20 -7.4243314752471426e-04 -6.3524860303507949e-04 -1.8675576627109685e-04
21 -1.1415041486189516e+01 -1.3016363071591645e+01 3.6007276733401099e+01 21 -1.1415041486189516e+01 -1.3016363071591645e+01 3.6007276733401099e+01
22 -1.7227422089792942e+01 -4.1746638094950628e+00 -2.7029162034499002e+01 22 -1.7227422090167991e+01 -4.1746638096674396e+00 -2.7029162034658786e+01
23 2.8642463575982458e+01 1.7191026881086707e+01 -8.9781146989020968e+00 23 2.8642463576357507e+01 1.7191026881259084e+01 -8.9781146987423170e+00
24 5.8150644491939154e+00 -3.3774314134628064e+01 1.7867788752379695e+01 24 5.8150644495067958e+00 -3.3774314132751599e+01 1.7867788754173183e+01
25 -2.3666545027773044e+01 3.8106021846559952e+00 -1.9896269873584632e+01 25 -2.3666545028156815e+01 3.8106021844353424e+00 -1.9896269873795571e+01
26 1.7843812244577855e+01 2.9960339884741117e+01 2.0167430316952100e+00 26 1.7843812244961626e+01 2.9960339884961773e+01 2.0167430319061479e+00
27 8.2825859209946024e+00 -3.6194570066818969e+01 1.4492694351988913e+01 27 8.2825859198612442e+00 -3.6194570067428131e+01 1.4492694352248694e+01
28 -2.8773892796642542e+01 1.2366374307374247e+01 -1.9468877181285176e+01 28 -2.8773892797077618e+01 1.2366374307246014e+01 -1.9468877181493664e+01
29 2.0497044211022661e+01 2.3831279505404666e+01 4.9748677441078746e+00 29 2.0497044211457737e+01 2.3831279505532898e+01 4.9748677443163629e+00
run_vdwl: 164.176727313193 run_vdwl: 164.176727343123
run_coul: 0 run_coul: 0
run_stress: ! |2- run_stress: ! |2-
4.5696669905868288e+02 4.7904134871830234e+02 1.0582153129928076e+03 -2.0794233475912671e+02 -2.0206236780739086e+00 1.5948889983617320e+02 4.5696669908578690e+02 4.7904134876274537e+02 1.0582153130184231e+03 -2.0794233476691440e+02 -2.0206236712071677e+00 1.5948889982266843e+02
run_forces: ! |2 run_forces: ! |2
1 9.2483453892728740e+00 7.5945844239974676e+01 6.1367289784738311e+01 1 9.2483453882111135e+00 7.5945844242532630e+01 6.1367289784064468e+01
2 2.2732852134554197e+01 1.8737493288759875e+01 -2.6669600068587577e+01 2 2.2732852134599622e+01 1.8737493288836678e+01 -2.6669600068598680e+01
3 -3.1964986246503269e+01 -9.3734416340601200e+01 -3.4671255449722295e+01 3 -3.1964986246546172e+01 -9.3734416340673434e+01 -3.4671255449709932e+01
4 -5.5001585694065913e-01 1.5070776418957152e-01 -3.9275226164990551e-01 4 -5.5001585741105219e-01 1.5070776450703974e-01 -3.9275226164162935e-01
5 -1.7969915375876547e-01 -3.0362292678324765e-01 8.6793365916706400e-01 5 -1.7969915441332354e-01 -3.0362292635308630e-01 8.6793365971015668e-01
6 -1.9796806590003533e+02 2.2990162655369829e+02 2.7501790745035373e+02 6 -1.9796806590165264e+02 2.2990162655270728e+02 2.7501790745023669e+02
7 7.9109888171861886e+00 -7.6524641847864288e+01 -4.2032627555539375e+02 7 7.9109888151537930e+00 -7.6524641845739424e+01 -4.2032627555780203e+02
8 4.5992175211879918e+01 -1.9873148355500923e+01 1.3922798515578887e+02 8 4.5992175210254615e+01 -1.9873148358014198e+01 1.3922798515443196e+02
9 1.1403303962503362e+01 1.2231165268449960e+01 5.0409090529604853e+01 9 1.1403303962567387e+01 1.2231165268357689e+01 5.0409090529541324e+01
10 1.3047784096912977e+02 -1.4556513307073578e+02 -4.4756481385998420e+01 10 1.3047784096840681e+02 -1.4556513307130629e+02 -4.4756481385556341e+01
11 -1.6346238349194633e-01 -4.0922088697872150e-01 -6.7303941060221573e-01 11 -1.6346238363922402e-01 -4.0922088693014880e-01 -6.7303941090576869e-01
12 1.7689668574627495e+00 6.3166692380469824e-01 -8.3233385524693426e-01 12 1.7689668582530251e+00 6.3166692115802814e-01 -8.3233385373913271e-01
13 5.6480104331071312e-01 -2.2395151872256039e-01 -9.5790993973179691e-03 13 5.6480104514758445e-01 -2.2395151881989347e-01 -9.5791000096411370e-03
14 -2.4030540495346994e-01 4.5179188229734012e-02 -6.0304369153488313e-01 14 -2.4030540489083232e-01 4.5179188324553629e-02 -6.0304369170302041e-01
15 -2.3220111317111478e-02 5.9408611078423390e-01 2.1676726911830960e-01 15 -2.3220111102749796e-02 5.9408611123108124e-01 2.1676726905655996e-01
16 1.0963039832302356e+02 -7.7096357855469549e+01 -2.8842624961188653e+02 16 1.0963039832494319e+02 -7.7096357857185083e+01 -2.8842624960984045e+02
17 -1.0862142117090467e+02 7.5521002117836630e+01 2.9024023117219969e+02 17 -1.0862142116718823e+02 7.5521002114629724e+01 2.9024023117422536e+02
18 -1.6565212275867415e-02 -2.7990268691876326e-02 2.6458602067006932e-02 18 -1.6565211729771653e-02 -2.7990264151000276e-02 2.6458599205916818e-02
19 7.1473709241289744e-04 6.1248437925700357e-04 1.7889258733572757e-04 19 7.1473670261989000e-04 6.1248404522910675e-04 1.7889248977386897e-04
20 -7.1473709241289744e-04 -6.1248437925700357e-04 -1.7889258733572757e-04 20 -7.1473670261989000e-04 -6.1248404522910675e-04 -1.7889248977386897e-04
21 -1.1536904971247665e+01 -1.3021625993962397e+01 3.6108894191673429e+01 21 -1.1536904971247079e+01 -1.3021625993961788e+01 3.6108894191671638e+01
22 -1.7333879764643559e+01 -4.2314763344327275e+00 -2.7103019136756011e+01 22 -1.7333879765020438e+01 -4.2314763346057394e+00 -2.7103019136915339e+01
23 2.8870784735891224e+01 1.7253102328395126e+01 -9.0058750549174214e+00 23 2.8870784736267517e+01 1.7253102328567529e+01 -9.0058750547563005e+00
24 6.1437425316795213e+00 -3.4297207023632204e+01 1.8296742414004438e+01 24 6.1437425319896413e+00 -3.4297207021755547e+01 1.8296742415795482e+01
25 -2.4276461284621075e+01 3.8560435260643189e+00 -2.0415720860228767e+01 25 -2.4276461285000558e+01 3.8560435258438814e+00 -2.0415720860437396e+01
26 1.8125049871956215e+01 3.0437790982988908e+01 2.1072387594169975e+00 26 1.8125049872338021e+01 3.0437790983209180e+01 2.1072387596271578e+00
27 8.4078124309265192e+00 -3.6323119973255714e+01 1.4505938075037919e+01 27 8.4078124297937791e+00 -3.6323119973862774e+01 1.4505938075297303e+01
28 -2.8937319168272772e+01 1.2421253477627801e+01 -1.9540501416079319e+01 28 -2.8937319168706676e+01 1.2421253477499173e+01 -1.9540501416287395e+01
29 2.0535244350189391e+01 2.3904950625827414e+01 5.0332497948309163e+00 29 2.0535244350622666e+01 2.3904950625953887e+01 5.0332497950390769e+00
... ...