commit JT 020818

- adding files for doc and reorg.
This commit is contained in:
julient31
2018-02-08 11:15:01 -07:00
parent 10b38cda93
commit 3d18f55155
36 changed files with 19629 additions and 0 deletions

Binary file not shown.

After

Width:  |  Height:  |  Size: 60 KiB

View File

@ -0,0 +1,40 @@
\documentclass[preview]{standalone}
\usepackage{varwidth}
\usepackage[utf8x]{inputenc}
\usepackage{amsmath,amssymb,amsthm,bm,tikz}
\usetikzlibrary{automata,arrows,shapes,snakes}
\begin{document}
\begin{varwidth}{50in}
\begin{tikzpicture}
%Global
\node (v1) at (0,6.0) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] { $\bm{v} \leftarrow \bm{v}+L_v.\Delta t/2$ };
\node (s1) at (0,4.5) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] { $\bm{s} \leftarrow \bm{s}+L_s.\Delta t/2$ };
\node (r) at (0,3.0) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] { $\bm{r} \leftarrow \bm{r}+L_r.\Delta t$ };
\node (s2) at (0,1.5) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] { $\bm{s} \leftarrow \bm{s}+L_s.\Delta t/2$ };
\node (v2) at (0,0.0) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] { $\bm{v} \leftarrow \bm{v}+L_v.\Delta t/2$ };
\draw[line width=2pt, ->] (v1) -- (s1);
\draw[line width=2pt, ->] (s1) -- (r);
\draw[line width=2pt, ->] (r) -- (s2);
\draw[line width=2pt, ->] (s2) -- (v2);
%Spin
\node (s01) at (6,6.0) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] {$\bm{s}_0 \leftarrow \bm{s}_0+L_{s_0}.\Delta t/4$ };
\node (sN1) at (6,4.5) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] {$\bm{s}_{\rm N-1}\leftarrow\bm{s}_{\rm N-1}+L_{s_{\rm N-1}}.\Delta t/4$};
\node (sN) at (6,3.0) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] {$\bm{s}_{\rm N} \leftarrow \bm{s}_{\rm N}+L_{s_{\rm N}}.\Delta t/2$ };
\node (sN2) at (6,1.5) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] {$\bm{s}_{\rm N-1}\leftarrow\bm{s}_{\rm N-1}+L_{s_{\rm N-1}}.\Delta t/4$};
\node (s02) at (6,0.0) [draw,thick,minimum width=0.2cm,minimum height=0.2cm] {$\bm{s}_0 \leftarrow \bm{s}_0+L_{s_0}.\Delta t/4$ };
\draw[line width=2pt,dashed, ->] (s01) -- (sN1);
\draw[line width=2pt, ->] (sN1) -- (sN);
\draw[line width=2pt, ->] (sN) -- (sN2);
\draw[line width=2pt,dashed, ->] (sN2) -- (s02);
%from Global to Spin
\draw[line width=2pt, dashed, ->] (s1) -- (s01.west);
\draw[line width=2pt, dashed, ->] (s1) -- (s02.west);
\end{tikzpicture}
\end{varwidth}
\end{document}

Binary file not shown.

After

Width:  |  Height:  |  Size: 11 KiB

View File

@ -0,0 +1,14 @@
\documentclass[preview]{standalone}
\usepackage{varwidth}
\usepackage[utf8x]{inputenc}
\usepackage{amsmath, amssymb, graphics, setspace}
\begin{document}
\begin{varwidth}{50in}
\begin{equation}
\frac{d \vec{s}_{i}}{dt} = \frac{1}{\left(1+\lambda^2 \right)} \left( \left(
\vec{\omega}_{i} +\vec{\eta} \right) \times \vec{s}_{i} + \lambda\, \vec{s}_{i}
\times\left( \vec{\omega}_{i} \times\vec{s}_{i} \right) \right), \nonumber
\end{equation}
\end{varwidth}
\end{document}

Binary file not shown.

After

Width:  |  Height:  |  Size: 7.9 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 7.3 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 16 KiB

View File

@ -0,0 +1,14 @@
\documentclass[preview]{standalone}
\usepackage{varwidth}
\usepackage[utf8x]{inputenc}
\usepackage{amsmath,amssymb,amsthm,bm}
\begin{document}
\begin{varwidth}{50in}
\begin{equation}
\vec{F}^{i} = \sum_{j}^{Neighbor} \frac{\partial {J} \left(r_{ij} \right)}{
\partial r_{ij}} \left( \vec{s}_{i}\cdot \vec{s}_{j} \right) \vec{r}_{ij}
~~{\rm and}~~ \vec{\omega}^{i} = \frac{1}{\hbar} \sum_{j}^{Neighbor} {J}
\left(r_{ij} \right)\,\vec{s}_{j} \nonumber
\end{equation}
\end{varwidth}
\end{document}

Binary file not shown.

After

Width:  |  Height:  |  Size: 13 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 7.0 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 15 KiB

View File

