Commit JT 033020

- modified all spin pairs (match nve)
- correct doc min_modify
- correct code max norm (square values)
- added draft nvt validation
This commit is contained in:
julient31 2020-03-30 08:09:11 -06:00
parent a739b8c6b7
commit 51e3f9dcda
28 changed files with 30374 additions and 74 deletions

View File

@ -91,7 +91,7 @@ The choice of a norm can be modified for the min styles *cg*\ , *sd*\
the 2-norm (Euclidean length) of the global force vector:
.. math::
|| \vec{F} ||_{2} = \sqrt{\vec{F}_1+ \cdots + \vec{F}_N}
|| \vec{F} ||_{2} = \sqrt{\vec{F}_1^2+ \cdots + \vec{F}_N^2}
The *max* norm computes the length of the 3-vector force
for each atom (2-norm), and takes the maximum value of those across

View File

@ -65,6 +65,6 @@ for t in range (0,N):
# calc. average magnetization
Sm = (S1+S2)*0.5
# calc. energy
en = -J0*(np.dot(S1,S2))
en = -2.0*J0*(np.dot(S1,S2))
# print res. in ps for comparison with LAMMPS
print(t*dt/1000.0,Sm[0],Sm[1],Sm[2],en)

View File

@ -13,7 +13,7 @@ en="$(echo "$en-$in" | bc -l)"
tail -n +$in log.lammps | head -n $en > res_lammps.dat
# compute Langevin
python3 -m llg_precession.py > res_llg.dat
python3 llg_precession.py > res_llg.dat
# plot results
python3 -m plot_precession.py res_lammps.dat res_llg.dat
python3 plot_precession.py res_lammps.dat res_llg.dat

View File

@ -3,7 +3,7 @@
units metal
atom_style spin
atom_modify map array
boundary p p p
boundary f f f
# read_data singlespin.data

View File

@ -1,12 +1,8 @@
#!/usr/bin/env python3
#Program fitting the exchange interaction
#Model curve: Bethe-Slater function
import numpy as np, pylab, tkinter
import matplotlib.pyplot as plt
import mpmath as mp
# from scipy.optimize import curve_fit
# from decimal import *
mub=5.78901e-5 # Bohr magneton (eV/T)
kb=8.617333262145e-5 # Boltzman constant (eV/K)
@ -24,4 +20,3 @@ tf=20.0
for i in range (0,npoints):
temp=ti+i*(tf-ti)/npoints
print('%lf %lf %lf' % (temp,func(temp),-g*mub*Hz*func(temp)))

View File

@ -20,7 +20,7 @@ do
done
# compute Langevin
python3 -m langevin.py > res_langevin.dat
python3 langevin.py > res_langevin.dat
# plot results
python3 -m plot_precession.py res_lammps.dat res_langevin.dat
python3 plot_precession.py res_lammps.dat res_langevin.dat

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,49 @@
# bcc iron in a 3d periodic box
clear
units metal
atom_style spin
dimension 3
boundary p p p
# necessary for the serial algorithm (sametag)
atom_modify map array
lattice bcc 2.8665
region box block 0.0 3.0 0.0 3.0 0.0 3.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for bcc iron
mass 1 55.845
set group all spin/random 31 2.2
velocity all create 100 4928459 rot yes dist gaussian
pair_style hybrid/overlay eam/alloy spin/exchange 3.5
pair_coeff * * eam/alloy Fe_Mishin2006.eam.alloy Fe
pair_coeff * * spin/exchange exchange 3.4 0.02726 0.2171 1.841
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all precession/spin zeeman 0.0 0.0 0.0 1.0
fix 2 all langevin/spin 0.0 0.00 21
fix 3 all nve/spin lattice moving
timestep 0.0001
# compute and output options
compute out_mag all spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable emag equal c_out_mag[5]
variable tmag equal c_out_mag[6]
thermo_style custom step time v_tmag temp v_emag ke pe etotal
thermo 200
run 100000

View File

