Commit 19aaf294 authored by julient31's avatar julient31
Browse files

Commit JT 092718 

- renamed pair/spin/long functions
- started to work on debugging ewald_dipole (force errors)
parent 6b4303c4

examples/SPIN/pppm_spin/in.spin.ewald_spin

+14 −13
Original line number Diff line number Diff line
@@ -19,30 +19,30 @@ create_atoms 1 box 


mass		1 58.93



#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

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/long 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/dipolar/long 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/long/qsymp 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0

pair_coeff 	* * eam/alloy ../examples/SPIN/pppm_spin/Co_PurjaPun_2012.eam.alloy Co

pair_coeff 	* * spin/exchange exchange 4.0 0.3593 1.135028015e-05 1.064568567

#pair_coeff 	* * eam/alloy ../examples/SPIN/pppm_spin/Co_PurjaPun_2012.eam.alloy Co

#pair_coeff 	* * spin/exchange exchange 4.0 0.3593 1.135028015e-05 1.064568567

#pair_coeff 	* * spin/long/qsymp long 8.0

pair_coeff 	* * spin/long long 8.0

pair_style	spin/dipolar/long 8.0

pair_coeff 	* * long 8.0



neighbor 	0.1 bin

neigh_modify 	every 10 check yes delay 20



#kspace_style pppm/dipole/spin 1.0e-4

#kspace_style ewald/dipole/spin 1.0e-4

kspace_style 	ewald/dipole/spin 1.0e-4

#kspace_modify 	mesh 32 32 32



#fix 		1 all precession/spin zeeman 1.0 0.0 0.0 1.0

fix 		1 all precession/spin zeeman 0.0 0.0 0.0 1.0

fix 		2 all langevin/spin 0.0 0.0 21

fix 		3 all nve/spin lattice yes

#fix 		3 all nve/spin lattice yes

fix 		3 all nve/spin lattice no



timestep	0.0001
@@ -57,7 +57,8 @@ variable magnorm equal c_out_mag[4] 
variable 	emag      equal c_out_mag[5]

variable 	tmag      equal c_out_mag[6]



thermo_style    custom step time v_magnorm v_emag temp etotal

thermo_style    custom step time v_magnorm v_tmag temp v_emag ke pe etotal

#thermo_style    custom step time v_magnorm v_emag temp etotal

thermo          10



compute 	outsp all property/atom spx spy spz sp fmx fmy fmz

examples/SPIN/pppm_spin/in.spin.pppm_spin

+8 −9
Original line number Diff line number Diff line
@@ -23,19 +23,18 @@ mass 1 58.93 
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/long 8.0

pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/long/qsymp 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0

pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/dipolar/long 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/dipolar/long/qsymp 8.0

pair_coeff 	* * eam/alloy ../examples/SPIN/pppm_spin/Co_PurjaPun_2012.eam.alloy Co

pair_coeff 	* * spin/exchange exchange 4.0 0.3593 1.135028015e-05 1.064568567

pair_coeff 	* * spin/long/qsymp long 8.0

#pair_coeff 	* * spin/long long 8.0

#pair_coeff 	* * spin/dipolar/long/qsymp long 8.0

pair_coeff 	* * spin/dipolar/long long 8.0



neighbor 	0.1 bin

neigh_modify 	every 10 check yes delay 20



kspace_style 	pppm/dipole/spin 1.0e-4

kspace_modify mesh 32 32 32

#kspace_modify 	mesh 32 32 32



#fix 		1 all precession/spin zeeman 1.0 0.0 0.0 1.0

fix 		1 all precession/spin zeeman 0.0 0.0 0.0 1.0
@@ -61,5 +60,5 @@ thermo 10 
compute 	outsp all property/atom spx spy spz sp fmx fmy fmz

dump 		100 all custom 1 dump_cobalt_hcp.lammpstrj type x y z c_outsp[1] c_outsp[2] c_outsp[3]



#run 		20000

run 		1
 
run 		20000

#run 		1

examples/SPIN/pppm_spin/in.spin.spin_dipolar_cut

+10 −8
Original line number Diff line number Diff line
@@ -19,16 +19,17 @@ create_atoms 1 box 


mass		1 58.93



#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

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/long 8.0

pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/dipolar/cut 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0 spin/dipolar/cut 8.0

pair_style 	spin/dipolar/cut 8.0

#pair_style 	hybrid/overlay eam/alloy spin/exchange 4.0

pair_coeff 	* * eam/alloy ../examples/SPIN/pppm_spin/Co_PurjaPun_2012.eam.alloy Co

pair_coeff 	* * spin/exchange exchange 4.0 0.3593 1.135028015e-05 1.064568567

pair_coeff 	* * spin/dipolar/cut long 8.0

#pair_coeff 	* * eam/alloy ../examples/SPIN/pppm_spin/Co_PurjaPun_2012.eam.alloy Co