@ -0,0 +1,13 @@
\documentclass[preview]{standalone}
\usepackage{varwidth}
\usepackage[utf8x]{inputenc}
\usepackage{amsmath,amssymb,amsthm,bm}
\begin{document}
\begin{varwidth}{50in}
\begin{equation}
\vec{F}^{i} = -\sum_{j}^{Neighbor} \left( \vec{s}_{i}\times \vec{s}_{j} \right)
\times \vec{E} ~~{\rm and}~~ \vec{\omega}^{i} = -\frac{1}{\hbar}
\sum_{j}^{Neighbor} \vec{s}_j \times \left(\vec{E}\times r_{ij} \right),\nonumber
\end{equation}
\end{varwidth}
\end{document}

Binary file not shown.

After

Width:  |  Height:  |  Size: 17 KiB

View File

@ -0,0 +1,12 @@
\documentclass[preview]{standalone}
\usepackage{varwidth}
\usepackage[utf8x]{inputenc}
\usepackage{amsmath,amssymb,amsthm,bm}
\begin{document}
\begin{varwidth}{50in}
\begin{equation}
\vec{\omega}_i = -\frac{1}{\hbar} \sum_{j}^{Neighb} \vec{s}_{j}\times\vec{D}(r_{ij}) ~~{\rm and}~~
\vec{F}_i = -\sum_{j}^{Neighb} \frac{\partial D(r_{ij})}{\partial r_{ij}} \left(\vec{s}_{i}\times \vec{s}_{j} \right) \cdot \vec{r}_{ij}, \nonumber
\end{equation}
\end{varwidth}
\end{document}

Binary file not shown.

After

Width:  |  Height:  |  Size: 7.9 KiB

View File

@ -0,0 +1,11 @@
\documentclass[preview]{standalone}
\usepackage{varwidth}
\usepackage[utf8x]{inputenc}
\usepackage{amsmath,amssymb,amsthm,bm}
\begin{document}
\begin{varwidth}{50in}
\begin{equation}
\bm{H}_{dm} = \sum_{{ i,j}=1,i\neq j}^{N} \vec{D}\left(r_{ij}\right) \cdot\left(\vec{s}_{i}\times \vec{s}_{j}\right), \nonumber
\end{equation}
\end{varwidth}
\end{document}

73
doc/src/pair_spin_me.txt Normal file
View File

@ -0,0 +1,73 @@
"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
pair_style pair/spin/me command :h3
[Syntax:]
pair_style pair/spin/me cutoff :pre
cutoff = global cutoff pair (distance in metal units) :ulb,l
:ule
[Examples:]
pair_style pair/spin/me 4.5
pair_coeff * * pair/spin/me me 4.5 0.00109 1.0 1.0 1.0 :pre
[Description:]
Style {pair/spin/me} computes a magneto-electric interaction between
pairs of magnetic spins. According to the derivation reported in
"(Katsura)"_#Katsura1, this interaction is defined as:
:c,image(Eqs/pair_spin_me_interaction.jpg)
where si and sj are neighboring magnetic spins of two particles,
rij = ri - rj is the normalized separation vector between the
two particles, and E is an electric polarization vector.
The norm and direction of E are giving the intensity and the
direction of a screened dielectric atomic polarization (in eV).
From this magneto-electric interaction, each spin i will be submitted
to a magnetic torque omega, and its associated atom can be submitted to a
force F for spin-lattice calculations (see
"fix_integration_spin"_fix_integration_spin.html), such as:
:c,image(Eqs/pair_spin_me_forces.jpg)
with h the Planck constant (in metal units).
More details about the derivation of these torques/forces are reported in
"(Tranchida)"_#Tranchida1.
:line
[Restrictions:]
All the {pair/spin} styles are part of the SPIN package.
These styles are only enabled if LAMMPS was built with this package, and
if the atom_style "spin" was declared.
See the "Making LAMMPS"_Section_start.html#start_3 section for more info.
[Related commands:]
"atom_style spin"_atom_style.html, "pair_coeff"_pair_coeff.html,
"pair_spin_exchange"_pair_spin_exchange.html, "pair_eam"_pair_eam.html,
[Default:] none
:line
:link(Katsura1)
[(Katsura)] H. Katsura, N. Nagaosa, A.V. Balatsky. Phys. Rev. Lett., 95(5), 057205. (2005)
:link(Tranchida1)
[(Tranchida)] J. Tranchida, S.J. Plimpton, P. Thibaudeau, A.P. Thompson.
arXiv preprint arXiv:1801.10233. (2018)

View File

