Hello,
I am trying to reproduce epitaxial growth of germanium on silicon as discussed in PRB vol. 34 pp 5621 (1986) by Chris Van de Walle.
the paper talks about superlattice structures of Si and Ge layers. In the limit when a thin layer of Ge is deposited on a thick layer of Silicon the “in-plane” lattice constant of Ge would be equal to lattice constant of Si (5.43A). this corresponds to a ~4% strain compared to its bulk lattice constant (a_ge = 5.65A)
As a consequence of this in-plane stress, the thin Ge layer distorts tetragonally, i.e. out of plane lattice constant become 5.82A. this is about 7.2% tensile strain referenced to silicon lattice spacing.
The paper gives analytic formulas to calculate the tetragonal distortion in terms of elastic constants of silicon and germanium.
I am trying to simulate this scenario with LAMMPS with molecular statics using Tersoff potential.
-
I calculate elastic constant using Adian Thompson’s ELASTIC input script.
-
I wrote a simple code which extracts strain tensor from LAMMPS xyz files using intial and final “minimized” geometry.
-
I set up the problem with a two layer system. bottom layer is Silicon; top layer is Germanium. I choose periodic boundary conditions in the x, y directions and “shrink-wrapped” in the z-direction. additionally, since I do not want to include a lot of silicon atoms, I use setforce to zero forces in x-y-z directions for all silicon atoms (idea is to simulate to a very thick silicon which experience negligible distortion from top Germanium layer).
LAMMPS output from MS simulations indicate that the Germanium atoms do show negligible in-plane strain (i.e. they move very little from the intialized positions in the x-y plane) and a tetragonally distorted out-plane lattice spacing. the strain calculated from MD simulations is slightly lower than 7.2% value I quoted above.
I notice that the strain tensor elements using my method depend on the choice of atoms I “pick” for calculating strain. If my code selects atoms at the edge of x-y boundaries then the strain calculation give absurd answers. once I leave out the 1-2 atomic layers at the x-y edges of the simulation box, I get sensible answers. I have repeated simulations for different thickness of Ge layer and results consistently show the wisdom of leaving out the atoms located at the x-y edge of Germanium layer.
one may hypothesize possible reasons for disagreement between LAMMPS and analytic formulation (as discussed in the paper). Tersoff potential could be a source of problem including other sources errors. elastic constant calculated using Tersoff potential do not agree very well with expt.
however, I repeated above simulations with one modification – I used “setforce” command to zero out forces in the x-y plane (in-plane) for the Germanium atoms. what surprised me was that I get better agreement with the analytic formulas.
so I am left wondering if there is an issue with using periodic boundary conditions in this heterogenous system. I asked a friend who has expertise in continuum stress simulations and he told me that they have something called “reflecting boundary conditions” which imply that the displacement vectors at the x-y edges (if he were to setup analogous continnuum simulation) would be zero.
periodic boundary conditions with setforce 0.0 0.0 NULL seem to be emulating the “reflecting boundary conditions” situation (and hence better agreement with continuum like formulation). another interesting bit is that with this setup, I do not have to disregard atoms at the edge of x-y boundaries for my strain calculator.
But, I am perplexed why periodic boundary conditions by themselves are not sufficient ?
I know this is a very long read for anyone but if you have inputs please let me know.
thank you.
regards