@ -0,0 +1,46 @@
#!/usr/bin/env python3
import numpy as np, pylab, tkinter
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from decimal import *
import sys, string, os
argv = sys.argv
if len(argv) != 3:
print("Syntax: ./plot_precession.py res_lammps.dat")
sys.exit()
dirname = os.path.join(os.getcwd(), "Feb_07")
lammps_file = sys.argv[1]
llg_file = sys.argv[2]
t,tmag,temp,e_mag,e_kin,e_pot,e_tot = np.loadtxt(lammps_file,skiprows=0, usecols=(1,2,3,4,5,6,7),unpack=True)
fig = plt.figure(figsize=(8,8))
ax1 = plt.subplot(3,1,1)
ax2 = plt.subplot(3,1,2)
ax3 = plt.subplot(3,1,3)
ax1.plot(t, e_tot, 'k-', label='Total energy')
ax1.plot(t, e_pot, 'r-', label='Potential energy')
ax1.set_ylabel("E (eV)")
ax1.legend(loc=3)
ax2.plot(t, e_kin, 'b-', label='Kinetic energy')
ax2.plot(t, e_mag, 'g-', label='Magnetic energy')
ax2.set_ylabel("E (eV)")
ax2.legend(loc=3)
ax3.plot(t, temp, 'b--', label='Latt. temperature')
ax3.plot(t, tmag, 'r--', label='Spin temperature')
ax3.set_ylabel("T (K)")
ax3.legend(loc=3)
plt.xlabel('Time (in ps)')
plt.legend()
plt.show()
fig.savefig(os.path.join(os.getcwd(), "nve_spin_lattice.pdf"), bbox_inches="tight")
plt.close(fig)

View File

@ -0,0 +1,16 @@
#!/bin/bash
# clean old res
rm res_*.dat
# compute Lammps
./../../../../src/lmp_serial \
-in in.spin.iron-nve
in="$(grep -n Step log.lammps | awk -F ':' '{print $1}')"
en="$(grep -n Loop log.lammps | awk -F ':' '{print $1}')"
in="$(echo "$in+1" | bc -l)"
en="$(echo "$en-$in" | bc -l)"
tail -n +$in log.lammps | head -n $en > res_lammps.dat
# plot results
python3 -m plot_nve.py res_lammps.dat res_llg.dat

View File

@ -0,0 +1,39 @@
#LAMMPS in.run
units metal
atom_style spin
atom_modify map array
boundary f f f
read_data two_spins.data
pair_style spin/exchange 3.1
pair_coeff * * exchange 3.1 11.254 0.0 1.0
group bead type 1
variable H equal 0.0
variable Kan equal 0.0
variable Temperature equal 0.0
variable RUN equal 30000
fix 1 all nve/spin lattice no
fix 2 all precession/spin zeeman ${H} 0.0 0.0 1.0 anisotropy ${Kan} 0.0 0.0 1.0
fix_modify 2 energy yes
fix 3 all langevin/spin ${Temperature} 0.01 12345
compute out_mag all spin
compute out_pe all pe
variable magx equal c_out_mag[1]
variable magy equal c_out_mag[2]
variable magz equal c_out_mag[3]
variable magnorm equal c_out_mag[4]
variable emag equal c_out_mag[5]
thermo_style custom step time v_magx v_magy v_magz v_emag pe etotal
thermo 10
timestep 0.0001
run ${RUN}

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,51 @@
# bcc iron in a 3d periodic box
clear
units metal
atom_style spin
dimension 3
boundary p p p
# necessary for the serial algorithm (sametag)
atom_modify map array
lattice bcc 2.8665
region box block 0.0 3.0 0.0 3.0 0.0 3.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for bcc iron
mass 1 55.845
set group all spin 2.2 0.0 0.0 1.0
velocity all create 400 4928459 rot yes dist gaussian
pair_style hybrid/overlay eam/alloy spin/exchange 4.0 spin/neel 4.0
pair_coeff * * eam/alloy Fe_Mishin2006.eam.alloy Fe
pair_coeff * * spin/exchange exchange 3.4 0.1 0.2171 1.841
pair_coeff * * spin/neel neel 4.0 0.02 0.0 1.841 0.0 0.0 1.0
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all precession/spin zeeman 0.0 0.0 0.0 1.0
fix 2 all langevin 200.0 200.0 10.0 48279
fix 3 all langevin/spin 0.0 0.00001 321
fix 4 all nve/spin lattice moving
timestep 0.001
# compute and output options
compute out_mag all spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable emag equal c_out_mag[5]
variable tmag equal c_out_mag[6]
thermo_style custom step time v_tmag temp v_emag ke pe etotal
thermo 200
run 200000