@ -0,0 +1,82 @@
"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
pair_style pair/spin/soc/dmi command :h3
[Syntax:]
pair_style pair/spin/soc/dmi cutoff :pre
cutoff = global cutoff pair (distance in metal units) :ulb,l
:ule
[Examples:]
pair_style pair/spin/soc/dmi 4.0
pair_coeff * * dmi 2.6 0.001 1.0 0.0 0.0
pair_coeff 1 2 dmi 4.0 0.00109 0.0 0.0 1.0 :pre
[Description:]
Style {pair/spin/soc/dmi} computes the Dzyaloshinskii-Moriya (DM) interaction
between pairs of magnetic spins:
:c,image(Eqs/pair_spin_soc_dmi_interaction.jpg)
where si and sj are two neighboring magnetic spins of two particles,
rij = ri - rj is the inter-atomic distance between the two particles,
and D(rij) is the DM vector defining the intensity and the
sign of the exchange interaction.
This function is defined as:
:c,image(Eqs/pair_spin_exchange_function.jpg)
where a, b and d are the three constant coefficients defined in the associated
"pair_coeff" command.
The coefficients a, b, and c need to be fitted so that the function above matches with
the value of the DM interaction for the N neighbor shells taken
into account.
Examples and more explanations about this function and its parametrization are reported
in "(Tranchida)"_#Tranchida1.
From this DM interaction, each spin i will be submitted to a magnetic torque
omega and its associated atom to a force F (for spin-lattice calculations only),
such as:
:c,image(Eqs/pair_spin_soc_dmi_forces.jpg)
with h the Planck constant (in metal units).
More details about the derivation of these torques/forces are reported in
"(Tranchida)"_#Tranchida1.
:line
[Restrictions:]
All the {pair/spin} styles are part of the SPIN package.
These styles are only enabled if LAMMPS was built with this package, and
if the atom_style "spin" was declared.
See the "Making LAMMPS"_Section_start.html#start_3 section for more info.
[Related commands:]
"atom_style spin"_atom_style.html, "pair_coeff"_pair_coeff.html,
"pair_eam"_pair_eam.html,
[Default:] none
:line
:link(Tranchida1)
[(Tranchida)]https://arxiv.org/abs/1801.10233

View File

@ -0,0 +1,55 @@
# sc iron atoms in bismuth oxide
clear
units metal
atom_style spin
dimension 3
boundary p p p
# check why?
atom_modify map array
lattice sc 3.96
region box block 0.0 34.0 0.0 34.0 0.0 5.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for bfo
mass 1 1.0
set group all spin/random 11 2.50
pair_style hybrid/overlay pair/spin/exchange 6.0 pair/spin/me 4.5
pair_coeff * * pair/spin/exchange exchange 6.0 -0.01575 0.0 1.965
pair_coeff * * pair/spin/me me 4.5 0.000109 1.0 1.0 1.0
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all force/spin anisotropy 0.0000035 0.0 0.0 1.0
fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 3 all integration/spin lattice no
timestep 0.0001
compute out_mag all compute/spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
thermo_style custom step time v_magnorm v_emag temp etotal
thermo 10
#dump 1 all custom 100 dump_spin_BFO.lammpstrj type x y z spx spy spz
run 10000

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,59 @@
# fcc cobalt in a 3d periodic box
clear
units metal
atom_style spin
dimension 3
boundary p p p
# check why?
atom_modify map array
lattice fcc 3.54
region box block 0.0 5.0 0.0 5.0 0.0 5.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for cobalt
mass 1 58.93
set group all spin/random 31 1.72
velocity all create 100 4928459 rot yes dist gaussian
#pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0 pair/spin/soc/neel 4.0
pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0
pair_coeff * * eam/alloy ../examples/SPIN/cobalt/Co_PurjaPun_2012.eam.alloy Co
pair_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
#pair_coeff * * pair/spin/soc/neel neel 4.0 0.003330282 0.864159 2.13731
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all force/spin zeeman 0.1 0.0 0.0 1.0
fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 3 all integration/spin lattice yes
timestep 0.0001
compute out_mag all compute/spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
thermo_style custom step time v_magnorm v_emag temp etotal
thermo 10
#dump 1 all custom 50 dump_cobalt.lammpstrj type x y z spx spy spz
run 1000

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,95 @@
#!/usr/bin/env python
# This file reads in the file log.lammps generated by the script ELASTIC/in.elastic
# It prints out the 6x6 tensor of elastic constants Cij
# followed by the 6x6 tensor of compliance constants Sij
# It uses the same conventions as described in:
# Sprik, Impey and Klein PRB (1984).
# The units of Cij are whatever was used in log.lammps (usually GPa)
# The units of Sij are the inverse of that (usually 1/GPa)
from numpy import zeros
from numpy.linalg import inv
# define logfile layout
nvals = 21
valpos = 4
valstr = '\nElastic Constant C'
# define order of Cij in logfile
cindices = [0]*nvals
cindices[0] = (0,0)
cindices[1] = (1,1)
cindices[2] = (2,2)
cindices[3] = (0,1)
cindices[4] = (0,2)
cindices[5] = (1,2)
cindices[6] = (3,3)
cindices[7] = (4,4)
cindices[8] = (5,5)
cindices[9] = (0,3)
cindices[10] = (0,4)
cindices[11] = (0,5)
cindices[12] = (1,3)
cindices[13] = (1,4)
cindices[14] = (1,5)
cindices[15] = (2,3)
cindices[16] = (2,4)
cindices[17] = (2,5)
cindices[18] = (3,4)
cindices[19] = (3,5)
cindices[20] = (4,5)
# open logfile
logfile = open("log.lammps",'r')
txt = logfile.read()
# search for 21 elastic constants
c = zeros((6,6))
s2 = 0
for ival in range(nvals):
s1 = txt.find(valstr,s2)
if (s1 == -1):
print "Failed to find elastic constants in log file"
exit(1)
s1 += 1
s2 = txt.find("\n",s1)
line = txt[s1:s2]
# print line
words = line.split()
(i1,i2) = cindices[ival]
c[i1,i2] = float(words[valpos])
c[i2,i1] = c[i1,i2]
print "C tensor [GPa]"
for i in range(6):
for j in range(6):
print "%10.8g " % c[i][j],
print
# apply factor of 2 to columns of off-diagonal elements
for i in range(6):
for j in range(3,6):
c[i][j] *= 2.0
s = inv(c)
# apply factor of 1/2 to columns of off-diagonal elements
for i in range(6):
for j in range(3,6):
s[i][j] *= 0.5
print "S tensor [1/GPa]"
for i in range(6):
for j in range(6):
print "%10.8g " % s[i][j],
print

