The problem with magnetic spins

LAMMPS LAMMPS General Discussion
1

Dear LAMMPS users,

I try to study the spins dynamics of alpha-Fe2O3.

The lattice potential term is Buckingham potential. Since the spin atom style does not support “pair_style buck/coul/long”, the atom style is set as “atom_style hybrid charge spin” .

But the total energy is about -2.155138e+11 eV, and lammps always reports “Out of range atoms - cannot compute PPPM”.

However, if I only consider the lattice dynamics without any spins, the program runs smoothly without errors, the crystal also maintains a stable structure。 The total energy is about -16146.742 eV.

On the one hand, if I only consider the spins dynamics and shield the movement of lattice, it still works. The total energy is about 83 eV.

Only when I put them together will things go wrong.

# Fe2O3 crystal in a 3d periodic box

**clear **
variable T equal 300
units metal
atom_style hybrid charge spin

dimension 3
boundary p p p

**atom_modify map array **

box tilt large

read_data “system.data”

#Fe
set type 1 spin/random 100 5.0

set type 1 charge 1.418

#O
set type 2 spin/random 100 0.0000001

set type 2 charge -0.945333333

velocity all create ${T} 4928459 rot yes dist gaussian

pair_style hybrid/overlay buck/coul/long 10.0 spin/exchange 4.0

pair_coeff 1 1 spin/exchange exchange 4.0 0.19014831 0.2881564913 0.317272971 offset yes

pair_coeff 1 2 spin/exchange exchange 4.0 0.0 0.2171 1.841 offset yes

pair_coeff 2 2 spin/exchange exchange 4.0 0.0 0.2171 1.841 offset yes

pair_coeff 1 1 buck/coul/long 8020.28853 0.228 42.546225

pair_coeff 1 2 buck/coul/long 8521.55602 0.190 0.0

pair_coeff 2 2 buck/coul/long 9022.82351 0.265 85.0924495

kspace_style pppm 1.0e-4

neighbor 0.1 bin

neigh_modify every 10 check yes delay 20

fix 1 all nve/spin lattice moving

fix 2 all langevin/spin ${T} 0.01 21

timestep 0.0001

# compute and output options
#

thermo 10000

compute outsp all property/atom spx spy spz sp fmx fmy fmz
dump 1 all custom 10000 dump_iron.lammpstrj type x y z c_outsp[1] c_outsp[2] c_outsp[3]

run 1000000

clear

I am looking forward to your reply,

Best wishes,

Wei

2

This is a question for the author of the SPIN package @julient

3

Hello @singlebook Wei and @akohlmey Axel!

Ah, this is a very interesting question Wei.
Considering the very large value of your energy, I’m guessing the system ‘exploded’ because of the mechanical forces generated by the spin Hamiltonian.
What is the initial pressure in your system when you have both spins and charges? I tried to run your script, but I do not have your input configuration (“system.data”). But I would assume your initial pressure was very wrong.

Here is my guess about the problem:
As you probably know, when you are running a spin-lattice calculation, your magnetic Hamiltonian also generates mechanical forces. So that in principle, the magnetic Hamiltonian needs to be fitted on magnetic energies and on some mechanical information.
For example, in this manuscript:
https://www.nature.com/articles/s41524-021-00617-2
You can see on Fig. 5 that the spin Hamiltonian is fitted on the energy and pressure fluctuations provided by spin-spiral DFT calculations.

I assume this was not the case in your calculation: your magnetic Hamiltonian is probably predicting decent magnetic energies, but very large (and incorrect) mechanical forces. You can probe that by outputting the mechanical forces in a dump file (so adding fx fy and fz to your dump custom command) when running a purely spin dynamics calculation.

Also, what is happening when you start your simulation with an equilibrium spin configuration (so not starting with a ‘spin random’), and without spin temperature(so setting your T variable to 0.0)? As you have the “offset yes” command in the definition of your magnetic pair_coeff, this should improve the situation. Is it the case?

Let me know when you have those elements of answer. But I need to have those confirmations before we can chat more of potential solutions to your problem.

All the best,
Julien.

4

Hello, Julien.

Thank you very much for your reply!

As a newly registered user of the forum, I have no permission to upload attachments. May I send the simulation results to your email？

As you suggested, I checked the force and pressure for hybrid atom style. The pressure is -3.1684489e+13 bars and the max force is 2.5195862e+11 eV/A.

However, as for the hybrid atom style, even if I set the spin of atom or exchange energy to zero, the values of force and pressure do not change.

I performed two other simulations separately. The atomic style of one is charge, and the atomic style of the other is spin. They both run smoothly! The max force for charge style is about 2 eV/A and pressure is
21796.224 bars. The crystal maintains a stable structure. And the force and pressure for spin style are 0 eV/A and 4102 bars.

Best wishes,

Wei

5

Hello Wei,

My apologies for the late reply.

Sure, feel free to send me some of your simulation results (and, if you can, your input data file) at my email: [email protected] It will allow me to try to run some tests of my own.
But for the benefit of other users, it is also important we keep exchanging here.

Hmm, the last result is very surprising (spin style only). The forces are probably computed but set to zero but I would assume the pressure to be very large.
For the hybrid situation, do you obtain that pressure and force after some steps? Or before you run any step?

Best,
Julien.