View File

@ -0,0 +1,49 @@
# bcc iron in a 3d periodic box
clear
units metal
atom_style spin
dimension 3
boundary p p p
# necessary for the serial algorithm (sametag)
atom_modify map array
lattice bcc 2.8665
region box block 0.0 3.0 0.0 3.0 0.0 3.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for bcc iron
mass 1 55.845
set group all spin 2.2 0.0 0.0 1.0
velocity all create 0 4928459 rot yes dist gaussian
pair_style hybrid/overlay eam/alloy spin/exchange 3.5
pair_coeff * * eam/alloy Fe_Mishin2006.eam.alloy Fe
pair_coeff * * spin/exchange exchange 3.4 0.1 0.2171 1.841
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all precession/spin zeeman 0.0 0.0 0.0 1.0
fix 2 all langevin/spin 200.0 0.1 321
fix 3 all nve/spin lattice moving
timestep 0.001
# compute and output options
compute out_mag all spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable emag equal c_out_mag[5]
variable tmag equal c_out_mag[6]
thermo_style custom step time v_tmag temp v_emag ke pe etotal
thermo 200
run 200000

View File

@ -0,0 +1,43 @@
#!/usr/bin/env python3
import numpy as np, pylab, tkinter
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from decimal import *
import sys, string, os
argv = sys.argv
if len(argv) != 3:
print("Syntax: ./plot_precession.py res_nvt_spin.dat res_nvt_lattice.dat")
sys.exit()
dirname = os.path.join(os.getcwd(), "Feb_07")
nvtspin_file = sys.argv[1]
nvtlatt_file = sys.argv[2]
ts,tmags,temps = np.loadtxt(nvtspin_file,skiprows=0, usecols=(1,2,3),unpack=True)
tl,tmagl,templ = np.loadtxt(nvtlatt_file,skiprows=0, usecols=(1,2,3),unpack=True)
fig = plt.figure(figsize=(8,8))
ax1 = plt.subplot(2,1,1)
ax2 = plt.subplot(2,1,2)
ax1.plot(ts, tmags, 'r-', label='Spin temp. (thermostat)')
ax1.plot(ts, temps, 'g-', label='Lattice temp.')
ax1.set_yscale("log")
ax1.set_ylabel("T (K)")
ax1.legend(loc=3)
ax2.plot(tl, tmagl, 'r-', label='Spin temp.')
ax2.plot(tl, templ, 'g-', label='Lattice temp. (thermostat)')
ax2.set_yscale("log")
ax2.set_ylabel("T (K)")
ax2.legend(loc=3)
plt.xlabel('Time (in ps)')
plt.legend()
plt.show()
fig.savefig(os.path.join(os.getcwd(), "nve_spin_lattice.pdf"), bbox_inches="tight")
plt.close(fig)

View File