View File

@ -0,0 +1,142 @@
# NOTE: This script should not need to be
# modified. See in.elastic for more info.
#
# Find which reference length to use
if "${dir} == 1" then &
"variable len0 equal ${lx0}"
if "${dir} == 2" then &
"variable len0 equal ${ly0}"
if "${dir} == 3" then &
"variable len0 equal ${lz0}"
if "${dir} == 4" then &
"variable len0 equal ${lz0}"
if "${dir} == 5" then &
"variable len0 equal ${lz0}"
if "${dir} == 6" then &
"variable len0 equal ${ly0}"
# Reset box and simulation parameters
clear
box tilt large
read_restart restart.equil
include potential.mod
# Negative deformation
variable delta equal -${up}*${len0}
variable deltaxy equal -${up}*xy
variable deltaxz equal -${up}*xz
variable deltayz equal -${up}*yz
if "${dir} == 1" then &
"change_box all x delta 0 ${delta} xy delta ${deltaxy} xz delta ${deltaxz} remap units box"
if "${dir} == 2" then &
"change_box all y delta 0 ${delta} yz delta ${deltayz} remap units box"
if "${dir} == 3" then &
"change_box all z delta 0 ${delta} remap units box"
if "${dir} == 4" then &
"change_box all yz delta ${delta} remap units box"
if "${dir} == 5" then &
"change_box all xz delta ${delta} remap units box"
if "${dir} == 6" then &
"change_box all xy delta ${delta} remap units box"
# Relax atoms positions
minimize ${etol} ${ftol} ${maxiter} ${maxeval}
# Obtain new stress tensor
variable tmp equal pxx
variable pxx1 equal ${tmp}
variable tmp equal pyy
variable pyy1 equal ${tmp}
variable tmp equal pzz
variable pzz1 equal ${tmp}
variable tmp equal pxy
variable pxy1 equal ${tmp}
variable tmp equal pxz
variable pxz1 equal ${tmp}
variable tmp equal pyz
variable pyz1 equal ${tmp}
# Compute elastic constant from pressure tensor
variable C1neg equal ${d1}
variable C2neg equal ${d2}
variable C3neg equal ${d3}
variable C4neg equal ${d4}
variable C5neg equal ${d5}
variable C6neg equal ${d6}
# Reset box and simulation parameters
clear
box tilt large
read_restart restart.equil
include potential.mod
# Positive deformation
variable delta equal ${up}*${len0}
variable deltaxy equal ${up}*xy
variable deltaxz equal ${up}*xz
variable deltayz equal ${up}*yz
if "${dir} == 1" then &
"change_box all x delta 0 ${delta} xy delta ${deltaxy} xz delta ${deltaxz} remap units box"
if "${dir} == 2" then &
"change_box all y delta 0 ${delta} yz delta ${deltayz} remap units box"
if "${dir} == 3" then &
"change_box all z delta 0 ${delta} remap units box"
if "${dir} == 4" then &
"change_box all yz delta ${delta} remap units box"
if "${dir} == 5" then &
"change_box all xz delta ${delta} remap units box"
if "${dir} == 6" then &
"change_box all xy delta ${delta} remap units box"
# Relax atoms positions
minimize ${etol} ${ftol} ${maxiter} ${maxeval}
# Obtain new stress tensor
variable tmp equal pe
variable e1 equal ${tmp}
variable tmp equal press
variable p1 equal ${tmp}
variable tmp equal pxx
variable pxx1 equal ${tmp}
variable tmp equal pyy
variable pyy1 equal ${tmp}
variable tmp equal pzz
variable pzz1 equal ${tmp}
variable tmp equal pxy
variable pxy1 equal ${tmp}
variable tmp equal pxz
variable pxz1 equal ${tmp}
variable tmp equal pyz
variable pyz1 equal ${tmp}
# Compute elastic constant from pressure tensor
variable C1pos equal ${d1}
variable C2pos equal ${d2}
variable C3pos equal ${d3}
variable C4pos equal ${d4}
variable C5pos equal ${d5}
variable C6pos equal ${d6}
# Combine positive and negative
variable C1${dir} equal 0.5*(${C1neg}+${C1pos})
variable C2${dir} equal 0.5*(${C2neg}+${C2pos})
variable C3${dir} equal 0.5*(${C3neg}+${C3pos})
variable C4${dir} equal 0.5*(${C4neg}+${C4pos})
variable C5${dir} equal 0.5*(${C5neg}+${C5pos})
variable C6${dir} equal 0.5*(${C6neg}+${C6pos})
# Delete dir to make sure it is not reused
variable dir delete

