ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Commit1 JT 052918 modifs files

Browse Source
This commit is contained in:
julient31
2018-05-29 14:51:21 -06:00
parent 5f0e6d0aa7
commit 6457e5eedb
11 changed files with 80 additions and 88 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
examples/SPIN/cobalt_fcc/in.spin.cobalt_fcc
View File

@ -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

3
examples/SPIN/cobalt_hcp/in.spin.cobalt_hcp
View File

@ -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

2
examples/SPIN/iron/in.spin.iron
View File

@ -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

4
src/SPIN/atom_vec_spin.cpp
View File

@ -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);
}

8
src/SPIN/fix_nve_spin.cpp
View File

@ -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++) {

6
src/SPIN/fix_precession_spin.cpp
View File

@ -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;

2
src/SPIN/fix_precession_spin.h
View File

@ -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

14
src/SPIN/pair_spin_dmi.cpp
View File

@ -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];

17
src/SPIN/pair_spin_exchange.cpp
View File

@ -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]);
}

97
src/SPIN/pair_spin_magelec.cpp
View File

@ -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]);
// 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++;
}
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) {

14
src/SPIN/pair_spin_neel.cpp
View File

@ -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];