@ -0,0 +1,23 @@
#!/bin/bash
# clean old res
rm res_*.dat
# compute NVT Spin -> Lattice
./../../../../src/lmp_serial -in in.spin.nvt_spin
in="$(grep -n Step log.lammps | awk -F ':' '{print $1}')"
en="$(grep -n Loop log.lammps | awk -F ':' '{print $1}')"
in="$(echo "$in+1" | bc -l)"
en="$(echo "$en-$in" | bc -l)"
tail -n +$in log.lammps | head -n $en > res_nvt_spin.dat
# compute NVT Lattice -> Spin
./../../../../src/lmp_serial -in in.spin.nvt_lattice
in="$(grep -n Step log.lammps | awk -F ':' '{print $1}')"
en="$(grep -n Loop log.lammps | awk -F ':' '{print $1}')"
in="$(echo "$in+1" | bc -l)"
en="$(echo "$en-$in" | bc -l)"
tail -n +$in log.lammps | head -n $en > res_nvt_lattice.dat
# plot results
python3 plot_nvt.py res_nvt_spin.dat res_nvt_lattice.dat

View File

@ -48,7 +48,6 @@ ComputeSpin::ComputeSpin(LAMMPS *lmp, int narg, char **arg) :
vector_flag = 1;
size_vector = 6;
extvector = 0;
// npairs = npairspin = 0;
// initialize the magnetic interaction flags
@ -180,15 +179,11 @@ void ComputeSpin::compute_vector()
mag[0] += sp[i][0];
mag[1] += sp[i][1];
mag[2] += sp[i][2];
// magenergy -= (sp[i][0]*fm[i][0] + sp[i][1]*fm[i][1] + sp[i][2]*fm[i][2]);
// update magnetic precession energies
if (precession_spin_flag) {
magenergy += lockprecessionspin->emag[i];
// magenergy -= lockprecessionspin->compute_zeeman_energy(sp[i]);
// magenergy -= lockprecessionspin->compute_anisotropy_energy(sp[i]);
// magenergy -= lockprecessionspin->compute_cubic_energy(sp[i]);
}
// update magnetic pair interactions
@ -222,13 +217,13 @@ void ComputeSpin::compute_vector()
magtot[2] *= scale;
magtot[3] = sqrt((magtot[0]*magtot[0])+(magtot[1]*magtot[1])+(magtot[2]*magtot[2]));
spintemperature = hbar*tempnumtot;
spintemperature /= (2.0*kb*tempdenomtot);
// spintemperature /= (2.0*kb*tempdenomtot);
spintemperature /= (kb*tempdenomtot);
vector[0] = magtot[0];
vector[1] = magtot[1];
vector[2] = magtot[2];
vector[3] = magtot[3];
// vector[4] = magenergytot*hbar;
vector[4] = magenergytot;
vector[5] = spintemperature;

View File