View File

@ -0,0 +1,59 @@
# fcc cobalt in a 3d periodic box
clear
units metal
atom_style spin
dimension 3
boundary p p p
# check why?
atom_modify map array
lattice fcc 3.54
region box block 0.0 5.0 0.0 5.0 0.0 5.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for cobalt
mass 1 58.93
set group all spin/random 31 1.72
velocity all create 100 4928459 rot yes dist gaussian
#pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0 pair/spin/soc/neel 4.0
pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0
pair_coeff * * eam/alloy ../examples/SPIN/cobalt/Co_PurjaPun_2012.eam.alloy Co
pair_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
#pair_coeff * * pair/spin/soc/neel neel 4.0 0.003330282 0.864159 2.13731
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all force/spin zeeman 0.1 0.0 0.0 1.0
fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 3 all integration/spin lattice yes
timestep 0.0001
compute out_mag all compute/spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
thermo_style custom step time v_magnorm v_emag temp etotal
thermo 10
#dump 1 all custom 50 dump_cobalt.lammpstrj type x y z spx spy spz
run 1000

View File

@ -0,0 +1,204 @@
# Compute elastic constant tensor for a crystal
#
# Written by Aidan Thompson (Sandia, athomps@sandia.gov)
#
# This script uses the following three include files.
#
# init.mod (must be modified for different crystal structures)
# Define units, deformation parameters and initial
# configuration of the atoms and simulation cell.
#
#
# potential.mod (must be modified for different pair styles)
# Define pair style and other attributes
# not stored in restart file
#
#
# displace.mod (displace.mod should not need to be modified)
# Perform positive and negative box displacements
# in direction ${dir} and size ${up}.
# It uses the resultant changes
# in stress to compute one
# row of the elastic stiffness tensor
#
# Inputs variables:
# dir = the Voigt deformation component
# (1,2,3,4,5,6)
# Global constants:
# up = the deformation magnitude (strain units)
# cfac = conversion from LAMMPS pressure units to
# output units for elastic constants
#
#
# To run this on a different system, it should only be necessary to
# modify the files init.mod and potential.mod. In order to calculate
# the elastic constants correctly, care must be taken to specify
# the correct units in init.mod (units, cfac and cunits). It is also
# important to verify that the minimization of energy w.r.t atom
# positions in the deformed cell is fully converged.
# One indication of this is that the elastic constants are insensitive
# to the choice of the variable ${up} in init.mod. Another is to check
# the final max and two-norm forces reported in the log file. If you know
# that minimization is not required, you can set maxiter = 0.0 in
# init.mod.
#
include init.mod
include potential.mod
# Compute initial state
fix 3 all box/relax aniso 0.0
minimize ${etol} ${ftol} ${maxiter} ${maxeval}
variable tmp equal pxx
variable pxx0 equal ${tmp}
variable tmp equal pyy
variable pyy0 equal ${tmp}
variable tmp equal pzz
variable pzz0 equal ${tmp}
variable tmp equal pyz
variable pyz0 equal ${tmp}
variable tmp equal pxz
variable pxz0 equal ${tmp}
variable tmp equal pxy
variable pxy0 equal ${tmp}
variable tmp equal lx
variable lx0 equal ${tmp}
variable tmp equal ly
variable ly0 equal ${tmp}
variable tmp equal lz
variable lz0 equal ${tmp}
# These formulas define the derivatives w.r.t. strain components
# Constants uses $, variables use v_
variable d1 equal -(v_pxx1-${pxx0})/(v_delta/v_len0)*${cfac}
variable d2 equal -(v_pyy1-${pyy0})/(v_delta/v_len0)*${cfac}
variable d3 equal -(v_pzz1-${pzz0})/(v_delta/v_len0)*${cfac}
variable d4 equal -(v_pyz1-${pyz0})/(v_delta/v_len0)*${cfac}
variable d5 equal -(v_pxz1-${pxz0})/(v_delta/v_len0)*${cfac}
variable d6 equal -(v_pxy1-${pxy0})/(v_delta/v_len0)*${cfac}
displace_atoms all random ${atomjiggle} ${atomjiggle} ${atomjiggle} 87287 units box
# Write restart
unfix 3
write_restart restart.equil
# uxx Perturbation
variable dir equal 1
include displace.mod
# uyy Perturbation
variable dir equal 2
include displace.mod
# uzz Perturbation
variable dir equal 3
include displace.mod
# uyz Perturbation
variable dir equal 4
include displace.mod
# uxz Perturbation
variable dir equal 5
include displace.mod
# uxy Perturbation
variable dir equal 6
include displace.mod
# Output final values
variable C11all equal ${C11}
variable C22all equal ${C22}
variable C33all equal ${C33}
variable C12all equal 0.5*(${C12}+${C21})
variable C13all equal 0.5*(${C13}+${C31})
variable C23all equal 0.5*(${C23}+${C32})
variable C44all equal ${C44}
variable C55all equal ${C55}
variable C66all equal ${C66}
variable C14all equal 0.5*(${C14}+${C41})
variable C15all equal 0.5*(${C15}+${C51})
variable C16all equal 0.5*(${C16}+${C61})
variable C24all equal 0.5*(${C24}+${C42})
variable C25all equal 0.5*(${C25}+${C52})
variable C26all equal 0.5*(${C26}+${C62})
variable C34all equal 0.5*(${C34}+${C43})
variable C35all equal 0.5*(${C35}+${C53})
variable C36all equal 0.5*(${C36}+${C63})
variable C45all equal 0.5*(${C45}+${C54})
variable C46all equal 0.5*(${C46}+${C64})
variable C56all equal 0.5*(${C56}+${C65})
# Average moduli for cubic crystals
variable C11cubic equal (${C11all}+${C22all}+${C33all})/3.0
variable C12cubic equal (${C12all}+${C13all}+${C23all})/3.0
variable C44cubic equal (${C44all}+${C55all}+${C66all})/3.0
variable bulkmodulus equal (${C11cubic}+2*${C12cubic})/3.0
variable shearmodulus1 equal ${C44cubic}
variable shearmodulus2 equal (${C11cubic}-${C12cubic})/2.0
variable poissonratio equal 1.0/(1.0+${C11cubic}/${C12cubic})
# For Stillinger-Weber silicon, the analytical results
# are known to be (E. R. Cowley, 1988):
# C11 = 151.4 GPa
# C12 = 76.4 GPa
# C44 = 56.4 GPa
print "========================================="
print "Components of the Elastic Constant Tensor"
print "========================================="
print "Elastic Constant C11all = ${C11all} ${cunits}"
print "Elastic Constant C22all = ${C22all} ${cunits}"
print "Elastic Constant C33all = ${C33all} ${cunits}"
print "Elastic Constant C12all = ${C12all} ${cunits}"
print "Elastic Constant C13all = ${C13all} ${cunits}"
print "Elastic Constant C23all = ${C23all} ${cunits}"
print "Elastic Constant C44all = ${C44all} ${cunits}"
print "Elastic Constant C55all = ${C55all} ${cunits}"
print "Elastic Constant C66all = ${C66all} ${cunits}"
print "Elastic Constant C14all = ${C14all} ${cunits}"
print "Elastic Constant C15all = ${C15all} ${cunits}"
print "Elastic Constant C16all = ${C16all} ${cunits}"
print "Elastic Constant C24all = ${C24all} ${cunits}"
print "Elastic Constant C25all = ${C25all} ${cunits}"
print "Elastic Constant C26all = ${C26all} ${cunits}"
print "Elastic Constant C34all = ${C34all} ${cunits}"
print "Elastic Constant C35all = ${C35all} ${cunits}"
print "Elastic Constant C36all = ${C36all} ${cunits}"
print "Elastic Constant C45all = ${C45all} ${cunits}"
print "Elastic Constant C46all = ${C46all} ${cunits}"
print "Elastic Constant C56all = ${C56all} ${cunits}"
print "========================================="
print "Average properties for a cubic crystal"
print "========================================="
print "Bulk Modulus = ${bulkmodulus} ${cunits}"
print "Shear Modulus 1 = ${shearmodulus1} ${cunits}"
print "Shear Modulus 2 = ${shearmodulus2} ${cunits}"
print "Poisson Ratio = ${poissonratio}"

