diff --git a/src/fix_adapt.cpp b/src/fix_adapt.cpp
index 65397b189466e57fd774377b56868b2e295225b2..f9a5ba5899e20cfdfec57f4113e8e832afb0a67e 100644
--- a/src/fix_adapt.cpp
+++ b/src/fix_adapt.cpp
@@ -11,6 +11,7 @@
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
+#include "math.h"
#include "string.h"
#include "stdlib.h"
#include "fix_adapt.h"
@@ -148,14 +149,17 @@ void FixAdapt::init()
error->all("Fix adapt pair style does not exist");
pairindex[m] =
pairptr[m]->pre_adapt(param[m],ilo[m],ihi[m],jlo[m],jhi[m]);
- if (pairindex[m] < 0)
- error->all("Fix adapt pair style is not compatible");
+ if (pairindex[m] == -1)
+ error->all("Fix adapt pair parameter is not recognized");
+ if (pairindex[m] == -2)
+ error->all("Fix adapt pair types are not valid");
+
} else if (which[m] == ATOM) {
if (strcmp(param[m],"diameter") == 0) {
awhich[m] = DIAMETER;
if (!atom->radius_flag)
error->all("Fix adapt requires atom attribute radius");
- } else error->all("Invalid fix adapt atom attribute");
+ } else error->all("Fix adapt atom attribute is not recognized");
}
ivar[m] = input->variable->find(var[m]);
@@ -182,14 +186,26 @@ void FixAdapt::pre_force(int vflag)
pairptr[m]->adapt(pairindex[m],ilo[m],ihi[m],jlo[m],jhi[m],value);
else if (which[m] == ATOM) {
- double *radius = atom->radius;
- int *mask = atom->mask;
- int nlocal = atom->nlocal;
+
+ // set radius from diameter
+ // set rmass if both rmass and density are defined
if (awhich[m] == DIAMETER) {
+ int mflag = 0;
+ if (atom->rmass_flag && atom->density_flag) mflag = 1;
+ double PI = 4.0*atan(1.0);
+
+ double *radius = atom->radius;
+ double *rmass = atom->rmass;
+ double *density = atom->density;
+ int *mask = atom->mask;
+ int nlocal = atom->nlocal;
+
for (int i = 0; i < nlocal; i++)
- if (mask[i] & groupbit)
+ if (mask[i] & groupbit) {
radius[i] = 0.5*value;
+ rmass[i] = 4.0*PI/3.0 * radius[i]*radius[i]*radius[i] * density[i];
+ }
}
}
}
diff --git a/src/pair.cpp b/src/pair.cpp
index d72566b4e24841bbf6708848f8dbe137df2e4cdd..e108c712a1aa14a0f06804ec80ae9467b247db81 100644
--- a/src/pair.cpp
+++ b/src/pair.cpp
@@ -23,7 +23,6 @@
#include "stdlib.h"
#include "string.h"
#include "pair.h"
-#include "pair_soft.h"
#include "atom.h"
#include "neighbor.h"
#include "neigh_list.h"
@@ -895,13 +894,6 @@ void Pair::write_file(int narg, char **arg)
force->init();
- // if pair style = soft, set prefactor to final value
-
- Pair *spair = force->pair_match("soft",1);
- if (spair)
- ((PairSoft *) spair)->prefactor[itype][jtype] =
- ((PairSoft *) spair)->prestop[itype][jtype];
-
// if pair style = any of EAM, swap in dummy fp vector
double eamfp[2];
diff --git a/src/pair_hybrid.h b/src/pair_hybrid.h
index 6b66ccd5ec14fcbdd64db718d78cba43d52af3a3..9554d0595d4999a6c367f35bfeab60517bce9d68 100644
--- a/src/pair_hybrid.h
+++ b/src/pair_hybrid.h
@@ -32,7 +32,7 @@ class PairHybrid : public Pair {
char **keywords; // style name of each Pair style
PairHybrid(class LAMMPS *);
- ~PairHybrid();
+ virtual ~PairHybrid();
void compute(int, int);
void settings(int, char **);
virtual void coeff(int, char **);
diff --git a/src/pair_hybrid_overlay.h b/src/pair_hybrid_overlay.h
index 770e0e066ccb6c16b333f0632e3d75827158bd9c..dea84bae4c28f5ddd58f0857a22834e396bb65bd 100644
--- a/src/pair_hybrid_overlay.h
+++ b/src/pair_hybrid_overlay.h
@@ -27,6 +27,7 @@ namespace LAMMPS_NS {
class PairHybridOverlay : public PairHybrid {
public:
PairHybridOverlay(class LAMMPS *);
+ ~PairHybridOverlay() {}
void coeff(int, char **);
private:
diff --git a/src/pair_soft.cpp b/src/pair_soft.cpp
index 92fab52cf6cedcefad868b68be08549956111cc0..8add30e67811bb4dc522adcee4dfc546af749ddd 100644
--- a/src/pair_soft.cpp
+++ b/src/pair_soft.cpp
@@ -14,6 +14,7 @@
#include "math.h"
#include "stdio.h"
#include "stdlib.h"
+#include "string.h"
#include "pair_soft.h"
#include "atom.h"
#include "comm.h"
@@ -43,8 +44,6 @@ PairSoft::~PairSoft()
memory->destroy_2d_int_array(setflag);
memory->destroy_2d_double_array(cutsq);
- memory->destroy_2d_double_array(prestart);
- memory->destroy_2d_double_array(prestop);
memory->destroy_2d_double_array(prefactor);
memory->destroy_2d_double_array(cut);
}
@@ -63,21 +62,6 @@ void PairSoft::compute(int eflag, int vflag)
if (eflag || vflag) ev_setup(eflag,vflag);
else evflag = vflag_fdotr = 0;
- // set current prefactor
- // for minimization, set to prestop
- // for dynamics, ramp from prestart to prestop
- // for 0-step dynamics, set to prestart
-
- double delta = update->ntimestep - update->beginstep;
- if (update->whichflag == 2) delta = 1.0;
- else if (update->nsteps) delta /= update->endstep - update->beginstep;
- else delta = 0.0;
- int ntypes = atom->ntypes;
- for (i = 1; i <= ntypes; i++)
- for (j = 1; j <= ntypes; j++)
- prefactor[i][j] = prestart[i][j] +
- delta * (prestop[i][j] - prestart[i][j]);
-
double **x = atom->x;
double **f = atom->f;
int *type = atom->type;
@@ -161,20 +145,8 @@ void PairSoft::allocate()
cutsq = memory->create_2d_double_array(n+1,n+1,"pair:cutsq");
- prestart = memory->create_2d_double_array(n+1,n+1,"pair:prestart");
- prestop = memory->create_2d_double_array(n+1,n+1,"pair:prestop");
prefactor = memory->create_2d_double_array(n+1,n+1,"pair:prefactor");
cut = memory->create_2d_double_array(n+1,n+1,"pair:cut");
-
- // init prestart and prestop to 0.0
- // since pair_hybrid can use all types even if pair_soft sub-class
- // never sets them
-
- for (int i = 1; i <= n; i++)
- for (int j = i; j <= n; j++) {
- prestart[i][j] = 0.0;
- prestop[i][j] = 0.0;
- }
}
/* ----------------------------------------------------------------------
@@ -203,24 +175,22 @@ void PairSoft::settings(int narg, char **arg)
void PairSoft::coeff(int narg, char **arg)
{
- if (narg < 4 || narg > 5) error->all("Incorrect args for pair coefficients");
+ if (narg < 3 || narg > 4) error->all("Incorrect args for pair coefficients");
if (!allocated) allocate();
int ilo,ihi,jlo,jhi;
force->bounds(arg[0],atom->ntypes,ilo,ihi);
force->bounds(arg[1],atom->ntypes,jlo,jhi);
- double prestart_one = force->numeric(arg[2]);
- double prestop_one = force->numeric(arg[3]);
+ double prefactor_one = force->numeric(arg[2]);
double cut_one = cut_global;
- if (narg == 5) cut_one = force->numeric(arg[4]);
+ if (narg == 4) cut_one = force->numeric(arg[3]);
int count = 0;
for (int i = ilo; i <= ihi; i++) {
for (int j = MAX(jlo,i); j <= jhi; j++) {
- prestart[i][j] = prestart_one;
- prestop[i][j] = prestop_one;
+ prefactor[i][j] = prefactor_one;
cut[i][j] = cut_one;
setflag[i][j] = 1;
count++;
@@ -239,13 +209,11 @@ double PairSoft::init_one(int i, int j)
// always mix prefactors geometrically
if (setflag[i][j] == 0) {
- prestart[i][j] = sqrt(prestart[i][i]*prestart[j][j]);
- prestop[i][j] = sqrt(prestop[i][i]*prestop[j][j]);
+ prefactor[i][j] = sqrt(prefactor[i][i]*prefactor[j][j]);
cut[i][j] = mix_distance(cut[i][i],cut[j][j]);
}
- prestart[j][i] = prestart[i][j];
- prestop[j][i] = prestop[i][j];
+ prefactor[j][i] = prefactor[i][j];
cut[j][i] = cut[i][j];
return cut[i][j];
@@ -264,8 +232,7 @@ void PairSoft::write_restart(FILE *fp)
for (j = i; j <= atom->ntypes; j++) {
fwrite(&setflag[i][j],sizeof(int),1,fp);
if (setflag[i][j]) {
- fwrite(&prestart[i][j],sizeof(double),1,fp);
- fwrite(&prestop[i][j],sizeof(double),1,fp);
+ fwrite(&prefactor[i][j],sizeof(double),1,fp);
fwrite(&cut[i][j],sizeof(double),1,fp);
}
}
@@ -289,12 +256,10 @@ void PairSoft::read_restart(FILE *fp)
MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world);
if (setflag[i][j]) {
if (me == 0) {
- fread(&prestart[i][j],sizeof(double),1,fp);
- fread(&prestop[i][j],sizeof(double),1,fp);
+ fread(&prefactor[i][j],sizeof(double),1,fp);
fread(&cut[i][j],sizeof(double),1,fp);
}
- MPI_Bcast(&prestart[i][j],1,MPI_DOUBLE,0,world);
- MPI_Bcast(&prestop[i][j],1,MPI_DOUBLE,0,world);
+ MPI_Bcast(&prefactor[i][j],1,MPI_DOUBLE,0,world);
MPI_Bcast(&cut[i][j],1,MPI_DOUBLE,0,world);
}
}
@@ -324,6 +289,38 @@ void PairSoft::read_restart_settings(FILE *fp)
MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
}
+/* ----------------------------------------------------------------------
+ check if name is recognized, return integer index for that name
+ if name not recognized, return -1
+ if type pair setting, return -2 if no type pairs are set
+------------------------------------------------------------------------- */
+
+int PairSoft::pre_adapt(char *name, int ilo, int ihi, int jlo, int jhi)
+{
+ int count = 0;
+ for (int i = ilo; i <= ihi; i++)
+ for (int j = MAX(jlo,i); j <= jhi; j++)
+ count++;
+ if (count == 0) return -2;
+
+ if (strcmp(name,"a") == 0) return 0;
+ return -1;
+}
+
+/* ----------------------------------------------------------------------
+ adapt parameter indexed by which
+ change all pair variables affected by the reset parameter
+ if type pair setting, set I-J and J-I coeffs
+------------------------------------------------------------------------- */
+
+void PairSoft::adapt(int which, int ilo, int ihi, int jlo, int jhi,
+ double value)
+{
+ for (int i = ilo; i <= ihi; i++)
+ for (int j = MAX(jlo,i); j <= jhi; j++)
+ prefactor[i][j] = prefactor[j][i] = value;
+}
+
/* ---------------------------------------------------------------------- */
double PairSoft::single(int i, int j, int itype, int jtype, double rsq,
diff --git a/src/pair_soft.h b/src/pair_soft.h
index 026c351f6403be0d554455d4002ceddadbbd4a99..eceb98ec6f6451a805b2da5a37244a4c583c1409 100644
--- a/src/pair_soft.h
+++ b/src/pair_soft.h
@@ -38,12 +38,13 @@ class PairSoft : public Pair {
void read_restart(FILE *);
void write_restart_settings(FILE *);
void read_restart_settings(FILE *);
+ int pre_adapt(char *, int, int, int, int);
+ void adapt(int, int, int, int, int, double);
double single(int, int, int, int, double, double, double, double &);
private:
double PI;
double cut_global;
- double **prestart,**prestop;
double **prefactor;
double **cut;