@ -155,8 +155,6 @@ FixPrecessionSpin::FixPrecessionSpin(LAMMPS *lmp, int narg, char **arg) : Fix(lm
FixPrecessionSpin::~FixPrecessionSpin()
{
delete [] magstr;
// test emag list storing mag energies
memory->destroy(emag);
}
@ -217,7 +215,6 @@ void FixPrecessionSpin::init()
if (varflag == CONSTANT) set_magneticprecession();
// test emag list storing mag energies
// init. size of energy stacking lists
nlocal_max = atom->nlocal;
@ -263,8 +260,8 @@ void FixPrecessionSpin::post_force(int /* vflag */)
const int nlocal = atom->nlocal;
double spi[4], fmi[3], epreci;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
@ -273,10 +270,7 @@ void FixPrecessionSpin::post_force(int /* vflag */)
eflag = 0;
eprec = 0.0;
for (int i = 0; i < nlocal; i++) {
// test emag list storing mag energies
emag[i] = 0.0;
if (mask[i] & groupbit) {
epreci = 0.0;
spi[0] = sp[i][0];

View File

@ -57,8 +57,9 @@ class FixPrecessionSpin : public Fix {
void compute_cubic(double *, double *);
double compute_cubic_energy(double *);
// test emag list storing mag energies
int nlocal_max; // max value of nlocal (for size of lists)
// storing magnetic energies
int nlocal_max; // max nlocal (for list size)
double *emag; // energy list
protected:

View File

@ -99,7 +99,6 @@ void PairSpin::init_style()
if (ifix >=0)
lattice_flag = ((FixNVESpin *) modify->fix[ifix])->lattice_flag;
// test emag list storing mag energies
// init. size of energy stacking lists
nlocal_max = atom->nlocal;

View File

@ -32,8 +32,9 @@ friend class FixNVESpin;
virtual void compute(int, int) {}
virtual void compute_single_pair(int, double *) {}
// test emag list storing mag energies
int nlocal_max; // max value of nlocal (for size of lists)
// storing magnetic energies
int nlocal_max; // max nlocal (for list size)
double *emag; // energy list
protected:

View File

@ -64,8 +64,6 @@ PairSpinDipoleCut::~PairSpinDipoleCut()
memory->destroy(setflag);
memory->destroy(cut_spin_long);
memory->destroy(cutsq);
// test emag list storing mag energies
memory->destroy(emag);
}
}
@ -188,8 +186,8 @@ void PairSpinDipoleCut::compute(int eflag, int vflag)
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
@ -209,6 +207,7 @@ void PairSpinDipoleCut::compute(int eflag, int vflag)
spi[1] = sp[i][1];
spi[2] = sp[i][2];
spi[3] = sp[i][3];
emag[i] = 0.0;
itype = type[i];
for (jj = 0; jj < jnum; jj++) {
@ -262,7 +261,9 @@ void PairSpinDipoleCut::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 *= 0.5*hbar;
// evdwl *= 0.5*hbar;
evdwl *= hbar;
emag[i] += evdwl;
}
} else evdwl = 0.0;

View File

@ -69,8 +69,6 @@ PairSpinDipoleLong::~PairSpinDipoleLong()
memory->destroy(setflag);
memory->destroy(cut_spin_long);
memory->destroy(cutsq);
// test emag list storing mag energies
memory->destroy(emag);
}
}
@ -215,8 +213,8 @@ void PairSpinDipoleLong::compute(int eflag, int vflag)
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
@ -242,8 +240,6 @@ void PairSpinDipoleLong::compute(int eflag, int vflag)
spi[1] = sp[i][1];
spi[2] = sp[i][2];
spi[3] = sp[i][3];
// test emag list storing mag energies
emag[i] = 0.0;
for (jj = 0; jj < jnum; jj++) {
@ -309,7 +305,8 @@ void PairSpinDipoleLong::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 *= 0.5*hbar;
// evdwl *= 0.5*hbar;
evdwl *= hbar;
emag[i] += evdwl;
}
} else evdwl = 0.0;

View File