View File

@ -0,0 +1,55 @@
# NOTE: This script can be modified for different atomic structures,
# units, etc. See in.spin.curie_temperature for more info.
#
# Define the finite deformation size. Try several values of this
# variable to verify that results do not depend on it.
variable up equal 1.0e-6
# Define the amount of random jiggle for atoms
# This prevents atoms from staying on saddle points
variable atomjiggle equal 1.0e-5
# Uncomment one of these blocks, depending on what units
# you are using in LAMMPS and for output
# metal units, elastic constants in eV/A^3
#units metal
#variable cfac equal 6.2414e-7
#variable cunits string eV/A^3
# metal units, elastic constants in GPa
units metal
variable cfac equal 1.0e-4
variable cunits string GPa
# real units, elastic constants in GPa
#units real
#variable cfac equal 1.01325e-4
#variable cunits string GPa
# Define minimization parameters
variable etol equal 0.0
variable ftol equal 1.0e-10
variable maxiter equal 100
variable maxeval equal 1000
variable dmax equal 1.0e-2
# generate the box and atom positions using a diamond lattice
variable a equal 5.43
boundary p p p
lattice fcc 3.54
region box block 0.0 5.0 0.0 5.0 0.0 5.0
create_box 1 box
create_atoms 1 box
#lattice diamond $a
#region box prism 0 2.0 0 3.0 0 4.0 0.0 0.0 0.0
#create_box 1 box
#create_atoms 1 box
# Need to set mass (even for fixed lattice calculation, just to satisfy LAMMPS)
mass 1 58.93