#pair_coeff 	* * spin/exchange exchange 4.0 0.3593 1.135028015e-05 1.064568567

pair_coeff 	* * long 8.0

#pair_coeff 	* * spin/long long 8.0



neighbor 	0.1 bin
@@ -37,7 +38,8 @@ neigh_modify every 10 check yes delay 20 
#fix 		1 all precession/spin zeeman 1.0 0.0 0.0 1.0

fix 		1 all precession/spin zeeman 0.0 0.0 0.0 1.0

fix 		2 all langevin/spin 0.0 0.0 21

fix 		3 all nve/spin lattice yes

#fix 		3 all nve/spin lattice yes

fix 		3 all nve/spin lattice no



timestep	0.0001

src/KSPACE/ewald_dipole.cpp

+14 −4
Original line number Diff line number Diff line
@@ -449,9 +449,16 @@ void EwaldDipole::compute(int eflag, int vflag) 
      exprl = cs[kx][0][i]*cypz - sn[kx][0][i]*sypz;

      expim = sn[kx][0][i]*cypz + cs[kx][0][i]*sypz;



      // mu dot k product



      //muik = mu[i][0]*kx + mu[i][1]*ky + mu[i][2]*kz;





      // taking im of struct_fact x exp(i*k*ri) (for force calc.)



      partial = (muk[k][i])*(expim*sfacrl_all[k] + exprl*sfacim_all[k]);

      //partial = (muk[k][i])*(expim*sfacrl_all[k] - exprl*sfacim_all[k]);

      //partial = muik * (expim*sfacrl_all[k] + exprl*sfacim_all[k]);

      ek[i][0] += partial * eg[k][0];

      ek[i][1] += partial * eg[k][1];

      ek[i][2] += partial * eg[k][2];
@@ -493,6 +500,9 @@ void EwaldDipole::compute(int eflag, int vflag) 
    f[i][0] += muscale * ek[i][0];

    f[i][1] += muscale * ek[i][1];

    if (slabflag != 2) f[i][2] += muscale * ek[i][2];

    //f[i][0] -= muscale * ek[i][0];

    //f[i][1] -= muscale * ek[i][1];

    //if (slabflag != 2) f[i][2] -= muscale * ek[i][2];

    t[i][0] += -mu[i][1]*tk[i][2] + mu[i][2]*tk[i][1];

    t[i][1] += -mu[i][2]*tk[i][0] + mu[i][0]*tk[i][2];

    if (slabflag != 2) t[i][2] += -mu[i][0]*tk[i][1] + mu[i][1]*tk[i][0];
@@ -545,7 +555,7 @@ void EwaldDipole::compute(int eflag, int vflag) 
}



/* ----------------------------------------------------------------------

   compute the 

   compute the struc. factors and mu dot k products

------------------------------------------------------------------------- */



void EwaldDipole::eik_dot_r()

src/KSPACE/ewald_dipole_spin.cpp

+11 −5
Original line number Diff line number Diff line
@@ -53,7 +53,7 @@ EwaldDipoleSpin::EwaldDipoleSpin(LAMMPS *lmp, int narg, char **arg) : 
  hbar = force->hplanck/MY_2PI;         	// eV/(rad.THz)

  mub = 5.78901e-5;                     	// in eV/T

  mu_0 = 1.2566370614e-6;               	// in T.m/A

  mub2mu0 = mub * mub * mu_0;           // in eV

  mub2mu0 = mub * mub * mu_0 / (4.0*MY_PI);	// in eV

  mub2mu0hbinv = mub2mu0 / hbar;        	// in rad.THz

}
@@ -500,6 +500,9 @@ void EwaldDipoleSpin::compute(int eflag, int vflag) 
  const double spscale2 = mub2mu0hbinv * scale;

  //const double muscale = qqrd2e * scale;



  printf("test ek: %g %g %g \n",ek[0][0],ek[0][1],ek[0][2]);

  printf("test tk: %g %g %g \n",tk[0][0],tk[0][1],tk[0][2]);



  for (i = 0; i < nlocal; i++) {

    f[i][0] += spscale * ek[i][0];

    f[i][1] += spscale * ek[i][1];
@@ -509,6 +512,9 @@ void EwaldDipoleSpin::compute(int eflag, int vflag) 
    if (slabflag != 2) fm_long[i][2] += spscale2 * tk[i][3];

  }



  printf("test f_l: %g %g %g \n",f[0][0],f[0][1],f[0][2]);

  printf("test fm_l: %g %g %g \n",fm_long[0][0],fm_long[0][1],fm_long[0][2]);

  

  // sum global energy across Kspace vevs and add in volume-dependent term

  // taking the re-part of struct_fact_i x struct_fact_j

  // substracting self energy and scaling

CRA Git | Maintained and supported by SUSTech CRA and CCSE