@ -53,8 +53,6 @@ PairSpinDmi::~PairSpinDmi()
memory->destroy(vmech_dmy);
memory->destroy(vmech_dmz);
memory->destroy(cutsq);
// test emag list storing mag energies
memory->destroy(emag);
}
}
@ -194,8 +192,8 @@ void PairSpinDmi::compute(int eflag, int vflag)
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
@ -216,8 +214,6 @@ void PairSpinDmi::compute(int eflag, int vflag)
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
// test emag list storing mag energies
emag[i] = 0.0;
// loop on neighbors
@ -271,7 +267,8 @@ void PairSpinDmi::compute(int eflag, int vflag)
if (eflag) {
evdwl -= (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
evdwl *= 0.5*hbar;
// evdwl *= 0.5*hbar;
evdwl *= hbar;
emag[i] += evdwl;
} else evdwl = 0.0;

View File

@ -50,8 +50,6 @@ PairSpinExchange::~PairSpinExchange()
memory->destroy(J2);
memory->destroy(J3);
memory->destroy(cutsq); // to be implemented
// test emag list storing mag energies
memory->destroy(emag);
}
}
@ -179,9 +177,9 @@ void PairSpinExchange::compute(int eflag, int vflag)
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
}
@ -201,8 +199,6 @@ void PairSpinExchange::compute(int eflag, int vflag)
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
// test emag list storing mag energies
emag[i] = 0.0;
// loop on neighbors
@ -255,12 +251,9 @@ void PairSpinExchange::compute(int eflag, int vflag)
if (eflag) {
evdwl -= (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
evdwl *= 0.5*hbar;
// printf("test ex energy: %g \n",evdwl);
// evdwl = -0.5*compute_energy(i,j,rsq,spi,spj);
// printf("test ex energy: %g \n",evdwl);
// evdwl *= 0.5*hbar;
evdwl *= hbar;
emag[i] += evdwl;
// evdwl *= hbar;
} else evdwl = 0.0;
if (evflag) ev_tally_xyz(i,j,nlocal,newton_pair,

View File

@ -51,8 +51,6 @@ PairSpinMagelec::~PairSpinMagelec()
memory->destroy(v_mey);
memory->destroy(v_mez);
memory->destroy(cutsq); // to be deteled
// test emag list storing mag energies
memory->destroy(emag);
}
}
@ -188,8 +186,8 @@ void PairSpinMagelec::compute(int eflag, int vflag)
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
@ -210,8 +208,6 @@ void PairSpinMagelec::compute(int eflag, int vflag)
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
// test emag list storing mag energies
emag[i] = 0.0;
// loop on neighbors
@ -264,7 +260,8 @@ void PairSpinMagelec::compute(int eflag, int vflag)
if (eflag) {
evdwl -= (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
evdwl *= 0.5*hbar;
// evdwl *= 0.5*hbar;
evdwl *= hbar;
emag[i] += evdwl;
} else evdwl = 0.0;

View File

@ -54,8 +54,6 @@ PairSpinNeel::~PairSpinNeel()
memory->destroy(q2);
memory->destroy(q3);
memory->destroy(cutsq); // to be deleted
// test emag list storing mag energies
memory->destroy(emag);
}
}
@ -193,8 +191,8 @@ void PairSpinNeel::compute(int eflag, int vflag)
numneigh = list->numneigh;
firstneigh = list->firstneigh;
// test emag list storing mag energies
// checking size of emag
if (nlocal_max < nlocal) { // grow emag lists if necessary
nlocal_max = nlocal;
memory->grow(emag,nlocal_max,"pair/spin:emag");
@ -215,8 +213,6 @@ void PairSpinNeel::compute(int eflag, int vflag)
spi[0] = sp[i][0];
spi[1] = sp[i][1];
spi[2] = sp[i][2];
// test emag list storing mag energies
emag[i] = 0.0;
// loop on neighbors
@ -272,9 +268,9 @@ void PairSpinNeel::compute(int eflag, int vflag)
}
if (eflag) {
// evdwl = (spi[0]*fmi[0] + spi[1]*fmi[1] + spi[2]*fmi[2]);
evdwl = compute_neel_energy(i,j,rsq,eij,spi,spj);
evdwl *= 0.5*hbar;
// evdwl *= 0.5*hbar;
evdwl *= hbar;
emag[i] += evdwl;
} else evdwl = 0.0;

View File

@ -902,13 +902,13 @@ double Min::fnorm_inf()
double local_norm_inf = 0.0;
for (i = 0; i < nvec; i++)
local_norm_inf = MAX(fabs(fvec[i]),local_norm_inf);
local_norm_inf = MAX(fvec[i]*fvec[i],local_norm_inf);
if (nextra_atom) {
for (int m = 0; m < nextra_atom; m++) {
fatom = fextra_atom[m];
n = extra_nlen[m];
for (i = 0; i < n; i++)
local_norm_inf = MAX(fabs(fatom[i]),local_norm_inf);
local_norm_inf = MAX(fatom[i]*fatom[i],local_norm_inf);
}
}
@ -917,7 +917,7 @@ double Min::fnorm_inf()
if (nextra_global)
for (i = 0; i < nextra_global; i++)
norm_inf = MAX(fabs(fextra[i]),norm_inf);
norm_inf = MAX(fextra[i]*fextra[i],norm_inf);
return norm_inf;
}