View File

@ -0,0 +1,30 @@
# NOTE: This script can be modified for different pair styles
# See in.elastic for more info.
# Choose potentials (spin or spin-lattice)
pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0
pair_coeff * * eam/alloy ../examples/SPIN/cobalt/Co_PurjaPun_2012.eam.alloy Co
air_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
# Setup neighbor style
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
# Setup minimization style
min_style cg
min_modify dmax ${dmax} line quadratic
# Setup output
compute out_mag all compute/spin
compute out_temp all temp
variable magnorm equal c_out_mag[5]
variable tmag equal c_out_mag[7]
thermo_style custom step time v_magnorm v_tmag temp
thermo 1
#thermo 1
#thermo_style custom step temp pe press pxx pyy pzz pxy pxz pyz lx ly lz vol
#thermo_modify norm no

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,125 @@
###################
#######Init########
###################
clear
units metal
dimension 3
#boundary p p p
boundary p p p
#newton off
#setting atom_style to spin:
atom_style spin
#Define sort for paramagnetic simulations (if no pair interaction)
#atom_modify sort 1000 4.0
atom_modify map array
###########################
#######Create atoms########
###########################
lattice fcc 3.54
region box block 0.0 5.0 0.0 5.0 0.0 5.0
create_box 1 box
create_atoms 1 box
#######################
#######Settings########
#######################
#isolating 1 atom into a group
group single_spin id 10
#Setting one or more properties of one or more atoms.
mass 1 58.93
#Setting spins orientation and moment
set group all spin/random 31 1.72
#set group all spin 1.72 0.0 0.0 1.0
#set group single_spin spin/random 11 1.72
#velocity all create 200 4928459 rot yes dist gaussian
velocity all create 200 4928459 rot yes dist gaussian
#Magneto-mechanic interactions for bulk fcc Cobalt
#pair_style pair/spin/exchange 4.0
#pair_style eam/alloy
#pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0 pair/spin/soc/neel 4.0
pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0
pair_coeff * * eam/alloy ../examples/SPIN/Co_PurjaPun_2012.eam.alloy Co
#pair_coeff * * ../Co_PurjaPun_2012.eam.alloy Co
#pair_style pair/spin/exchange 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
#pair_coeff * * exchange 4.0 0.0 0.003496 1.4885
pair_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
#pair_coeff * * pair/spin/exchange exchange 2.5 0.0446928 0.003496 1.4885
#pair_coeff * * pair/spin/exchange exchange 4.0 0.0 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
#type i and j | interaction type | cutoff | K1 (eV) | K2 (adim) | K3 (Ang) (for SOC)
#pair_coeff * * pair/spin/soc/neel neel 4.0 0.003330282 0.864159 2.13731
#pair_coeff * * pair/spin/soc/neel neel 4.0 0.0 0.864159 2.13731
#Define a skin distance, update neigh list every
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
#neighbor 1.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 0.0 0.0 0.0 1.0
#fix 1 all force/spin anisotropy 0.01 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.0 0.0 21
#Magnetic integration fix
fix 3 all integration/spin lattice yes
#fix 3 all integration/spin lattice no
#compute real time, total magnetization, magnetic energy, and spin temperature
#Iteration | Time | Mx | My | Mz | |M| | Em | Tm
#Setting the timestep for the simulation (in ps)
timestep 0.0001
##################
#######run########
##################
compute out_mag all compute/spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
#variable test equal etotal-0.5*c_out_mag[6]
thermo 10
thermo_style custom step time v_magnorm v_emag temp etotal
thermo_modify format float %20.15g
#Dump the positions and spin directions of magnetic particles (vmd format)
dump 1 all custom 50 dump_cobalt.lammpstrj type x y z spx spy spz
#Running the simulations for N timesteps
run 10
#run 1000

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,88 @@
###################
#######Init########
###################
clear
units metal
dimension 3
boundary p p p
#boundary f f f
#setting atom_style to spin:
atom_style spin
atom_modify map array
###########################
#######Create atoms########
###########################
read_data ../examples/SPIN/Norm_randXY_8x8x32_test.data
#######################
#######Settings########
#######################
#Setting one or more properties of one or more atoms.
mass 1 58.93
#velocity all create 200 4928459 rot yes dist gaussian
#Magneto-mechanic interactions for bulk fcc Cobalt
pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0 pair/spin/soc/neel 4.0
# cobalt eam potential
pair_coeff * * eam/alloy ../examples/SPIN/Co_PurjaPun_2012.eam.alloy Co
#type i and j | interaction type | cutoff | J1 (eV) | J2 (adim) | J3 (Ang) (for Exchange)
pair_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
#type i and j | interaction type | cutoff | K1 (eV) | K2 (adim) | K3 (Ang) (for SOC)
pair_coeff * * pair/spin/soc/neel neel 4.0 0.003330282 0.864159 2.13731
#Define a skin distance, update neigh list every
#neighbor 1.0 bin
#neigh_modify every 10 check yes delay 20
neighbor 1.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 0.0 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.0 0.0 21
#Magnetic integration fix
fix 3 all integration/spin serial
#Setting the timestep for the simulation (in ps)
timestep 0.0001
##################
#######run########
##################
compute out_mag all compute/spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
thermo 10
thermo_style custom step time v_magnorm v_emag v_tmag temp etotal
thermo_modify format float %20.15g
#Dump the positions and spin directions of magnetic particles (vmd format)
dump 1 all custom 1 dump.lammpstrj type x y z spx spy spz
#Running the simulations for N timesteps
run 10
#run 10000

