Transcript maws2 5436

A Level-Set Method for Modeling Epitaxial
Growth and Self-Organization of Quantum Dots
Christian Ratsch, UCLA, Department of Mathematics
Collaborators:
Outline:
•Xiaobin Niu
•The level-set method for epitaxial growth
•Raffaele Vardavas
•Results for irreversible and reversible
aggregation
•Russel Caflisch
•Spatially varying diffusion can be used for
self-organization of islands (quantum dots)
• Coupling of level-set formalism with an
elastic model.
Modeling thin film growth
Methods used
(Atomistic) KMC simulations:
• Completely stochastic method
• Rate parameters can be obtained from DFT
(a)
(a)
9750-00-444
(i)
(c)
(g)
(h)
(f)
(e)
(b)
(d)
Continuum equations (PDEs):
• essentially deterministic
• no microscopic details.
• parameters can be obtained from atomistic
model (but difficult)
New Method
Time scale
Length scale
Atomic
motion
~ 1013
seconds
Ångstroms
Islands/
devices
Seconds hours
Microns and
larger
Level set method:
• PDE - based, (almost) deterministic
• atomistic details can be included
• microscopic parameters can be
obtained from DFT
Idea of the level set appproach
Island dynamics
Atomistic picture
(i.e., kinetic Monte Carlo)
(a)
F
(a)
9750-00-444
v
(i)
(c)
(g)
(h)
(f)
(e)
(b)
(d)
D
•Describe motion of island boundaries by a level-set function
•Adatoms are described in a mean-field approach with a diffusion equation
The level set method: schematic
Level set function j
Surface morphology
j=0
j=0
t
j=0
j=1
j=0
• Level set function is continuous in plane, but has discrete height resolution
• Adatoms are treated in a mean field picture
The level set method: the basic formalism
• Governing Equation:
j
 vn | j | 0
t
j=0
• Diffusion equation for the adatom density (x,t):
• Boundary condition:

dN
 F  D2   2
t
dt

   eq ( Ddet , x)
• Velocity:
vn  D(n     n    )
• Nucleation Rate:
dN
 D  ( x, t ) 2
dt
Seeding position chosen stochastically
(weighted with local value of 2)
• Stochastic break-up of islands (depends on:
Ddet , size, local environment)
Numerical details
Level set function
• 3rd order essentially non-oscillatory (ENO) scheme for spatial part of
levelset function
• 3rd order Runge-Kutta for temporal part
Diffusion equation
• Implicit scheme to solve diffusion equation (Backward Euler)
• Use ghost-fluid method to make matrix symmetric
• Use PCG Solver (Preconditioned Conjugate Gradient)
S. Chen, M. Kang, B. Merriman, R.E. Caflisch, C. Ratsch, R. Fedkiw, M.F. Gyure, and S. Osher, JCP (2001)
A typical level set simulation
Fluctuations need to be included in nucleation of islands
Nucleation rate:
Validation:
dN
 D  ( x, t ) 2
dt
max
Probabilistic seeding
weight by local 2
Scaling of island densities

Nucleation Theory:
N ~ (D/F)-1/3
Scaled island size distribution
C. Ratsch et al., Phys. Rev. B 61, R10598 (2000)
Detachment of adatoms and breakup of islands
• Detachment of atoms (from boundary) is accounted for by boundary condition:
   eq ( Ddet , x)
• The numerical timestep remains unchanged. Thus, no increase in CPU time!
• Stochastic element is needed for breakup of islands
• For “small” islands, calculate probability of island
break-up.
• This probability is related to Ddet, and local
environment
• Pick random number to decide break-up
• If island is removed, atoms are distributed uniformly
in an area that corresponds to the diffusion length
Validation: Scaling and sharpening of island size distribution
Experimental Data for
Fe/Fe(001),
Stroscio and Pierce,
Phys. Rev. B 49 (1994)
Petersen, Ratsch, Caflisch, Zangwill, Phys. Rev. E 64, 061602 (2001).
Computational efficiency
• Fast events can be included without decreasing the numerical timestep (due to
mean-field treatment of adatoms)
Modeling self-organization of quantum dots
• Ultimate goal: Solve elastic equations at every timestep, and couple the strain field
to the simulation parameters (i.e., D, Ddet).
• This is possible because the simulation timestep can be kept rather large.
• Needed: Spatially varying, anisotropic diffusion and detachment rates.
Modifications to the code will be discussed!
• So far: We assume simple variation of potential energy surface.
• Next (with some preliminary results): couple with elastic code of Caflisch,
Connell, Luo, Lee
Vertical alignment of stacked quantum dots
Stacked InAs quantum dots on GaAs(001)
•Islands nucleate “on top” of lower
islands
•Size and separation becomes more
uniform
•Interpretation: buried islands lead
to strain (there is a 7% misfit)
Spatially varying potential
energy surface
Spatially varying nucleation
probability
B. Lita et al. (Goldman group), APL 74, 2824 (1999)
Aligned islands due to buried dislocation lines
Ge on relaxed SiGe buffer layers
• Islands align along lines
• Dislocation lines are buried underneath
• Interpretation: buried dislocation lines
lead to strain
Spatially varying potential
energy surface
Spatially varying nucleation
probability
H. J. Kim, Z. M. Zhao, Y. H. Xie, PRB 68,
205312 (2003).
Level Set formalism is ideally suited to
incorporate anisotropic, spatially varying
diffusion and thus nucleation without
extra computational cost
Modifications to the level set formalism for non-constant diffusion
0 
 Dxx (x)


• Replace diffusion constant by matrix: D  D(x) 
 0
Dyy (x) 

Diffusion in x-direction
• Diffusion equation:

dN
 F    D( )  2
 drift
