diff --git a/examples/SPIN/cobalt_fcc/in.spin.cobalt_fcc b/examples/SPIN/cobalt_fcc/in.spin.cobalt_fcc
index c1d7e4f9038b38253236e357d4dca6b3a6b4831f..fd5eeb38b8fd248ecf31085dc7cf03bc77040f97 100644
--- a/examples/SPIN/cobalt_fcc/in.spin.cobalt_fcc
+++ b/examples/SPIN/cobalt_fcc/in.spin.cobalt_fcc
@@ -61,4 +61,3 @@ compute outsp all property/atom spx spy spz sp fmx fmy fmz
dump 100 all custom 1 dump_cobalt_fcc.lammpstrj type x y z c_outsp[1] c_outsp[2] c_outsp[3]
run 1000
-
diff --git a/examples/SPIN/cobalt_hcp/in.spin.cobalt_hcp b/examples/SPIN/cobalt_hcp/in.spin.cobalt_hcp
index e712e511be48a40257ef59d9a762a4dc68e27327..76d41f724d42e157579413c9e056d51799eba4b0 100644
--- a/examples/SPIN/cobalt_hcp/in.spin.cobalt_hcp
+++ b/examples/SPIN/cobalt_hcp/in.spin.cobalt_hcp
@@ -19,7 +19,8 @@ create_atoms 1 box
mass 1 58.93
-set group all spin/random 31 1.72
+#set group all spin/random 31 1.72
+set group all spin 1.72 0.0 0.0 1.0
velocity all create 100 4928459 rot yes dist gaussian
#pair_style hybrid/overlay eam/alloy spin/exchange 4.0 spin/neel 4.0
diff --git a/examples/SPIN/iron/in.spin.iron b/examples/SPIN/iron/in.spin.iron
index 0f40a5af216c66bde4e050e2eb707a01b12ba33a..758cc6c5899fbac5bc2da0bb6fc40a683262409f 100644
--- a/examples/SPIN/iron/in.spin.iron
+++ b/examples/SPIN/iron/in.spin.iron
@@ -53,4 +53,4 @@ thermo 50
compute outsp all property/atom spx spy spz sp fmx fmy fmz
dump 100 all custom 1 dump_iron.lammpstrj type x y z c_outsp[1] c_outsp[2] c_outsp[3]
-run 100000
+run 50000
diff --git a/src/SPIN/atom_vec_spin.cpp b/src/SPIN/atom_vec_spin.cpp
index 369ce8671eb92ce2d699b0e2683b9949d95aa497..22cd78036b79e6c90ae5e18412e61f9d9f5f9209 100644
--- a/src/SPIN/atom_vec_spin.cpp
+++ b/src/SPIN/atom_vec_spin.cpp
@@ -889,12 +889,10 @@ void AtomVecSpin::pack_data(double **buf)
int AtomVecSpin::pack_data_hybrid(int i, double *buf)
{
-
buf[0] = sp[i][0];
buf[1] = sp[i][1];
buf[2] = sp[i][2];
buf[3] = sp[i][3];
-
return 4;
}
@@ -910,7 +908,7 @@ void AtomVecSpin::write_data(FILE *fp, int n, double **buf)
"%-1.16e %d %d %d\n",
(tagint) ubuf(buf[i][0]).i,(int) ubuf(buf[i][1]).i,
buf[i][2],buf[i][3],buf[i][4],
- buf[i][5],buf[i][6],buf[i][7],buf[i][8],buf[i][9],
+ buf[i][5],buf[i][6],buf[i][7],buf[i][8],
(int) ubuf(buf[i][10]).i,(int) ubuf(buf[i][11]).i,
(int) ubuf(buf[i][12]).i);
}
diff --git a/src/SPIN/fix_nve_spin.cpp b/src/SPIN/fix_nve_spin.cpp
index bf45b8cb416af744cff65c4f5f3dc4b46d32c88d..699426de9c22ae242923931546fa961be2f87873 100644
--- a/src/SPIN/fix_nve_spin.cpp
+++ b/src/SPIN/fix_nve_spin.cpp
@@ -69,9 +69,9 @@ enum{NONE};
FixNVESpin::FixNVESpin(LAMMPS *lmp, int narg, char **arg) :
Fix(lmp, narg, arg),
+ pair(NULL), spin_pairs(NULL),
rsec(NULL), stack_head(NULL), stack_foot(NULL),
- backward_stacks(NULL), forward_stacks(NULL),
- pair(NULL), spin_pairs(NULL)
+ backward_stacks(NULL), forward_stacks(NULL)
{
if (lmp->citeme) lmp->citeme->add(cite_fix_nve_spin);
@@ -248,8 +248,8 @@ void FixNVESpin::init()
nlocal_max = atom->nlocal;
stack_head = memory->grow(stack_head,nsectors,"NVE/spin:stack_head");
stack_foot = memory->grow(stack_foot,nsectors,"NVE/spin:stack_foot");
- forward_stacks = memory->grow(forward_stacks,nlocal_max,"NVE/spin:forward_stacks");
backward_stacks = memory->grow(backward_stacks,nlocal_max,"NVE/spin:backward_stacks");
+ forward_stacks = memory->grow(forward_stacks,nlocal_max,"NVE/spin:forward_stacks");
if (nlocal_max == 0)
error->all(FLERR,"Incorrect value of nlocal_max");
@@ -385,8 +385,8 @@ void FixNVESpin::pre_neighbor()
if (nlocal_max < nlocal) { // grow linked lists if necessary
nlocal_max = nlocal;
- forward_stacks = memory->grow(forward_stacks,nlocal_max,"NVE/spin:forward_stacks");
backward_stacks = memory->grow(backward_stacks,nlocal_max,"NVE/spin:backward_stacks");
+ forward_stacks = memory->grow(forward_stacks,nlocal_max,"NVE/spin:forward_stacks");
}
for (int j = 0; j < nsectors; j++) {
diff --git a/src/SPIN/fix_precession_spin.cpp b/src/SPIN/fix_precession_spin.cpp
index 4ee5e87540e89e0e6aae364828f1cb1891743bae..412b9d903b7a04b9ddaf0a9f209326db02d95a19 100644
--- a/src/SPIN/fix_precession_spin.cpp
+++ b/src/SPIN/fix_precession_spin.cpp
@@ -200,7 +200,7 @@ void FixPrecessionSpin::post_force(int vflag)
}
if (aniso_flag) { // compute magnetic anisotropy
- compute_anisotropy(i,spi,fmi);
+ compute_anisotropy(spi,fmi);
emag -= (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
emag *= hbar;
}
@@ -219,7 +219,7 @@ void FixPrecessionSpin::compute_single_precession(int i, double spi[3], double f
compute_zeeman(i,fmi);
}
if (aniso_flag) {
- compute_anisotropy(i,spi,fmi);
+ compute_anisotropy(spi,fmi);
}
}
@@ -236,7 +236,7 @@ void FixPrecessionSpin::compute_zeeman(int i, double fmi[3])
/* ---------------------------------------------------------------------- */
-void FixPrecessionSpin::compute_anisotropy(int i, double spi[3], double fmi[3])
+void FixPrecessionSpin::compute_anisotropy(double spi[3], double fmi[3])
{
double scalar = nax*spi[0] + nay*spi[1] + naz*spi[2];
fmi[0] += scalar*Kax;
diff --git a/src/SPIN/fix_precession_spin.h b/src/SPIN/fix_precession_spin.h
index 2a616b61f0842a673965facc70003d1c354d7c89..9e7bd1a830c943992c4a96c521a21e87c7beff60 100644
--- a/src/SPIN/fix_precession_spin.h
+++ b/src/SPIN/fix_precession_spin.h
@@ -40,7 +40,7 @@ class FixPrecessionSpin : public Fix {
int zeeman_flag, aniso_flag;
void compute_single_precession(int, double *, double *);
void compute_zeeman(int, double *);
- void compute_anisotropy(int, double *, double *);
+ void compute_anisotropy(double *, double *);
protected:
int style; // style of the magnetic precession
diff --git a/src/SPIN/pair_spin_dmi.cpp b/src/SPIN/pair_spin_dmi.cpp
index c89ddd085a1b530ceb4f47fc83dda527645fb8bf..5f735cd1a8900b3dd45e95ae1aafcfd0f71c94ce 100755
--- a/src/SPIN/pair_spin_dmi.cpp
+++ b/src/SPIN/pair_spin_dmi.cpp
@@ -46,7 +46,8 @@ using namespace MathConst;
/* ---------------------------------------------------------------------- */
-PairSpinDmi::PairSpinDmi(LAMMPS *lmp) : PairSpin(lmp)
+PairSpinDmi::PairSpinDmi(LAMMPS *lmp) : PairSpin(lmp),
+lockfixnvespin(NULL)
{
single_enable = 0;
no_virial_fdotr_compute = 1;
@@ -103,9 +104,10 @@ void PairSpinDmi::settings(int narg, char **arg)
void PairSpinDmi::coeff(int narg, char **arg)
{
-
if (!allocated) allocate();
+ // check if args correct
+
if (strcmp(arg[2],"dmi") != 0)
error->all(FLERR,"Incorrect args in pair_style command");
if (narg != 8)
@@ -239,16 +241,15 @@ void PairSpinDmi::compute(int eflag, int vflag)
i = ilist[ii];
itype = type[i];
+ jlist = firstneigh[i];
+ jnum = numneigh[i];
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
-
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
- jlist = firstneigh[i];
- jnum = numneigh[i];
// loop on neighbors
@@ -280,7 +281,7 @@ void PairSpinDmi::compute(int eflag, int vflag)
local_cut2 = cut_spin_dmi[itype][jtype]*cut_spin_dmi[itype][jtype];
- // compute magnetic and mechanical components of soc_dmi
+ // compute dmi interaction
if (rsq <= local_cut2) {
compute_dmi(i,j,eij,fmi,spj);
@@ -296,7 +297,6 @@ void PairSpinDmi::compute(int eflag, int vflag)
fm[i][1] += fmi[1];
fm[i][2] += fmi[2];
- // check newton pair => see if needs correction
if (newton_pair || j < nlocal) {
f[j][0] -= fi[0];
f[j][1] -= fi[1];
diff --git a/src/SPIN/pair_spin_exchange.cpp b/src/SPIN/pair_spin_exchange.cpp
index 36708293f7568c22b143e721f893e6e6b7c171d8..39c5823d0374e434fab44fa43663cdd15e8eaced 100755
--- a/src/SPIN/pair_spin_exchange.cpp
+++ b/src/SPIN/pair_spin_exchange.cpp
@@ -49,7 +49,6 @@ using namespace MathConst;
PairSpinExchange::PairSpinExchange(LAMMPS *lmp) : PairSpin(lmp),
lockfixnvespin(NULL)
{
- hbar = force->hplanck/MY_2PI;
single_enable = 0;
no_virial_fdotr_compute = 1;
lattice_flag = 0;
@@ -203,7 +202,7 @@ void *PairSpinExchange::extract(const char *str, int &dim)
void PairSpinExchange::compute(int eflag, int vflag)
{
int i,j,ii,jj,inum,jnum,itype,jtype;
- double evdwl,ecoul;
+ double evdwl, ecoul;
double xi[3], rij[3], eij[3];
double spi[3], spj[3];
double fi[3], fmi[3];
@@ -233,15 +232,16 @@ void PairSpinExchange::compute(int eflag, int vflag)
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
+ itype = type[i];
+
+ jlist = firstneigh[i];
+ jnum = numneigh[i];
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
- jlist = firstneigh[i];
- jnum = numneigh[i];
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
- itype = type[i];
// loop on neighbors
@@ -293,13 +293,10 @@ void PairSpinExchange::compute(int eflag, int vflag)
}
if (eflag) {
- if (rsq <= local_cut2) {
- evdwl -= (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
- evdwl *= hbar;
- }
+ evdwl -= (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
+ evdwl *= hbar;
} else evdwl = 0.0;
-
if (evflag) ev_tally_xyz(i,j,nlocal,newton_pair,
evdwl,ecoul,fi[0],fi[1],fi[2],rij[0],rij[1],rij[2]);
}
diff --git a/src/SPIN/pair_spin_magelec.cpp b/src/SPIN/pair_spin_magelec.cpp
index a4a17bd1db28b385b6559c3ecdbe0bddd0a9164e..7b6fd323334234938e4a6ad5eefa72fddbc7c34d 100755
--- a/src/SPIN/pair_spin_magelec.cpp
+++ b/src/SPIN/pair_spin_magelec.cpp
@@ -46,9 +46,9 @@ using namespace MathConst;
/* ---------------------------------------------------------------------- */
-PairSpinMagelec::PairSpinMagelec(LAMMPS *lmp) : PairSpin(lmp)
+PairSpinMagelec::PairSpinMagelec(LAMMPS *lmp) : PairSpin(lmp),
+lockfixnvespin(NULL)
{
- hbar = force->hplanck/MY_2PI;
single_enable = 0;
no_virial_fdotr_compute = 1;
lattice_flag = 0;
@@ -103,42 +103,45 @@ void PairSpinMagelec::settings(int narg, char **arg)
void PairSpinMagelec::coeff(int narg, char **arg)
{
- const double hbar = force->hplanck/MY_2PI;
-
if (!allocated) allocate();
- if (strcmp(arg[2],"magelec")==0) {
- if (narg != 8) error->all(FLERR,"Incorrect args in pair_style command");
- int ilo,ihi,jlo,jhi;
- force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi);
- force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi);
-
- const double rij = force->numeric(FLERR,arg[3]);
- const double magelec = (force->numeric(FLERR,arg[4]));
- double mex = force->numeric(FLERR,arg[5]);
- double mey = force->numeric(FLERR,arg[6]);
- double mez = force->numeric(FLERR,arg[7]);
-
- double inorm = 1.0/(mex*mex+mey*mey+mez*mez);
- mex *= inorm;
- mey *= inorm;
- mez *= inorm;
-
- int count = 0;
- for (int i = ilo; i <= ihi; i++) {
- for (int j = MAX(jlo,i); j <= jhi; j++) {
- cut_spin_magelec[i][j] = rij;
- ME[i][j] = magelec/hbar;
- ME_mech[i][j] = magelec;
- v_mex[i][j] = mex;
- v_mey[i][j] = mey;
- v_mez[i][j] = mez;
- setflag[i][j] = 1;
- count++;
- }
+ // check if args correct
+
+ if (strcmp(arg[2],"magelec") != 0)
+ error->all(FLERR,"Incorrect args in pair_style command");
+ if (narg != 8)
+ error->all(FLERR,"Incorrect args in pair_style command");
+
+ int ilo,ihi,jlo,jhi;
+ force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi);
+ force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi);
+
+ const double rij = force->numeric(FLERR,arg[3]);
+ const double magelec = (force->numeric(FLERR,arg[4]));
+ double mex = force->numeric(FLERR,arg[5]);
+ double mey = force->numeric(FLERR,arg[6]);
+ double mez = force->numeric(FLERR,arg[7]);
+
+ double inorm = 1.0/(mex*mex+mey*mey+mez*mez);
+ mex *= inorm;
+ mey *= inorm;
+ mez *= inorm;
+
+ int count = 0;
+ for (int i = ilo; i <= ihi; i++) {
+ for (int j = MAX(jlo,i); j <= jhi; j++) {
+ cut_spin_magelec[i][j] = rij;
+ ME[i][j] = magelec/hbar;
+ ME_mech[i][j] = magelec;
+ v_mex[i][j] = mex;
+ v_mey[i][j] = mey;
+ v_mez[i][j] = mez;
+ setflag[i][j] = 1;
+ count++;
}
- if (count == 0) error->all(FLERR,"Incorrect args in pair_style command");
- } else error->all(FLERR,"Incorrect args in pair_style command");
+ }
+ if (count == 0)
+ error->all(FLERR,"Incorrect args in pair_style command");
}
@@ -210,8 +213,7 @@ void PairSpinMagelec::compute(int eflag, int vflag)
double evdwl, ecoul;
double xi[3], rij[3], eij[3];
double spi[3], spj[3];
- double fi[3], fj[3];
- double fmi[3], fmj[3];
+ double fi[3], fmi[3];
double local_cut2;
double rsq, inorm;
int *ilist,*jlist,*numneigh,**firstneigh;
@@ -240,11 +242,11 @@ void PairSpinMagelec::compute(int eflag, int vflag)
i = ilist[ii];
itype = type[i];
+ jlist = firstneigh[i];
+ jnum = numneigh[i];
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
- jlist = firstneigh[i];
- jnum = numneigh[i];
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
@@ -254,6 +256,7 @@ void PairSpinMagelec::compute(int eflag, int vflag)
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
+ jtype = type[j];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
@@ -262,21 +265,16 @@ void PairSpinMagelec::compute(int eflag, int vflag)
evdwl = 0.0;
fi[0] = fi[1] = fi[2] = 0.0;
- fj[0] = fj[1] = fj[2] = 0.0;
fmi[0] = fmi[1] = fmi[2] = 0.0;
- fmj[0] = fmj[1] = fmj[2] = 0.0;
- rij[0] = rij[1] = rij[2] = 0.0;
rij[0] = x[j][0] - xi[0];
rij[1] = x[j][1] - xi[1];
rij[2] = x[j][2] - xi[2];
rsq = rij[0]*rij[0] + rij[1]*rij[1] + rij[2]*rij[2];
inorm = 1.0/sqrt(rsq);
- eij[0] *= inorm;
- eij[1] *= inorm;
- eij[2] *= inorm;
-
- jtype = type[j];
+ eij[0] = inorm*rij[0];
+ eij[1] = inorm*rij[1];
+ eij[2] = inorm*rij[2];
local_cut2 = cut_spin_magelec[itype][jtype]*cut_spin_magelec[itype][jtype];
@@ -296,11 +294,10 @@ void PairSpinMagelec::compute(int eflag, int vflag)
fm[i][1] += fmi[1];
fm[i][2] += fmi[2];
- // check newton pair => see if needs correction
if (newton_pair || j < nlocal) {
- f[j][0] -= fj[0];
- f[j][1] -= fj[1];
- f[j][2] -= fj[2];
+ f[j][0] -= fi[0];
+ f[j][1] -= fi[1];
+ f[j][2] -= fi[2];
}
if (eflag) {
diff --git a/src/SPIN/pair_spin_neel.cpp b/src/SPIN/pair_spin_neel.cpp
index 408522045fca764ade9366187dbfe3d747aba733..23328736b1e790d4fc9ae73e83e966005b62837d 100755
--- a/src/SPIN/pair_spin_neel.cpp
+++ b/src/SPIN/pair_spin_neel.cpp
@@ -47,7 +47,8 @@ using namespace MathConst;
/* ---------------------------------------------------------------------- */
-PairSpinNeel::PairSpinNeel(LAMMPS *lmp) : PairSpin(lmp)
+PairSpinNeel::PairSpinNeel(LAMMPS *lmp) : PairSpin(lmp),
+lockfixnvespin(NULL)
{
single_enable = 0;
no_virial_fdotr_compute = 1;
@@ -108,8 +109,6 @@ void PairSpinNeel::settings(int narg, char **arg)
void PairSpinNeel::coeff(int narg, char **arg)
{
- const double hbar = force->hplanck/MY_2PI;
-
if (!allocated) allocate();
// check if args correct
@@ -246,11 +245,13 @@ void PairSpinNeel::compute(int eflag, int vflag)
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
+ itype = type[i];
+
+ jlist = firstneigh[i];
+ jnum = numneigh[i];
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
- jlist = firstneigh[i];
- jnum = numneigh[i];
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
@@ -285,7 +286,7 @@ void PairSpinNeel::compute(int eflag, int vflag)
local_cut2 = cut_spin_neel[itype][jtype]*cut_spin_neel[itype][jtype];
- // compute magnetic and mechanical components of neel
+ // compute neel interaction
if (rsq <= local_cut2) {
compute_neel(i,j,rsq,eij,fmi,spi,spj);
@@ -301,7 +302,6 @@ void PairSpinNeel::compute(int eflag, int vflag)
fm[i][1] += fmi[1];
fm[i][2] += fmi[2];
- // check newton pair => see if needs correction
if (newton_pair || j < nlocal) {
f[j][0] -= fi[0];
f[j][1] -= fi[1];