diff --git a/doc/src/fix_neb.txt b/doc/src/fix_neb.txt index 7d57aebfd4..52d8a7df84 100644 --- a/doc/src/fix_neb.txt +++ b/doc/src/fix_neb.txt @@ -14,10 +14,10 @@ fix ID group-ID neb Kspring keyword value :pre ID, group-ID are documented in "fix"_fix.html command :ulb,l neb = style name of this fix command :l -Kspring = parallel spring constant (force/distance units or force units, see nudge keyword) :l +Kspring = spring constant for parallel nudging force (force/distance units or force units, see parallel keyword) :l zero or more keyword/value pairs may be appended :l -keyword = {nudge} or {perp} or {ends} :l - {nudge} value = {neigh} or {ideal} +keyword = {parallel} or {perp} or {end} :l + {parallel} value = {neigh} or {ideal} {neigh} = parallel nudging force based on distance to neighbor replicas (Kspring = force/distance units) {ideal} = parallel nudging force based on interpolated ideal position (Kspring = force units) {perp} value = {Kspring2} @@ -28,7 +28,7 @@ keyword = {nudge} or {perp} or {ends} :l {last} = apply force to last replica {last/efirst} = apply force to last replica and set its target energy to that of first replica {last/efirst/middle} = same as {last/efirst} plus prevent middle replicas having lower energy than first replica - {Kspring3} = spring constant for target energy term (1/distance units) :pre + {Kspring3} = spring constant for target energy term (1/distance units) :pre,ule [Examples:] @@ -59,22 +59,22 @@ interatomic force Fi = -Grad(V) for each replica I is altered. For all intermediate replicas (i.e. for 1 < I < N, except the climbing replica) the force vector becomes: -Fi = -Grad(V) + (Grad(V) dot T') T' + Fnudge_parallel + Fspring_perp :pre +Fi = -Grad(V) + (Grad(V) dot T') T' + Fnudge_parallel + Fnudge_perp :pre T' is the unit "tangent" vector for replica I and is a function of Ri, Ri-1, Ri+1, and the potential energy of the 3 replicas; it points roughly in the direction of (Ri+i - Ri-1); see the -"(Henkelman1)"_#Henkelman1 paper for details. Ri gives the atomic +"(Henkelman1)"_#Henkelman1 paper for details. Ri are the atomic coordinates of replica I; Ri-1 and Ri+1 are the coordinates of its neighbor replicas. The term (Grad(V) dot T') is used to remove the component of the gradient parallel to the path which would tend to distribute the replica unevenly along the path. Fnudge_parallel is an artificial nudging force which is applied only in the tangent direction and which maintains the equal spacing between replicas (see -below for more information). Fspring_perp is an optional artificial -spring which is applied only perpendicular to the tangent and which -prevent the paths from forming acute kinks (see below for more -information). +below for more information). Fnudge_perp is an optional artificial +spring which is applied in a direction perpendicular to the tangent +direction and which prevent the paths from forming acute kinks (see +below for more information). In the second stage of the NEB calculation, the interatomic force Fi for the climbing replica (the replica of highest energy after the @@ -86,7 +86,7 @@ and the relaxation procedure is continued to a new converged MEP. :line -The keyword {nudge} specifies how the parallel nudging force is +The keyword {parallel} specifies how the parallel nudging force is computed. With a value of {neigh}, the parallel nudging force is computed as in "(Henkelman1)"_#Henkelman1 by connecting each intermediate replica with the previous and the next image: @@ -98,7 +98,7 @@ units. With a value of {ideal}, the spring force is computed as suggested in "(WeinenE)"_#WeinenE : - + Fnudge_parallel = -{Kspring} * (RD-RDideal) / (2 * meanDist) :pre where RD is the "reaction coordinate" see "neb"_neb.html section, and @@ -113,15 +113,16 @@ keeping the replicas equally spaced. :line -The keyword {perp} adds a spring force perpendicular to the path in -order to prevent the path from becoming too kinky. It -can significantly improve the convergence of the NEB calculation when -the resolution is poor. I.e. when too few replicas are used; see +The keyword {perp} specifies if and how a perpendicual nudging force +is computed. It adds a spring force perpendicular to the path in +order to prevent the path from becoming too kinky. It can +significantly improve the convergence of the NEB calculation when the +resolution is poor. I.e. when few replicas are used; see "(Maras)"_#Maras1 for details. The perpendicular spring force is given by -Fspring_perp = {Kspring2} * F(Ri-1,Ri,Ri+1) (Ri+1 + Ri-1 - 2 Ri) :pre +Fnudge_perp = {Kspring2} * F(Ri-1,Ri,Ri+1) (Ri+1 + Ri-1 - 2 Ri) :pre where {Kspring2} is the specified value. F(Ri-1 Ri R+1) is a smooth scalar function of the angle Ri-1 Ri Ri+1. It is equal to 0.0 when @@ -133,11 +134,11 @@ force is added. :line -By default, no nudging forces act on the first and last replicas during -the NEB relaxation, so these replicas simply relax toward their -respective local minima. By using the key word {end}, additional forces -can be applied to the first or last replica, to enable them to relax -toward a MEP while constraining their energy. +By default, no additional forces act on the first and last replicas +during the NEB relaxation, so these replicas simply relax toward their +respective local minima. By using the key word {end}, additional +forces can be applied to the first and/or last replicas, to enable +them to relax toward a MEP while constraining their energy. The interatomic force Fi for the specified replica becomes: @@ -177,9 +178,10 @@ only be done if a particular intermediate replica has a lower energy than the first replica. This should effectively prevent the intermediate replicas from over-relaxing. -After converging a NEB calculation using an {estyle} of {last/efirst/middle}, -you should check that all intermediate replicas have a larger energy than the -first replica. If this is not the case, the path is probably not a MEP. +After converging a NEB calculation using an {estyle} of +{last/efirst/middle}, you should check that all intermediate replicas +have a larger energy than the first replica. If this is not the case, +the path is probably not a MEP. Finally, note that if the last replica converges toward a local minimum which has a larger energy than the energy of the first diff --git a/examples/neb/README b/examples/neb/README index d01ca35966..5ef32f2ba6 100644 --- a/examples/neb/README +++ b/examples/neb/README @@ -8,9 +8,7 @@ mpirun -np 3 lmp_g++ -partition 3x1 -in in.neb.sivac mpirun -np 8 lmp_g++ -partition 4x2 -in in.neb.hop1 mpirun -np 8 lmp_g++ -partition 4x2 -in in.neb.hop2 mpirun -np 8 lmp_g++ -partition 4x2 -in in.neb.hop1.end -mpirun -np 6 lmp_g++ -partition 3x2 -in in.neb.sivac -mpirun -np 9 lmp_g++ -partition 3x3 -in in.neb.sivac - +mpirun -np 8 lmp_g++ -partition 4x2 -in in.neb.sivac Note that more than 4 replicas should be used for a precise estimate of the activation energy corresponding to a transition. diff --git a/examples/neb/in.neb.hop1 b/examples/neb/in.neb.hop1 index 9b5dcb07ec..f26b52a28a 100644 --- a/examples/neb/in.neb.hop1 +++ b/examples/neb/in.neb.hop1 @@ -51,7 +51,7 @@ set group nebatoms type 3 group nonneb subtract all nebatoms fix 1 lower setforce 0.0 0.0 0.0 -fix 2 nebatoms neb 1.0 #nudge ideal +fix 2 nebatoms neb 1.0 parallel ideal fix 3 all enforce2d thermo 100 diff --git a/examples/neb/in.neb.hop1.end b/examples/neb/in.neb.hop1.end index 2f4ba526d8..81e5315306 100644 --- a/examples/neb/in.neb.hop1.end +++ b/examples/neb/in.neb.hop1.end @@ -41,7 +41,7 @@ set group nebatoms type 3 group nonneb subtract all nebatoms fix 1 lower setforce 0.0 0.0 0.0 -fix 2 nebatoms neb 1.0 nudge ideal end first 1.0 +fix 2 nebatoms neb 1.0 parallel ideal end first 1.0 fix 3 all enforce2d thermo 100 diff --git a/examples/neb/log.19June17.neb.hop1.end.g++.4 b/examples/neb/log.19Jun17.neb.hop1.end.g++.4 similarity index 100% rename from examples/neb/log.19June17.neb.hop1.end.g++.4 rename to examples/neb/log.19Jun17.neb.hop1.end.g++.4 diff --git a/examples/neb/log.19June17.neb.hop1.end.g++.8 b/examples/neb/log.19Jun17.neb.hop1.end.g++.8 similarity index 100% rename from examples/neb/log.19June17.neb.hop1.end.g++.8 rename to examples/neb/log.19Jun17.neb.hop1.end.g++.8 diff --git a/examples/neb/log.19Jun17.neb.hop1.g++.4 b/examples/neb/log.19Jun17.neb.hop1.g++.4 new file mode 100644 index 0000000000..e2984c031c --- /dev/null +++ b/examples/neb/log.19Jun17.neb.hop1.g++.4 @@ -0,0 +1,9 @@ +LAMMPS (19 May 2017) +Running on 4 partitions of processors +Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN +0 4327.2753 2746.3378 0.082169072 4.9967651 4514.5424 0.42933428 0.42323635 1.8941131 0 -3.0535948 0.33333333 -2.6242605 0.66666667 -2.7623811 1 -3.0474969 +87 0.095951502 0.052720903 0.005588927 0.065110105 0.12467831 0.0071014928 0.0022798007 2.3003372 0 -3.0535967 0.32435271 -3.0473127 0.62805027 -3.0464952 1 -3.048775 +Climbing replica = 3 +Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN +87 0.14137277 0.11108954 0.005588927 0.065110105 0.12467831 0.0071014928 0.0022798007 2.3003372 0 -3.0535967 0.32435271 -3.0473127 0.62805027 -3.0464952 1 -3.048775 +124 0.099583263 0.085936899 0.0044220372 0.023873795 0.091308308 0.0071061754 0.0022863931 2.308121 0 -3.0535968 0.32223905 -3.0473329 0.61673898 -3.0464906 1 -3.048777 diff --git a/examples/neb/log.19Jun17.neb.hop1.g++.8 b/examples/neb/log.19Jun17.neb.hop1.g++.8 new file mode 100644 index 0000000000..d1be1284fa --- /dev/null +++ b/examples/neb/log.19Jun17.neb.hop1.g++.8 @@ -0,0 +1,9 @@ +LAMMPS (19 May 2017) +Running on 4 partitions of processors +Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN +0 4327.2753 2746.3378 0.082169072 4.9967651 4514.5424 0.42933428 0.42323635 1.8941131 0 -3.0535948 0.33333333 -2.6242605 0.66666667 -2.7623811 1 -3.0474969 +87 0.095951792 0.052720902 0.0055889267 0.065110091 0.12467831 0.0071014928 0.0022798007 2.3003372 0 -3.0535967 0.32435271 -3.0473127 0.62805027 -3.0464952 1 -3.048775 +Climbing replica = 3 +Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN +87 0.14137297 0.11108954 0.0055889267 0.065110091 0.12467831 0.0071014928 0.0022798007 2.3003372 0 -3.0535967 0.32435271 -3.0473127 0.62805027 -3.0464952 1 -3.048775 +124 0.099582186 0.08593683 0.0044220345 0.023873731 0.091308197 0.0071061754 0.0022863931 2.3081211 0 -3.0535968 0.32223904 -3.0473329 0.61673896 -3.0464906 1 -3.048777 diff --git a/examples/neb/log.19June17.neb.hop2.g++.4 b/examples/neb/log.19Jun17.neb.hop2.g++.4 similarity index 100% rename from examples/neb/log.19June17.neb.hop2.g++.4 rename to examples/neb/log.19Jun17.neb.hop2.g++.4 diff --git a/examples/neb/log.19June17.neb.hop2.g++.8 b/examples/neb/log.19Jun17.neb.hop2.g++.8 similarity index 100% rename from examples/neb/log.19June17.neb.hop2.g++.8 rename to examples/neb/log.19Jun17.neb.hop2.g++.8 diff --git a/examples/neb/log.19June17.neb.sivac.g++.4 b/examples/neb/log.19Jun17.neb.sivac.g++.4 similarity index 100% rename from examples/neb/log.19June17.neb.sivac.g++.4 rename to examples/neb/log.19Jun17.neb.sivac.g++.4 diff --git a/examples/neb/log.19June17.neb.sivac.g++.8 b/examples/neb/log.19Jun17.neb.sivac.g++.8 similarity index 100% rename from examples/neb/log.19June17.neb.sivac.g++.8 rename to examples/neb/log.19Jun17.neb.sivac.g++.8 diff --git a/examples/neb/log.19June17.neb.hop1.g++.4 b/examples/neb/log.19June17.neb.hop1.g++.4 deleted file mode 100644 index ee1596cc0b..0000000000 --- a/examples/neb/log.19June17.neb.hop1.g++.4 +++ /dev/null @@ -1,10 +0,0 @@ -LAMMPS (19 May 2017) -Running on 4 partitions of processors -Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN -0 4327.2753 2746.3378 0.082169072 4.9967651 4514.5424 0.42933428 0.42323635 1.8941131 0 -3.0535948 0.33333333 -2.6242605 0.66666667 -2.7623811 1 -3.0474969 -100 0.10482184 0.085218486 0.0051952047 0.04785954 0.19041553 0.0070900402 0.0022691875 2.3031875 0 -3.0535967 0.31839181 -3.0473647 0.63987598 -3.0465067 1 -3.0487759 -111 0.096708467 0.07803707 0.0048656875 0.03613038 0.19671332 0.0070871172 0.0022668002 2.3052945 0 -3.0535968 0.31853431 -3.0473633 0.64178871 -3.0465096 1 -3.0487764 -Climbing replica = 3 -Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN -111 0.2023467 0.1777038 0.0048656875 0.03613038 0.19671332 0.0070871172 0.0022668002 2.3052945 0 -3.0535968 0.31853431 -3.0473633 0.64178871 -3.0465096 1 -3.0487764 -179 0.096874474 0.090676856 0.0034851031 0.0094134782 0.093630619 0.0071047642 0.0022856172 2.3122768 0 -3.0535969 0.31577311 -3.0473955 0.61798541 -3.0464922 1 -3.0487778 diff --git a/examples/neb/log.19June17.neb.hop1.g++.8 b/examples/neb/log.19June17.neb.hop1.g++.8 deleted file mode 100644 index 325ce73c14..0000000000 --- a/examples/neb/log.19June17.neb.hop1.g++.8 +++ /dev/null @@ -1,10 +0,0 @@ -LAMMPS (19 May 2017) -Running on 4 partitions of processors -Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN -0 4327.2753 2746.3378 0.082169072 4.9967651 4514.5424 0.42933428 0.42323635 1.8941131 0 -3.0535948 0.33333333 -2.6242605 0.66666667 -2.7623811 1 -3.0474969 -100 0.10482171 0.085218406 0.0051952008 0.047859379 0.1904156 0.0070900401 0.0022691875 2.3031875 0 -3.0535967 0.31839181 -3.0473647 0.639876 -3.0465067 1 -3.0487759 -111 0.096708718 0.078036984 0.0048656841 0.036130268 0.1967134 0.0070871172 0.0022668002 2.3052946 0 -3.0535968 0.31853431 -3.0473633 0.64178873 -3.0465096 1 -3.0487764 -Climbing replica = 3 -Step MaxReplicaForce MaxAtomForce GradV0 GradV1 GradVc EBF EBR RDT RD1 PE1 RD2 PE2 ... RDN PEN -111 0.20234693 0.17770387 0.0048656841 0.036130268 0.1967134 0.0070871172 0.0022668002 2.3052946 0 -3.0535968 0.31853431 -3.0473633 0.64178873 -3.0465096 1 -3.0487764 -178 0.09975409 0.093814031 0.0035463662 0.010006594 0.096949208 0.0071042931 0.0022851195 2.312004 0 -3.0535969 0.31607934 -3.0473923 0.618931 -3.0464926 1 -3.0487777 diff --git a/src/REPLICA/fix_neb.cpp b/src/REPLICA/fix_neb.cpp index bd21a4d2ca..6daaf94710 100644 --- a/src/REPLICA/fix_neb.cpp +++ b/src/REPLICA/fix_neb.cpp @@ -49,8 +49,6 @@ FixNEB::FixNEB(LAMMPS *lmp, int narg, char **arg) : displacements(NULL) { - - if (narg < 4) error->all(FLERR,"Illegal fix neb command"); kspring = force->numeric(FLERR,arg[3]); @@ -68,7 +66,7 @@ FixNEB::FixNEB(LAMMPS *lmp, int narg, char **arg) : int iarg = 4; while (iarg < narg) { - if (strcmp(arg[iarg],"nudge") == 0) { + if (strcmp(arg[iarg],"parallel") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal fix neb command"); if (strcmp(arg[iarg+1],"ideal") == 0) { NEBLongRange = true; diff --git a/src/dump.cpp b/src/dump.cpp index dde89024ff..44098298ba 100644 --- a/src/dump.cpp +++ b/src/dump.cpp @@ -30,7 +30,7 @@ using namespace LAMMPS_NS; -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) // allocate space for static class variable Dump *Dump::dumpptr; #else @@ -692,7 +692,7 @@ void Dump::sort() index[idsort[i]-idlo] = i; } -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) if (!reorderflag) { dumpptr = this; for (i = 0; i < nme; i++) index[i] = i; @@ -728,7 +728,8 @@ void Dump::sort() memcpy(&buf[i*size_one],&bufsort[index[i]*size_one],nbytes); } -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) + /* ---------------------------------------------------------------------- compare two atom IDs called via qsort() in sort() method @@ -795,6 +796,7 @@ int Dump::bufcompare_reverse(const void *pi, const void *pj) compare two atom IDs called via merge_sort() in sort() method ------------------------------------------------------------------------- */ + int Dump::idcompare(const int i, const int j, void *ptr) { tagint *idsort = ((Dump *)ptr)->idsort; diff --git a/src/irregular.cpp b/src/irregular.cpp index 74fdb6f772..6cd1b22c2f 100644 --- a/src/irregular.cpp +++ b/src/irregular.cpp @@ -23,17 +23,13 @@ using namespace LAMMPS_NS; -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) // allocate space for static class variable // prototype for non-class function - int *Irregular::proc_recv_copy; static int compare_standalone(const void *, const void *); - #else - #include "mergesort.h" - // prototype for non-class function static int compare_standalone(const int, const int, void *); #endif @@ -433,7 +429,7 @@ int Irregular::create_atom(int n, int *sizes, int *proclist, int sortflag) for (i = 0; i < nrecv_proc; i++) order[i] = i; -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) proc_recv_copy = proc_recv; qsort(order,nrecv_proc,sizeof(int),compare_standalone); #else @@ -464,7 +460,7 @@ int Irregular::create_atom(int n, int *sizes, int *proclist, int sortflag) return nrecvsize; } -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) /* ---------------------------------------------------------------------- comparison function invoked by qsort() @@ -705,12 +701,13 @@ int Irregular::create_data(int n, int *proclist, int sortflag) for (i = 0; i < nrecv_proc; i++) order[i] = i; -#if defined(LMP_USE_LIBC_QSORT) +#if defined(LMP_QSORT) proc_recv_copy = proc_recv; qsort(order,nrecv_proc,sizeof(int),compare_standalone); #else merge_sort(order,nrecv_proc,(void *)proc_recv,compare_standalone); #endif + int j; for (i = 0; i < nrecv_proc; i++) { j = order[i];