Transcript Molecular Dynamics - University of Washington
Molecular Dynamics Simulations Joo Chul Yoon with Prof. Scott T. Dunham Electrical Engineering University of Washington
Contents Introduction to MD Simulation Setup Integration Method Force Calculation and MD Potential MD Simulations of Silicon Recrystallization Simulation Preparation SW Potential Tersoff Potential
Introduction to Molecular Dynamics Calculate how a system of particles evolves in time Consider a set of atoms with positions /velocities and the potential energy function of the system Predict the next positions of particles over some short time interval by solving Newtonian mechanics
Basic MD Algorithm
v
i
(
t
0 ) Get new forces
F
i
(
r
i
) Solve the equations of motion numerically over a short step
r
i
v
i
(
t
(
t
) )
r
i
v
i
(
t
(
t
t
t
) )
t t
t
t t
t
Calculate results and finish
Simulation Setup Simulation Cell Boundary Condition Constructing neighboring cells Initial atom velocities MD Time step Temperature Control
Simulation Cell usually using orthogonal cells Open boundary for a molecule or nanocluster in vacuum not for a continuous medium Fixed boundary fixed boundary atoms completely unphysical Periodic boundary conditions obtaining bulk properties
Periodic boundary conditions An atom moving out of boundary comes back on the other side considered in force calculation
r cut
L
2
r cut
Constructing neighboring cells pair potential calculation .
2
O
( N 2 ) not necessary to search all atoms Verlet neighbor list containing all neighbor atoms within
r L
updating every time steps
L
where
r L
r cut
N L
v
t
2
r cut i r L
skin
Constructing neighbor cells Linked cell method divide MD cell into smaller subcells :
n
n
n l
L n
r L
instead where
N N
(
N
1 )!
C
N
/
n
3 reducing it to
O
( N )
L
: the length of MD cell
r L
26 skin cells
Simulation Setup Simulation Cell Boundary Condition Constructing neighboring cells Initial atom velocities MD Time step Temperature Control
Initial Velocities Maxwell-Boltzmann distribution The probability of finding a particle with speed
P
(
v
x
)
m
2
k B T
1 / 2 exp 1 2
m
v
2
x
/
k B T
Generate random initial atom velocities scaling T with equipartition theorem 3 2
k B T
1 2
m
v
2
MD Time Step
r
/
t
1/20 of the nearest atom distance
t
MD is limited to <~100 ns
Temperature Control Velocity Scaling Scale velocities to the target T Efficient, but limited by energy transfer Larger system takes longer to equilibrate Nose-Hoover thermostat Fictitious degree of freedom is added Produces canonical ensemble (NVT) Unwanted kinetic effects from T oscillation
Integration Method Finite difference method Numerical approximation of the integral over time Verlet Method Better long-tem energy conservation Not for forces depending on the velocities Predictor-Corrector Long-term energy drift (error is linear in time) Good local energy conservation (minimal fluctuation)
Verlet Method
a
(
r
) 1
m
F
(
r
(
t
))
r
i
v
i
(
t
) Obtain the positions and velocities at
t
t
r
(
t
t
)
r
(
t
)
v
(
t
)
t
1 2
a
(
r
)
t
2
a
(
t
v
(
t
t
)
t
/ 2 ) 1
m
a
(
r
(
t
t
))
v
(
t
)
t
1 2
a
(
r
)
t
v
(
t
t
)
v
(
t
t
/ 2 ) 1 2
a
(
t
t
)
t
Predictor-Corrector Method Predictor Step
v
i
(
t
)
a
(
r
) 1
m
F
(
r
(
t
))
v
(
t
t
)
r
P
(
t
v
P
(
t
t
)
t
)
r
(
t
)
v
(
t
)
t
v
(
t
)
a
(
t
)
t
a
(
t
) 2
t
2
a
P
(
t
t
)
a
(
t
)
r
iii
(
t
)
t
r
iii
: 3rd order derivatives
Predictor-Corrector Method Corrector Step get corrected acceleration
a
C
(
r
)
F
(
r
P
(
t m
t
)) using error in acceleration
a
(
t
t
)
a
C
(
t
t
)
a
P
(
t
t
) correct positions and velocities
r
(
t
v
(
t
t
)
t
)
r v
P P
(
t
(
t
t
t
) )
C
0
C
1
t
2
a
(
t t
2
a
(
t
t
t
) )
C n
: constants depending accuracy
Force Calculation The force on an atom is determined by F
i
U
( r )
j N
i
u
(
r ij
r ij
) rˆ
ij U
( r ) : potential function
N r ij
: number of atoms in the system : vector distance between atoms
i
and
j
MD Potential Classical Potential
U
i U
1 (
r
i
)
i
,
j U
2 (
r
i
,
r
j
)
i j
, ,
k U
3 (
r
i
,
r
j
,
r
k
) ...
U
1 : Single particle potential Ex) external electric field, zero if no external force
U
2 : Pair potential only depending on
F
ij
r
i U
(
r ij
)
U
3 : Three-body potential with an angular dependence
F
i
i
j
(
V ij
V ji
)
j k V jki
Using Classical Potential Born-Oppenheimer Approximation
m e M
0 Assume total wavefunction as (
R
i
,
r
) (
R
i
) (
r
,
R
i
) (
R
i
) : Nuclei wavefunction (
r
,
R
i
) : Electron wavefunction parametrically depending on
R
i
The equation of motion for nuclei is given by
H N
i P i
2 2
M i
U
(
R
i
) (approximated to classical motion)
MD Potential Models Empirical Potential functional form for the potential fitting the parameters to experimental data Ex) Lennard-Jones, Morse, Born-Mayer Semi-empirical Potential calculate the electronic wavefunction for fixed atomic positions from QM Ex) EAM, Glue Model, Tersoff Ab-initio MD direct QM calculation of electronic structure Ex) Car-Parrinello using plane-wave psuedopotential
Stillinger-Weber Potential works fine with crystalline and liquid silicon
U
i
,
j U
2 (
r
i
,
r
j
)
i
,
j U
,
k
3 (
r
i
,
r
j
,
r
k
)
U
2 (
r ij
)
f
2 (
r ij
/ )
U
3 (
r
i
,
r
j
,
r
k
)
f
3 (
r
i
/ ,
r
j
/ ,
r
k
/ ) , : energy and length units Pair potential function
f
2
(
r
)
A
(
Br
p
r
q
) exp[(
r
a
)
1
] ,
0 ,
r
a r
a
Stillinger-Weber Potential Three body potential function
f
3 (
r
i
,
r
j
,
r
k
)
h
(
r ij
,
r ik
,
jik
)
h
(
r ji
,
r jk
,
ijk
)
h
(
r ki
,
r kj
,
ikj
)
h
(
r ij
,
r ik
,
jik
)
exp[
(
r ij
a
)
1
(
r ik
a
)
1
]
(cos
jik
1 3 ) 2
Stillinger-Weber Potential Limited by the cosine term (cos
jik
1 3 ) 2 forces the ideal tetrahedral angle not for various equilibrium angles too low coordination in liquid silicon incorrect surface structures incorrect energy and structure for small clusters Bond-order potential for Si, Ge, C bond strength dependence on local environment Tersoff, Brenner
Tersoff Potential cluster-functional potential environment dependence without absolute minimum at the tetrahedral angle The more neighbors, the weaker bondings
U
U repulsive
(
r ij
)
b ijk U attractive
(
r ij
)
b ijk
: environment-dependent parameter weakening the pair interaction when coordination number increases
Tersoff Potential
U ij
f C
(
r ij
)[
a ij f R
(
r ij
)
b ij f A
(
r ij
)] where repulsive part attractive part
f R
(
r
)
Ae
1
r f A
(
r
)
Be
2
r f C
(
r
) potential cutoff function 1 , 1 2 1 2
r
sin
r
R
2
R
(
r
0 ,
D
D R D
) ,
R
D
r
R
D
Tersoff Potential
b ij
( 1
n
ij n
) 1 / 2
n
ij
g
( )
k
i
,
j f C
(
r ik
)
g
(
jik
) exp[
3 3
(
r ij
1
c
2
d
2
d
2 (
h c
2 cos ) 2
r ik
)
3
]
a ij
( 1
n
ij n
) 1 / 2
n
ij
k
i
,
j f C
(
r ik
) exp[
3 3
(
r ij
r ik
)
3
]
Contents Introduction to MD Simulation Setup Integration Method Force Calculation and MD Potential MD Simulations of Silicon Recrystallization Simulation Preparation SW Potential Tersoff Potential
4 x 4 x 13 cells MD Simulation Setup Initial Setup 5 TC layer 1 static layer
MD Simulation Setup System Preparation Ion Implantation(1 keV) Cooling to 0K
Recrystallization 1200 K for 0.5 ns
Recrystallization SW Potential 1200K Crystal Rate a/c interface displacement
5 x 5 x 13 cells MD Simulation Setup Initial Setup 6 TC layer 1 static layer
MD Simulation Setup System Preparation Ion Implantation(1 keV) Cooled to 0K
Recrystallization 1900 K for 0.85 ns
Recrystallization Tersoff Potential 1900K Crystal Rate a/c interface displacement
Recrystallization Crystal Rate SW Potential 1200K Tersoff Potential 1900K
Recrystallization a/c interface displacement SW Potential 1200K Tersoff Potential 1900K
2 x 2 x 13 cells MD Simulation Setup Initial Setup 6 TC layer 1 static layer
Recrystallization 1800 K for 20 ns
Tersoff Potential Melting temperature of Tersoff: about 2547K Potential energy per particle versus temperature: the system with a/c interface is heated by adding energy at a rate of 1000K/ns
Tersoff Potential As in recrystallized Si : o 2 0.82 in amorphized Si o A 2
Tersoff Potential As in recrystallized Si : o A 2 in amorphized Si o 2 0.20 in crystalline Si
Summary Review Molecular Dynamics MD simulation for recrystallization of Si with SW, Tersoff with As