View File

@ -0,0 +1,89 @@
# read a dump file for a sim. of magnetic cobalt
clear
units metal
atom_style spin
dimension 3
boundary p p p
# check why?
atom_modify map array
lattice fcc 3.54
region box block 0.0 5.0 0.0 5.0 0.0 5.0
create_box 1 box
read_dump ../examples/SPIN/Norm_randXY_8x8x32.dump 0 x y z box yes
create_atoms 1 box
# setting mass, mag. moments, and interactions
mass 1 58.93
#Setting spins orientation and moment
#set group all spin/random 31 1.72
set group all spin 1.72 0.0 0.0 1.0
set group single_spin spin/random 11 1.72
velocity all create 200 4928459 rot yes dist gaussian
#Magneto-mechanic interactions for bulk fcc Cobalt
pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0 pair/spin/soc/neel 4.0
# cobalt eam potential
pair_coeff * * eam/alloy ../examples/SPIN/Co_PurjaPun_2012.eam.alloy Co
#type i and j | interaction type | cutoff | J1 (eV) | J2 (adim) | J3 (Ang) (for Exchange)
pair_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
#type i and j | interaction type | cutoff | K1 (eV) | K2 (adim) | K3 (Ang) (for SOC)
pair_coeff * * pair/spin/soc/neel neel 4.0 0.003330282 0.864159 2.13731
#Define a skin distance, update neigh list every
#neighbor 1.0 bin
#neigh_modify every 10 check yes delay 20
neighbor 1.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 0.0 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.0 0.0 21
#Magnetic integration fix
fix 3 all integration/spin serial
#Setting the timestep for the simulation (in ps)
timestep 0.0001
##################
#######run########
##################
compute out_mag all compute/spin
compute out_pe all pe
compute out_ke all ke
compute out_temp all temp
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
thermo 10
thermo_style custom step time v_magnorm v_emag v_tmag temp etotal
thermo_modify format float %20.15g
#Dump the positions and spin directions of magnetic particles (vmd format)
dump 1 all custom 20 dump.lammpstrj type x y z spx spy spz
#Running the simulations for N timesteps
run 100
#run 10000

View File

@ -0,0 +1,55 @@
# fcc cobalt in a 3d box
clear
units metal
atom_style spin
dimension 3
boundary p p f
# check why?
atom_modify map array
lattice fcc 3.54
region box block 0.0 50.0 0.0 50.0 0.0 4.0
create_box 1 box
create_atoms 1 box
# setting mass, mag. moments, and interactions for cobalt
mass 1 58.93
set group all spin/random 31 1.72
#pair_style hybrid/overlay eam/alloy pair/spin/exchange 4.0 pair/spin/soc/neel 4.0
pair_style hybrid/overlay pair/spin/exchange 4.0 pair/spin/soc/dmi 2.6
pair_coeff * * pair/spin/exchange exchange 4.0 0.0446928 0.003496 1.4885
pair_coeff * * pair/spin/soc/dmi dmi 2.6 0.01 0.0 0.0 1.0
#pair_coeff * * pair/spin/soc/neel neel 4.0 0.003330282 0.864159 2.13731
neighbor 0.1 bin
neigh_modify every 10 check yes delay 20
fix 1 all force/spin anisotropy 0.0001 0.0 0.0 1.0
fix 2 all langevin/spin 0.0 0.1 0.0 21
fix 3 all integration/spin lattice no
timestep 0.0002
# define output and run
compute out_mag all compute/spin
variable magz equal c_out_mag[4]
variable magnorm equal c_out_mag[5]
variable emag equal c_out_mag[6]
variable tmag equal c_out_mag[7]
variable mag_force equal f_1
thermo_style custom step time v_magnorm v_emag etotal
thermo 50
dump 1 all custom 50 dump_skyrmion.lammpstrj type x y z spx spy spz
run 10000

View File

@ -0,0 +1,91 @@
#!/bin/bash
# example prepare_vmd.sh /home/jtranch/Documents/lammps/src/dump_VSRSV.lammpstrj
# you will get a return file
echo "vmd script for file $1 is preparing..."
timestamp(){
date +%s
}
TS=$(timestamp)
FILE=view_${TS}.vmd
cat >${FILE} <<EOF
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 2.0
set res 50
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 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 {$1} type {lammpstrj} autobonds off first 0 last -1 step 1 waitfor all]
scale by 0.5
animate style Loop
enable_trace
EOF
echo "$FILE is ready..."

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 2.0
set res 50
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.lammpstrj} type {lammpstrj} autobonds off first 0 last -1 step 1 waitfor all]
scale by 0.5
animate style Loop
enable_trace