t
dt
drift ~ Dxx  x Ead  D yy  y Ead
• Velocity:
vn  n  D( )   n  D(  )
• Nucleation Rate:
Dxx (x)  Dyy (x)
dN

 (x, t ) 2
dt
2
Diffusion in y-direction
Possible variations of
potential energy surface
no drift
drift
Isotropic diffusion with sinusoidal variation in x-direction
Dxx  D yy ~ sin( ax)
Only variation of transition
energy, and constant adsorption
energy
• Islands nucleate in regions
of fast diffusion
• Little subsequent nucleation
in regions of slow diffusion
fast diffusion
slow diffusion
Comparison with experimental results
Simulations
Results of Xie et al.
(UCLA, Materials Science Dept.)
Isotropic diffusion with sinusoidal variation in x- and y-direction
Dxx  D yy ~ sin( ax) sin( ay )
Anisotropic diffusion with variation of adsorption energy
What is the effect of thermodynamic drift ?
Etran
Ead
Spatially constant
adsorption and transition
energies, i.e., no drift
small amplitude
large amplitude
Regions of fast surface diffusion
Most nucleation does not occur in region of fast diffusion, but is dominated by drift
Transition from thermodynamically to kinetically controlled diffusion
Constant transition
energy
(thermodynamic drift)
Constant adsorption
energy
(no drift)
In all cases, diffusion constant D has the same form:
D
•No drift (right): nucleation dominated by fast diffusion
•Large Drift (left): nucleation dominated by drift
dN
 D  ( x, t ) 2
dt
x
Time evolution in the kinetic limit
• A properly modified PES (in the “kinetic limit”) leads to very regular, 1-D
structures
• Can this approach used to produce quantum wires?
Combination of island dynamics model with elastic code
•In contrast to an atomistic (KMC) simulation, the timestep is rather large, even
when we have a large detachment rate (high temperature).
•A typical timestep in our simulation is O(10-2 s); compare to typical atomistic
simulation, where it is O(10-6 s).
•This allows us to do an “expensive calculation” at every timestep.
•For example, we can solve the elastic equations at every timestep, and couple
the local value of the strain to the microscopic parameters.
•This work is currently in progress ….. but here are some initial results.
Our Elastic model
•Write down an atomistic energy density, that includes the following terms
(lattice statics) (this is work by Caflisch, Connell, Luo, Lee, et al.):
E  k ( S xx  S yy )
2
Nearest neighbor springs
2
Diagonal springs E  kdiag (S xx  2S xy  S yy ) 2  kdiag (S xx  2S xy  S yy ) 2
Bond bending terms
E  mSxy
2
•This can be related to (and interpreted as) continuum energy density
E   ( S xx  S yy )  S xy  S xx S yy
2
2
2
• Minimize energy with respect to all displacements: u E [u] = 0
Numerical Method
• PCG using Algebraic MultiGrid (poster by Young-Ju Lee)
• Artificial boundary conditions at top of substrate (poster by Young-Ju Lee)
• Additional physics, such as more realistic potential or geometry easily included
Couple elastic code to island dynamics model
Example:
Sxx
• Epilayer is 4% bigger than
substrate (I.e., Ge on Si)
• Choose elastic constants
representative for Ge, Si
• Deposit 0.2 monolayers
Syy
Modification of diffusion field
The dependence of D on strain can
be based on DFT results.
Example: Stain dependent
diffusion for Ag/Ag(111)
Ediff  Ediff , 0  const  ( S xx  S yy )
C. Ratsch, A.P. Seitsonen, and M. Scheffler
Phys. Rev. B 55, 6750-6753 (1997).
Results with strain-dependent detachment rate
Constant diffusion
Change diffusion as a function of
strain at every timestep
•It is not clear whether there is an effect on ordering
•More quantitative analysis needed
Modification of detachment rates
Edet  Edet,0  const  S xx2  S yy2
• The detachment rate has only
physical meaning at the island edge
(where it changes the boundary
condition eq)
• The model shown here indicates that
it is more likely to detach from a
bigger (more strained island) than from
a smaller one.
• Previous (KMC) work suggests that
this leads to more uniform island size
distribution.
Results with strain-dependent detachment rate
No change of Ddet
Strain induced change of Ddet at
every timestep
•Maybe fewer islands are close together upon strain induced increase of Ddet (?)
•Obviously, a more quantitative analysis is needed!
Conclusions
• We have developed a numerically stable and accurate level set method to
describe epitaxial growth.
• Only the relevant microscopic fluctuations are included.
• Fast events can be included without changing the timestep of the simulations.
• This framework is ideally suited to include anisotropic, spatially varying
diffusion.
• A properly modified potential energy surface can be exploited to obtain a high
regularity in the arrangement of islands.
• We have combined this model with a strain model, to modify the microscopic
parameters of the model according to the local value of the strain.
Essentially-Non-Oscillatory (ENO) Schemes
Need 4 points to discretize  j with third order accuracy
i-3
i-2
i-1
i
i+1
i+2
i+3
i+4
Set 1 Set 2 Set 3
This often leads to oscillations at the interface
Fix: pick the best four points out of a larger set of grid points to get rid of
oscillations (“essentially-non-oscillatory”)
Solution of Diffusion Equation
• Standard Discretization:
ik 1  ik
t
D

dN
 F  D 2   2
t
dt
ik11  2 ik 1  ik11
(x) 2
Aρ k 1  b
• Leads to a symmetric system of equations:
• Use preconditional conjugate gradient method
f  0
Problem at boundary:
  f   i    i   i 1 

  



x

x



(  xx ) i   1
1
1x  x 
2
i-2
Matrix not symmetric anymore ; replace by:
  g  i   i  i 1 

  


x

x



(  xx )i  
x
i
i-1
1x
g
 g : Ghost value at i
“ghost fluid method”
i+1