jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Commit JT 030419 

- correction of pair_spin calculations
- corrects an error between i and ii lists in single/pair calc.
This commit is contained in:
julient31 2019-03-04 08:04:11 -07:00
parent 849e52040a
commit ab0c35be93
4 changed files with 228 additions and 134 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

104
src/SPIN/pair_spin_dmi.cpp
View File

@ -329,7 +329,9 @@ void PairSpinDmi::compute(int eflag, int vflag)
}
/* ---------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
update the pair interactions fmi acting on the spin ii
------------------------------------------------------------------------- */
void PairSpinDmi::compute_single_pair(int ii, double fmi[3])
{
@ -341,52 +343,76 @@ void PairSpinDmi::compute_single_pair(int ii, double fmi[3])
double delx,dely,delz;
double spj[3];
int i,j,jnum,itype,jtype;
int *ilist,*jlist,*numneigh,**firstneigh;
int i,j,jnum,itype,jtype,ntypes;
int k,locflag;
int *jlist,*numneigh,**firstneigh;
double rsq, inorm;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// check if interaction applies to type of ii
i = ilist[ii];
itype = type[i];
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
jlist = firstneigh[i];
jnum = numneigh[i];
for (int jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
spj[2] = sp[j][2];
delx = xi[0] - x[j][0];
dely = xi[1] - x[j][1];
delz = xi[2] - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
inorm = 1.0/sqrt(rsq);
eij[0] = -inorm*delx;
eij[1] = -inorm*dely;
eij[2] = -inorm*delz;
local_cut2 = cut_spin_dmi[itype][jtype]*cut_spin_dmi[itype][jtype];
if (rsq <= local_cut2) {
compute_dmi(i,j,eij,fmi,spj);
}
itype = type[ii];
ntypes = atom->ntypes;
locflag = 0;
k = 1;
while (k <= ntypes) {
if (k <= itype) {
if (setflag[k][itype] == 1) {
locflag =1;
break;
}
k++;
} else if (k > itype) {
if (setflag[itype][k] == 1) {
locflag =1;
break;
}
k++;
} else error->all(FLERR,"Wrong type number");
}
// if interaction applies to type ii,
// locflag = 1 and compute pair interaction
//i = ilist[ii];
if (locflag == 1) {
xi[0] = x[ii][0];
xi[1] = x[ii][1];
xi[2] = x[ii][2];
jlist = firstneigh[ii];
jnum = numneigh[ii];
for (int jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
spj[2] = sp[j][2];
delx = xi[0] - x[j][0];
dely = xi[1] - x[j][1];
delz = xi[2] - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
inorm = 1.0/sqrt(rsq);
eij[0] = -inorm*delx;
eij[1] = -inorm*dely;
eij[2] = -inorm*delz;
local_cut2 = cut_spin_dmi[itype][jtype]*cut_spin_dmi[itype][jtype];
if (rsq <= local_cut2) {
compute_dmi(ii,j,eij,fmi,spj);
}
}
}
}
/* ----------------------------------------------------------------------

88
src/SPIN/pair_spin_exchange.cpp
View File

@ -318,7 +318,6 @@ void PairSpinExchange::compute(int eflag, int vflag)
void PairSpinExchange::compute_single_pair(int ii, double fmi[3])
{
int *type = atom->type;
double **x = atom->x;
double **sp = atom->sp;
@ -327,46 +326,70 @@ void PairSpinExchange::compute_single_pair(int ii, double fmi[3])
double delx,dely,delz;
double spj[3];
int i,j,jnum,itype,jtype;
int *ilist,*jlist,*numneigh,**firstneigh;
int i,j,jnum,itype,jtype,ntypes;
int k,locflag;
int *jlist,*numneigh,**firstneigh;
double rsq;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
i = ilist[ii];
itype = type[i];
// check if interaction applies to type of ii
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
jlist = firstneigh[i];
jnum = numneigh[i];
for (int jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
local_cut2 = cut_spin_exchange[itype][jtype]*cut_spin_exchange[itype][jtype];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
spj[2] = sp[j][2];
delx = xi[0] - x[j][0];
dely = xi[1] - x[j][1];
delz = xi[2] - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
if (rsq <= local_cut2) {
compute_exchange(i,j,rsq,fmi,spj);
}
itype = type[ii];
ntypes = atom->ntypes;
locflag = 0;
k = 1;
while (k <= ntypes) {
if (k <= itype) {
if (setflag[k][itype] == 1) {
locflag =1;
break;
}
k++;
} else if (k > itype) {
if (setflag[itype][k] == 1) {
locflag =1;
break;
}
k++;
} else error->all(FLERR,"Wrong type number");
}
// if interaction applies to type ii,
// locflag = 1 and compute pair interaction
if (locflag == 1) {
xi[0] = x[ii][0];
xi[1] = x[ii][1];
xi[2] = x[ii][2];
jlist = firstneigh[ii];
jnum = numneigh[ii];
for (int jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
local_cut2 = cut_spin_exchange[itype][jtype]*cut_spin_exchange[itype][jtype];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
spj[2] = sp[j][2];
delx = xi[0] - x[j][0];
dely = xi[1] - x[j][1];
delz = xi[2] - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
if (rsq <= local_cut2) {
compute_exchange(ii,j,rsq,fmi,spj);
}
}
}
}
/* ----------------------------------------------------------------------
@ -528,3 +551,4 @@ void PairSpinExchange::read_restart_settings(FILE *fp)
MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
}

56
src/SPIN/pair_spin_magelec.cpp
View File

@ -319,7 +319,9 @@ void PairSpinMagelec::compute(int eflag, int vflag)
}
/* ---------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
update the pair interactions fmi acting on the spin ii
------------------------------------------------------------------------- */
void PairSpinMagelec::compute_single_pair(int ii, double fmi[3])
{
@ -331,30 +333,48 @@ void PairSpinMagelec::compute_single_pair(int ii, double fmi[3])
double delx,dely,delz;
double spj[3];
int i,j,inum,jnum,itype,jtype;
int *ilist,*jlist,*numneigh,**firstneigh;
int i,j,jnum,itype,jtype,ntypes;
int k,locflag;
int *jlist,*numneigh,**firstneigh;
double rsq, inorm;
inum = list->inum;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// check if interaction applies to type of ii
// compute pair if
itype = type[ii];
ntypes = atom->ntypes;
locflag = 0;
k = 1;
while (k <= ntypes) {
if (k <= itype) {
if (setflag[k][itype] == 1) {
locflag =1;
break;
}
k++;
} else if (k > itype) {
if (setflag[itype][k] == 1) {
locflag =1;
break;
}
k++;
} else error->all(FLERR,"Wrong type number");
}
if (ii < inum) {
// if interaction applies to type ii,
// locflag = 1 and compute pair interaction
i = ilist[ii];
itype = type[i];
if (locflag == 1) {
xi[0] = xi[1] = xi[2] = 0.0;
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
xi[0] = x[ii][0];
xi[1] = x[ii][1];
xi[2] = x[ii][2];
jlist = firstneigh[i];
jnum = numneigh[i];
jlist = firstneigh[ii];
jnum = numneigh[ii];
for (int jj = 0; jj < jnum; jj++) {
@ -377,12 +397,10 @@ void PairSpinMagelec::compute_single_pair(int ii, double fmi[3])
eij[2] = -inorm*delz;
if (rsq <= local_cut2) {
compute_magelec(i,j,eij,fmi,spj);
compute_magelec(ii,j,eij,fmi,spj);
}
}
}
}
}
/* ---------------------------------------------------------------------- */

114
src/SPIN/pair_spin_neel.cpp
View File

@ -330,7 +330,9 @@ void PairSpinNeel::compute(int eflag, int vflag)
}
/* ---------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
update the pair interactions fmi acting on the spin ii
------------------------------------------------------------------------- */
void PairSpinNeel::compute_single_pair(int ii, double fmi[3])
{
@ -342,57 +344,81 @@ void PairSpinNeel::compute_single_pair(int ii, double fmi[3])
double xi[3], rij[3], eij[3];
double spi[3], spj[3];
int i,j,jnum,itype,jtype;
int *ilist,*jlist,*numneigh,**firstneigh;
int i,j,jnum,itype,jtype,ntypes;
int k,locflag;
int *jlist,*numneigh,**firstneigh;
double rsq, inorm;
ilist = list->ilist;
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// check if interaction applies to type of ii
i = ilist[ii];
itype = type[i];
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
xi[0] = x[i][0];
xi[1] = x[i][1];
xi[2] = x[i][2];
eij[0] = eij[1] = eij[2] = 0.0;
jlist = firstneigh[i];
jnum = numneigh[i];
for (int jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
local_cut2 = cut_spin_neel[itype][jtype]*cut_spin_neel[itype][jtype];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
spj[2] = sp[j][2];
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*rij[0];
eij[1] = inorm*rij[1];
eij[2] = inorm*rij[2];
if (rsq <= local_cut2) {
compute_neel(i,j,rsq,eij,fmi,spi,spj);
}
itype = type[ii];
ntypes = atom->ntypes;
locflag = 0;
k = 1;
while (k <= ntypes) {
if (k <= itype) {
if (setflag[k][itype] == 1) {
locflag =1;
break;
}
k++;
} else if (k > itype) {
if (setflag[itype][k] == 1) {
locflag =1;
break;
}
k++;
} else error->all(FLERR,"Wrong type number");
}
// if interaction applies to type ii,
// locflag = 1 and compute pair interaction
if (locflag == 1) {
spi[0] = sp[ii][0];
spi[1] = sp[ii][1];
spi[2] = sp[ii][2];
xi[0] = x[ii][0];
xi[1] = x[ii][1];
xi[2] = x[ii][2];
eij[0] = eij[1] = eij[2] = 0.0;
jlist = firstneigh[ii];
jnum = numneigh[ii];
for (int jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
local_cut2 = cut_spin_neel[itype][jtype]*cut_spin_neel[itype][jtype];
spj[0] = sp[j][0];
spj[1] = sp[j][1];
spj[2] = sp[j][2];
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*rij[0];
eij[1] = inorm*rij[1];
eij[2] = inorm*rij[2];
if (rsq <= local_cut2) {
compute_neel(ii,j,rsq,eij,fmi,spi,spj);
}
}
}
}
/* ---------------------------------------------------------------------- */