First version of the spin tutorial (2)

Examples (example/SPIN), for BFO and Co
This commit is contained in:
julient31 2017-08-16 09:32:35 -06:00
parent 023b018ed2
commit 45ea7b3cc7
5 changed files with 695 additions and 0 deletions

255
doc/src/tutorial_spin.txt Normal file
View File

@ -0,0 +1,255 @@
<script type="text/javascript"
src="https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML">
</script>
<script type="text/x-mathjax-config">
MathJax.Hub.Config({ TeX: { equationNumbers: {autoNumber: "AMS"} } });
</script>
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
:link(lws,http://lammps.sandia.gov)
:link(ld,Manual.html)
:link(lc,Section_commands.html#comm)
:line
Tutorial for the SPIN package in LAMMPS :h3
This tutorial explains how to use the spin dynamics in LAMMPS,
and to perform spin and spin--lattice simulations using the
SPIN package. As an illustration, input files are documented.
First of all, LAMMPS has to be compiled with the SPIN package
activated. Then, the data file and input scripts have to be
modified to include the magnetic spins and to handle them.
:line
[Overview of the spin and spin--lattice dynamics]
At the atomic scale, magnetic materials can be seen as ensemble of
atoms, each one of them carying a magnetic moment, refered to as its
atomic spin. In ref "Antropov"_#Antropov, a formalism allowing to
simulate approximate classical spins, and to couple them to lattice
vibrations was introduced. Each of these spins is simulated via a
classical vector, associated to each magnetic atom, and whose
trajectory is defined by an equation of motion.
Lattice vibrations are simulated by the usual equations of MD.
A mechanical potential (EAM, Finnis-Siclair, or Dudarev-Derlet) ensure
the cohesion of the particles. The coupling between the magnetic and
lattice degrees of freedom is performed via the inter-atomic dependence
of the magnetic interactions.
:ole
:line
[Preparation of the data file]
For the mechanical potentials, the data file is similar to a standard LAMMPS
data file for {atom_style full}. The DPs and the {harmonic bonds} connecting
them to their DC should appear in the data file as normal atoms and bonds.
For the magnetic interactions, no data file is necessary, every interaction
input will be define in the input file. :pre
:line
[Basic input file]
Up to know, spin simulations only accept the metal units.
The atom style should be set to {spin}, so that you
can define atomic spin vectors.
The set group command defines the Lande factor (the norm) of
the magnetic vectors in a given group, and their initial
oriantation. The command is:
set group A B C D E F :pre
with A, B, C, D, E, and F the following input parameters:
A is set to all, or to the number of a specific and pre-defined
group of atoms,
B is set to {spin} or {spin/random}, depending on if the spins of
the group have to be initialized in a particulat direction, or randomly.
If B is defined as {spin}, the C is the Lande factor for the spins of
the group, and [D,E,F] is the vector giving the direction for the
initialization.
If B is defined as {spin/random}, C is a number giving the seed for the random
difinition of the directions, and D is the Lande factor of the spins
in this group. E and F are not defined.
Examples:
set group all spin 1.72 1.0 0.0 0.0
set group all spin/random 11 1.72
setting the initial direction of all spins, with a Lande factor of 1.72,
and either in the x direction, or randomly. :l
The pair style has to be set to {pair/spin}. The command is
pair_style pair/spin A
with A a global radius cutoff.
The different pair interactions and their associated coefficients can then be
defined via the pair coeff command.
pair_coeff A B C D E F G H
where A and B are setting the pair concerned by this pair coeff command.
For example, A=1 and B=2 to set the pair coefficients between spins of
type 1 and spins of type 2. Use {*} for setting a pair coeffcient between
all pairs of spins.
C defines the type of the interaction. It can be set to {exchange} for an
exchange interaction, {dmi} for a Dzyaloshinskii-Moriya (DM) interaction, or to
{me} for a magneto-electric (ME) interaction.
If C is set to {exchange}, D is the radius cutoff (in Angtrom) associated
to the exchange interaction, and E, F and G are the three parameters of the
Bethe--Slater function (E in eV, F without dimension, and G in Angtrom).
H is not defined.
If C is set to {dmi}, D is the radius cutoff (in Angtrom) associated to the DM
interaction, E is the intensity of the interaction in eV (which is also the
norm of the DM vector), and [F, G, H] are giving the direction of the DM
vector.
If C is set to {me}, D is the radius cutoff (in Angtrom) associated to the ME
interaction, E is the intensity of the interaction in eV (which corresponds to
the intensity of the electric polarization), and [F, G, H] are giving the
direction of the polarization vector.
Examples:
pair_coeff * * exchange 4.0 0.0446928 0.003496 1.4885
pair_coeff * * dmi 2.6 0.001 0.0 0.0 1.0
are setting an exchange interaction between every type of spins, with a radius
cutoff of 4.0 Angtrom, and E=0.0446928 eV, F=0.003496, and G=1.4885 Angtrom as
coefficient for the associated Bethe--Slater function, and a DM interaction
with a radius cutoff of 2.6 Angtrom, an intensity of 0.001 eV and a DM vector
in the direction [0.0 0.0 1.0]. :l
A fix has to be set to {force/spin} to define local magnetic forces, like the Zeeman
interaction or an anisotropic interaction. The command is:
fix A B C D E F G H
with A the label of the associated fix, B defining to which group of spins the
force is applied to, C defined as {force/spin} for magnetic local fixes, and
D defining the type of this fix. D can be equal to {zeeman} for a Zeeman interaction,
or to {anisotropy} for an anisotropic interaction.
If D is set to {zeeman}, then E gives the intensity of the applied magnetic field (in
Tesla), and [F, G, H] the direction of this field.
If is set to {anisotropy}, hen E gives the intensity of the anisotropic field (in eV),
and [F, G, H] the direction of this anisotropy.
Examples:
fix 1 all force/spin zeeman 1.0 0.0 0.0 1.0
fix 2 all force/spin anisotropy 0.001 0.0 1.0 0.0
are setting two fixes, the first one, labelled 1, is applied to all spins and defines
a Zeeman interaction corresponding to a field of 1 Tesla in the z direction, whereas
the second fix, labelled 2, is also applied to all spins, and defines a uniaxial
anisotropy of intensity 0.001 eV, in the direction y. :l
To simulate the temperature effects, a fix has to be set to {langevin/spin}. The command
is
fix A B C D E F G
where A is the label of the associated fix, B defines to which spins the fix is applied,
and C is set equal to {langevin/spin}. Then, D defines the temperature of the random bath
(in K), E is the transverse magnetic damping (no dimension), F is a longitudinal magnetic
damping, and G the seed for the random variables.
Note: the transverse damping is not implemented into LAMMPS yet. It is necessary for
micromagnetic simulations only.
Examples:
fix 2 all langevin/spin 300.0 0.01 0.0 21
is setting a fix labelled as 2, which is connecting all spins to a random bath. The temperature
of this bath is set to 300 K, and the value of the transverse Gilbert damping is set to 0.01.
The seed is set to 21.
For LAMMPS to understand that the motion of spins has to be taken into account, one has to set
a fix to {nve/spin}. The command is:
fix A B C D
with A the label of this fix, B defining to which group of atoms this fix is applied to, C has
to be set to {nve/spin}, and D can be set to {serial} for a serial calculation, with one
processor only, or to {mpi} for a parallel calculation, with N processors.
Example:
fix 3 all nve/spin mpi
is setting a fix labelled 3, and applies it to all the particles. The calculation is running in
parallel as the option {mpi} is defined.
Two main outputs of the computed magnetic quntities can be performed.
The first one is via a compute and a fix options. The compute command is defined as:
compute A B C
with A the label of this compute, B to which group of particles it is applied to, and finally,
the option {compute/spin} has to be set.
This compute is associated to a fix to define the output frequency and format. A typical command is:
fix A B C D E F G H
where A is the label of the fix, B the group of particles it is applied to, C defines the type of the
output fix, for example we chose to use {ave/time}, which can output every N timesteps, and perform
a time average over those steps. D, E, and F are the usual command of {ave/time}. G stands for the
magnetic quantities that are computed by {compute/spin}. Those quantities are:
c_mag[1] Physical time (in ps)
c_mag[2] Magnetization along x (adim)
c_mag[3] Magnetization along y (adim)
c_mag[4] Magnetization along z (adim)
c_mag[5] Magnetization Norm (adim)
c_mag[6] Magnetic energy (in eV)
c_mag[7] Spin temperature (in K)
And H stands for the output format, and is defined as in other .
Example:
compute 1 all compute/spin
fix 3 all ave/time 1 1 10 c_mag[1] c_mag[2] c_mag[3] c_mag[4] c_mag[5] c_mag[6] c_mag[7]
file mag_output.dat format %20.16g
is defining a compute of the magnetic quantities applied to all spins. The fix then outputs
every magnetic quantities every 10 time steps without performing any time average. These
quantities are stored in a file called mag_output.dat, with a special format (defined by %20.16g)
It is also possible to output the evolution of the spin configuration. This can be done with
the dump command. The only difference with a regular dump command is that the velocities
vi are replaced by the spin components spi.
Example:
dump 1 all custom 100 dump_spin.lammpstrj type x y z spx spy spz
is dumping every 100 timesteps the spin configuration in a file called dump_spin.lammpstrj.
:line
:link(Antropov)
[(Antropov)] Antropov, Katsnelson, Harmon, Van Schilfgaarde, and Kusnezov, Phys Rev B, 54(2), 1019 (1996)

109
examples/SPIN/in.BFO Normal file
View File

@ -0,0 +1,109 @@
###################
#######Init########
###################
clear
#setting units to metal (Ang, picosecs, eV, ...):
units metal
#setting dimension of the system (N=2 or 3):
dimension 3
#setting boundary conditions. (p for periodic, f for fixed, ...)
boundary p p f
#setting atom_style to spin:
atom_style spin
#Define sort for paramagnetic simulations (if no pair interaction)
#atom_modify sort 1000 4.0
#why?
atom_modify map array
#newton off for pair spin in SEQNEI
#newton off off
###########################
#######Create atoms########
###########################
#Lattice constant of sc Iron atoms of BFO
lattice sc 3.96
#Defining a geometric region of space. Sets ID(user's choice), style(block, sphere, ...), then, args depends on the style chosen
#(for block, one has x0, xf, y0, yf, z0, zf, in distance units)
region box block 0.0 17.0 0.0 17.0 0.0 5.0
#Creating a simulation box based on the specified region. Entries: number of atom types and box ref.
create_box 1 box
#Creating atoms (or molecules) on a lattice, or a single atom (or molecule), ...
#Entries: atom type,
create_atoms 1 box
#Replicating NxNxN the entire set of atoms
#replicate 1 1 1
#######################
#######Settings########
#######################
#Setting one or more properties of one or more atoms.
#Setting mass
mass 1 1.0
#set group all mass 1.0
#Setting spins orientation and moment
set group all spin/random 11 2.50
#set group all spin 2.50 1.0 0.0 0.0
#Magnetic exchange interaction coefficient for bulk fcc Cobalt
pair_style pair/spin 5.7
#type i and j | interaction type | cutoff | J1 (eV) | J2 (adim) | J3 (Ang) (for Exchange)
pair_coeff * * exchange 5.7 -0.01575 0.0 1.965
#type i and j | interaction type | cutoff | Int (eV) | [dx,dy,dz] (for DMI and ME)
#pair_coeff * * dmi 2.6 0.001 0.0 0.0 1.0
#pair_coeff * * me 4.0 0.000109 1.0 1.0 1.0
#Mex10
pair_coeff * * me 4.0 0.00109 1.0 1.0 1.0
#Define a skin distance, update neigh list every
#neighbor 1.0 bin
#neigh_modify every 10 check yes delay 20
neighbor 0.0 bin
neigh_modify every 1 check no delay 0
#Magnetic field fix
#Type | Intensity (T or eV) | Direction
#fix 1 all force/spin zeeman 1.0 0.0 0.0 1.0
fix 1 all force/spin anisotropy 0.0000035 0.0 0.0 1.0
#Fix Langevin spins (merging damping and temperature)
#Temp | Alpha_trans | Alpha_long | Seed
#fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 2 all langevin/spin 0.0 0.1 0.0 21
#fix 2 all langevin/spin 0.0 0.0 0.0 21
#Magnetic integration fix
fix 3 all nve/spin mpi
#compute real time, total magnetization, magnetic energy, and spin temperature
#Iteration | Time | Mx | My | Mz | |M| | Em | Tm
compute mag all compute/spin
fix outmag all ave/time 1 1 500 c_mag[1] c_mag[2] c_mag[3] c_mag[4] c_mag[5] c_mag[6] c_mag[7] file mag_BFO.dat format %20.16g
#Setting the timestep for the simulation (in ps)
timestep 0.0001
##################
#######run########
##################
#Dump the positions and spin directions of magnetic particles (vmd format)
dump 1 all custom 5000 dump_spin_BFO.lammpstrj type x y z spx spy spz
#Running the simulations for N timesteps
run 30000
#run 1

126
examples/SPIN/in.cobalt_SD Normal file
View File

@ -0,0 +1,126 @@
###################
#######Init########
###################
clear
#setting units to metal (Ang, picosecs, eV, ...):
units metal
#setting dimension of the system (N=2 or 3):
dimension 3
#setting boundary conditions. (p for periodic, f for fixed, ...)
boundary p p p
#boundary f f f
#setting atom_style to spin:
atom_style spin
#Define sort for paramagnetic simulations (if no pair interaction)
#atom_modify sort 1000 4.0
#why?
atom_modify map array
#newton off for pair spin in SEQNEI
#newton off off
###########################
#######Create atoms########
###########################
#Lattice constant of fcc Cobalt
lattice fcc 3.54
#lattice sc 2.50
#Test Kagome
#variable a equal sqrt(3.0)/8.0
#variable b equal 3.0*sqrt(3.0)/8.0
#variable c equal sqrt(3.0)/4.0
#lattice custom 2.5 a1 1.0 0.0 0.0 &
# a2 0.0 1.0 0.0 &
# a3 0.0 0.0 1.0 &
# basis 0.0 $a 0.0 &
# basis 0.25 $a 0.0 &
# basis 0.375 0.0 0.0 &
# basis 0.25 $b 0.0 &
# basis 0.5 $b 0.0 &
# basis 0.625 $c 0.0
#Defining a geometric region of space. Sets ID(user's choice), style(block, sphere, ...), then, args depends on the style chosen
#(for block, one has x0, xf, y0, yf, z0, zf, in distance units)
region box block 0.0 5.0 0.0 5.0 0.0 5.0
#Creating a simulation box based on the specified region. Entries: number of atom types and box ref.
create_box 1 box
#Creating atoms (or molecules) on a lattice, or a single atom (or molecule), ...
#Entries: atom type,
create_atoms 1 box
#Replicating NxNxN the entire set of atoms
#replicate 1 1 1
#######################
#######Settings########
#######################
#Setting one or more properties of one or more atoms.
#Setting mass
mass 1 1.0
#set group all mass 1.0
#Setting spins orientation and moment
#set group all spin/random 11 1.72
set group all spin 1.72 1.0 0.0 0.0
#Magnetic exchange interaction coefficient for bulk fcc Cobalt
pair_style pair/spin 4.0
#type i and j | interaction type | cutoff | J1 (eV) | J2 (adim) | J3 (Ang) (for Exchange)
pair_coeff * * exchange 4.0 0.0446928 0.003496 1.4885
#type i and j | interaction type | cutoff | Int (eV) | [dx,dy,dz] (for DMI and ME)
#pair_coeff * * dmi 2.6 0.001 0.0 0.0 1.0
#pair_coeff * * me 2.6 0.01 1.0 1.0 1.0
#Define a skin distance, update neigh list every
#neighbor 1.0 bin
#neigh_modify every 10 check yes delay 20
neighbor 0.0 bin
neigh_modify every 1 check no delay 0
#Magnetic field fix
#Type | Intensity (T or eV) | Direction
fix 1 all force/spin zeeman 1.0 0.0 0.0 1.0
#fix 1 all force/spin anisotropy 0.001 0.0 0.0 1.0
#fix 1 all force/spin anisotropy -0.1 0.0 0.0 1.0
#fix 1 all force/spin anisotropy 0.1 0.0 1.0 0.0
#Fix Langevin spins (merging damping and temperature)
#Temp | Alpha_trans | Alpha_long | Seed
#fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 2 all langevin/spin 1.0 0.1 0.0 21
#fix 2 all langevin/spin 0.0 0.0 0.0 21
#Magnetic integration fix
fix 3 all nve/spin mpi
#compute real time, total magnetization, magnetic energy, and spin temperature
#Iteration | Time | Mx | My | Mz | |M| | Em | Tm
compute mag all compute/spin
fix outmag all ave/time 1 1 10 c_mag[1] c_mag[2] c_mag[3] c_mag[4] c_mag[5] c_mag[6] c_mag[7] file mag_Co_nodamp.dat format %20.16g
#Setting the timestep for the simulation (in ps)
timestep 0.0001
##################
#######run########
##################
#Dump the positions and spin directions of magnetic particles (vmd format)
dump 1 all custom 100 dump_spin_T100.lammpstrj type x y z spx spy spz
#Running the simulations for N timesteps
run 2000
#run 1

126
examples/SPIN/in.cobalt_dev Normal file
View File

@ -0,0 +1,126 @@
###################
#######Init########
###################
clear
#setting units to metal (Ang, picosecs, eV, ...):
units metal
#setting dimension of the system (N=2 or 3):
dimension 3
#setting boundary conditions. (p for periodic, f for fixed, ...)
boundary p p p
#boundary f f f
#setting atom_style to spin:
atom_style spin
#Define sort for paramagnetic simulations (if no pair interaction)
#atom_modify sort 1000 4.0
#why?
atom_modify map array
#newton off for pair spin in SEQNEI
#newton off off
###########################
#######Create atoms########
###########################
#Lattice constant of fcc Cobalt
lattice fcc 3.54
#lattice sc 2.50
#Test Kagome
#variable a equal sqrt(3.0)/8.0
#variable b equal 3.0*sqrt(3.0)/8.0
#variable c equal sqrt(3.0)/4.0
#lattice custom 2.5 a1 1.0 0.0 0.0 &
# a2 0.0 1.0 0.0 &
# a3 0.0 0.0 1.0 &
# basis 0.0 $a 0.0 &
# basis 0.25 $a 0.0 &
# basis 0.375 0.0 0.0 &
# basis 0.25 $b 0.0 &
# basis 0.5 $b 0.0 &
# basis 0.625 $c 0.0
#Defining a geometric region of space. Sets ID(user's choice), style(block, sphere, ...), then, args depends on the style chosen
#(for block, one has x0, xf, y0, yf, z0, zf, in distance units)
region box block 0.0 5.0 0.0 5.0 0.0 5.0
#Creating a simulation box based on the specified region. Entries: number of atom types and box ref.
create_box 1 box
#Creating atoms (or molecules) on a lattice, or a single atom (or molecule), ...
#Entries: atom type,
create_atoms 1 box
#Replicating NxNxN the entire set of atoms
#replicate 1 1 1
#######################
#######Settings########
#######################
#Setting one or more properties of one or more atoms.
#Setting mass
mass 1 1.0
#set group all mass 1.0
#Setting spins orientation and moment
#set group all spin/random 11 1.72
set group all spin 1.72 1.0 0.0 0.0
#Magnetic exchange interaction coefficient for bulk fcc Cobalt
pair_style pair/spin 4.0
#type i and j | interaction type | cutoff | J1 (eV) | J2 (adim) | J3 (Ang) (for Exchange)
pair_coeff * * exchange 4.0 0.0446928 0.003496 1.4885
#type i and j | interaction type | cutoff | Int (eV) | [dx,dy,dz] (for DMI and ME)
#pair_coeff * * dmi 2.6 0.001 0.0 0.0 1.0
#pair_coeff * * me 2.6 0.01 1.0 1.0 1.0
#Define a skin distance, update neigh list every
#neighbor 1.0 bin
#neigh_modify every 10 check yes delay 20
neighbor 0.0 bin
neigh_modify every 1 check no delay 0
#Magnetic field fix
#Type | Intensity (T or eV) | Direction
fix 1 all force/spin zeeman 1.0 0.0 0.0 1.0
#fix 1 all force/spin anisotropy 0.001 0.0 0.0 1.0
#fix 1 all force/spin anisotropy -0.1 0.0 0.0 1.0
#fix 1 all force/spin anisotropy 0.1 0.0 1.0 0.0
#Fix Langevin spins (merging damping and temperature)
#Temp | Alpha_trans | Alpha_long | Seed
#fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 2 all langevin/spin 1.0 0.1 0.0 21
#fix 2 all langevin/spin 0.0 0.0 0.0 21
#Magnetic integration fix
fix 3 all nve/spin mpi
#compute real time, total magnetization, magnetic energy, and spin temperature
#Iteration | Time | Mx | My | Mz | |M| | Em | Tm
compute mag all compute/spin
fix outmag all ave/time 1 1 10 c_mag[1] c_mag[2] c_mag[3] c_mag[4] c_mag[5] c_mag[6] c_mag[7] file mag_Co_nodamp.dat format %20.16g
#Setting the timestep for the simulation (in ps)
timestep 0.0001
##################
#######run########
##################
#Dump the positions and spin directions of magnetic particles (vmd format)
dump 1 all custom 100 dump_spin_T100.lammpstrj type x y z spx spy spz
#Running the simulations for N timesteps
run 2000
#run 1

View File

@ -0,0 +1,79 @@
proc vmd_draw_arrow {mol start end} {
set middle [vecadd $start [vecscale 0.9 [vecsub $end $start]]]
graphics $mol cylinder $start $middle radius 0.05
graphics $mol cone $middle $end radius 0.01 color 3
}
proc vmd_draw_vector {args} {
set usage {"draw vector {x1 y1 z1} {x2 y2 z2} [scale <s>] [resolution <res>] [radius <r>] [filled <yes/no>]"}
# defaults
set scale 3.0
set res 5
set radius 0.1
set filled yes
if {[llength $args] < 3} {
error "wrong # args: should be $usage"
}
set mol [lindex $args 0]
set center [lindex $args 1]
set vector [lindex $args 2]
if {[llength $center] != 3 || [llength $vector] != 3} {
error "wrong type of args: should be $usage"
}
foreach {flag value} [lrange $args 3 end] {
switch -glob $flag {
scale {set scale $value}
res* {set res $value}
rad* {set radius $value}
fill* {set filled $value}
default {error "unknown option '$flag': should be $usage" }
}
}
set vechalf [vecscale [expr $scale * 0.5] $vector]
return [list \
[graphics $mol color yellow]\
[graphics $mol cylinder [vecsub $center $vechalf]\
[vecadd $center [vecscale 0.7 $vechalf]] \
radius $radius resolution $res filled $filled] \
[graphics $mol color orange]\
[graphics $mol cone [vecadd $center [vecscale 0.6 $vechalf]] \
[vecadd $center $vechalf] radius [expr $radius * 2.5] \
resolution $res]]
}
proc vmd_draw_spin {args} {
global molid
graphics $molid delete all
set frame [molinfo $molid get frame]
set natoms [molinfo $molid get numatoms]
for {set i 0} {$i < $natoms} {incr i} {
set sel [atomselect top "index $i"]
# set sel [atomselect top "index 1200"]
set coords [lindex [$sel get {x y z}] $molid]
set velocities [lindex [$sel get {vx vy vz}] $molid]
draw vector $coords $velocities
set uvx [lindex [$sel get {vx}] $molid]
set uvy [lindex [$sel get {vy}] $molid]
set uvz [lindex [$sel get {vz}] $molid]
$sel set user [vecadd [vecadd [vecscale $uvy $uvy] [vecscale $uvz $uvz] ] [vecscale $uvx $uvx]]
$sel set user $uvy
#draw vector $coords {0.0 uvy 0.0}
}
#pbc box -color 3
}
proc enable_trace {} {
global vmd_frame
trace variable vmd_frame([molinfo top]) w vmd_draw_spin
}
set molid [mol addfile {/home/jtranch/Documents/lammps/src/dump_spin_T100.lammpstrj} type {lammpstrj} autobonds off first 0 last -1 step 1 waitfor all]
scale by 0.5
animate style Loop
enable_trace