# Threading Dislocations Piercing the Free Surface of an Anisotropic Hexagonal Crystal: Review of Theoretical Approaches

Authors: Salem Neily; Sami Dhouibi; Roland Bonnet
Journal: Advances in Condensed Matter Physics
 (2018)

---

## Abstract

Inclined threading dislocations (TDs) piercing the oriented free surface of a crystal are currently observed after growth of oriented thin films on substrates. Up to date the unique way to treat their anisotropic elastic properties nearby the free surface region is to use the integral formalism, which assumes no dislocation core size and needs numerical double integrations. In a first stage of the work, a new and alternative approach to the integral formalism is developed using double Fourier series and the concept of a finite core size, which is often observed in high-resolution transmission electron microscopy. In a second stage, the integral formalism and the Fourier series approaches are applied to the important case of a TD piercing the basal free surface of a hexagonal crystal. For this particular geometry, easy-to-use expressions are derived and compared to a third approach previously known for a plate-like crystal. Finally, the numerical interest and the convergence of these approaches are tested using the basal free surface of the GaN compound, in particular for TDs with Burgers vectorsc and (a +c).


---

## Body

## 1. Introduction

Epitaxial growth of semiconductor materials offers opportunities in combining and modifying structural and electronic parameters [1–5]. Unfortunately, the thin crystalline film deposited on the substrate often contains high densities of dislocations piercing the free surface, the so-called threading dislocations (TDs), which affect adversely device properties. At high temperature, a reduction of the density of some TDs can occur via fusion or annihilation reactions. It is therefore of importance to study their intrinsic elastic properties. The aim of the present work is first to investigate the elastic field of an inclined TD piercing the free surface of any anisotropic, semi-infinite crystal. The inclination angle with respect to the normal is denoted θ. For isotropic crystals, exact solutions are well known for θ<π/2 [6, 7] and θ=π/2 [8, 9]. We present a new, alternative approach to the TD elastic field to that developed previously by Lothe et al. [10] who uses the integral formalism and requires in general a set of double integrations. This new approach is developed in terms of double Fourier series and can take account of a finite dislocation core size, which is not possible with the integral formalism. These two approaches are then applied to the important case of a hexagonal crystal (parameters a and c) limited by a basal plane. The particular geometry of a TD parallel to thec axis enables simplified expressions to be derived. These latter are compared to a third approach developed previously for a plate of finite thickness [11] from stress potentials [12]. Since wurtzite gallium nitride (GaN) has attractive optoelectronic applications, such as light emitting diode and laser diodes [1, 3–5], this compound is used in the present work to perform numerical applications.


## 2. Methods

### 2.1. The Fourier Calculation Approach

A description of the elastic field of inclined dislocations lying in an isotropic plate was recently presented from a Fourier calculation method [20, 21]. Below, this method is generalized to treat any anisotropic semi-infinite crystal. For the calculation, we use two Cartesian coordinate systems with the same origin O, Figure 1(a). The system Ox1x2x3 is linked to the free surface and oriented with the axis Ox3 directed outside matter. The direction Ox1 is any direction taken in the surface. The second system denoted Ox1′x2′x3′, is attached to a TD placed and oriented along Ox3′. This latter system is obtained from Ox1x2x3 by a rotation θ around Ox1. The TD terminates in O, is oriented by the unit vector ξ//Ox3′, and has the Burgers vectorb.Figure 1
Cartesian frames and symbols attached to (a) an isolated inclined threading dislocation (TD) oriented byξ//Ox3′ in the plane Ox2x3 and (b), two infinite sets of periodically distributed TDs oriented along ξ. Each of these TDs has a cylindrical core, the section of which is a small square.
(a) (b)The first step of the calculation is to consider a biperiodic series of largely spaced and identical emerging TDs. Its stress field can be written as the sum(j,k=1–3)(1)σjk=σjk∞+σjkrel,in which σjk∞ is the stress field generated in an infinite crystal and σjkrel is the unknown relaxation stress field that cancels the forces applied along the plane chosen for the surface. To fulfil (1) the idea is to use partial results obtained in [20] for σjk∞ and in [22] for σjkrel.Let us first consider the calculation ofσjk∞.In the planex3′=0, Figure 1(b), the theoretical dislocation lines form a square pattern of black points. These lines are distant by the period T along Ox1′ and Ox2′. Each dislocation line is oriented parallel to Ox3′ and is associated with a finite core. According to some high-resolution transmission electron microscopy observations and numerical simulations, for example, [23–25], a TD core can have a variety of atomic structures. These latter have certainly an influence on the elastic properties very near the theoretical dislocation line [25]. However, beyond a few b values and thanks to the St. Venant principle, elasticity can be applied without accounting for true limiting boundary conditions around the core. In the present work, for simplicity, the core section is chosen as a square, limited by an infinite prism with a section e2, limited by portions of planes parallel to planes x1′=x2′=±e/2. Inside the core the dislocation density is taken as b/e2, while outside the prism, the dislocation density is zero. Each prism can therefore be described as an extended core dislocation. Around the axis Ox3′, when T increases and the core size e decreases, the elastic field tends gradually to that of a classical Volterra dislocation. Along the axes Ox1 and Ox2, the periodicity of the defects is T and T/cos⁡θ.According to [20], the displacement uk∞ and stress fields σjk∞ of the infinite dislocation set can be written in the Ox1x2x3 frame in the form of double Fourier series. The σjk∞ stress field can be expressed as (integers m,n≠0)(2)σjk∞=Σm,nFjkm,nexp⁡2iπG1x1+G2x2+G3x3,in which the preexponential terms Fjkm,n can be calculated explicitly with(3)G1=mT,G2=nTcos⁡θ,G3=-nTsin⁡θ.The terms Fjkm,n depend on the elastic constants Cijkl taken in Ox1x2x3, on the Burgers vectorb, on the inclination angle θ, and on the geometrical parameters e and T. However, since their expressions are much too long to be written extensively, they are not given in the present text.Let us now consider the calculation ofσjkrel.Since this stress field should also be periodic, it can be expressed formally as a double Fourier series. As shown in [22], the formal forms of the associated displacement field can be written as (m and n≠0)(4)ukrel=Σm,nΣα=1–3Pαλαkexp⁡2iπpαx3exp⁡2iπG1x1+G2x2,in which the constants Pα are three unknown complex constants attached to each harmonic term. They depend on the order (m,n) of each harmonic term and on the limiting boundary conditions of the problem. The constants pα and λαk are found from the nontrivial solutions of the equations(5)ψjk+pαϕjk+pα2Cj3k3λαk=0.For α=1–3, the constants pα are complex chosen with a negative complex part. The corresponding constants λαk are normalized solutions of the same equation. In addition(6)ψjk=Cj1k1G12+Cj2k2G22+Cj1k2+Cj2k1G1G2,ϕjk=Cj1k3+Cj3k1G1+Cj3k2+Cj2k3G2.For convergence towards x3<0, only the negative complex parts of the constants pα are retained. The stress field is finally calculated from derivation of (4) and application of the Hooke law:(7)σjkrel=2iπΣm,nΣα=1–3PαLαjkexp⁡2iπpαx3exp⁡2iπG1x1+G2x2,The constants Lαjk are(8)Lαjk=λαqG1Cjkq1+G2Cjkq1+pαCjkq3.The unknown constants Pα are found from the following condition: the addition of the stress fields (2) and (7) are such that zero stresses act on the plane x3=0, that is, (k=1–3),(9)σ3k∞+σ3krel=0forx3=0.For each k value, this expression represents in fact a double Fourier series, which must be identical to zero. The only solution is that all its terms (m,n) are zero. This leads, for each term (m,n), to a system of three linear equations, the solution of which gives the values of the three complex unknowns Pα. From the knowledge of the Pα values associated with each term (m,n), the total displacement field is obtained by the addition (u∞+urel), in which u∞ is calculated from expression (5) in [20] and urel from expression (4).The total stress field is finally obtained from the addition of expressions (2) and (7).As a consequence, three parameters are required to evaluate the elastic field around an isolated TD: the periodT should much larger than b, while the size e could be estimated to a few b. Of course, a sufficiently large number of terms (m,n) must be retained to obtain a sufficiently precise solution.For the case of a hexagonal crystal,θ=0 and Ox3//c, simplified expressions can be obtained, briefly listed below. The axis Ox1 can be chosen anywhere in the basal plane.Constantspα are the following:(10)p1=-iC11-C122C44G12+G221/2,p2=-iβ1-β22C33C44G12+G221/2,p3=-i-β2+β12C33C44G12+G221/2.Compact expressions are then obtained for the Pα constants because the terms Lαj1 and Lαj2 are not involved in the calculations. As a result, for brevity, denoting Lαj3(α=1–3,j=1–3) as Lαj and Fj3 as Fj, respectively, the Pα constants become(11)P1=iβ3F3L22L31-L21L32-F2L23L31-L21L33-F1L22L33-L23L32,P2=iβ3F3L11L32-L12L31+F2L13L31-L11L33+F1L12L33-L13L32,P3=iβ3F3L12L21-L11L22+F2L13L21-L11L23+F1L12L23-L13L22with(12)β1=C132-C11C33C13+2C442-C11C331/2,β2=C13C13+2C44-C11C33,β3=2πL13L22L31-L21L32-L12L23L31-L21L33+L11L23L32-L22L32.Similar expressions are found for the particular cases m=0 or n=0.The elastic field described by the above Fourier approach can be compared with that calculated with the integral formalism [10]. After a brief presentation of this latter approach, we present simplified expressions for a hexagonal crystal (parameters a and c) that has a basal free surface.


### 2.2. The Integral Formalism Approach

In the integral formalism, the relaxation stress field of a TD dislocation is calculated for convenience in the frameOXYZ chosen in [10] Figure 1(a). The components of the Burgers vector are denoted b=bX,bY,bZ. It is generated by a plane fan-shaped distribution of straight infinitesimal dislocations placed in the plane x3=0 of the unlimited crystal. This distribution is located in the plane Z=0 and has the Burgers vector density bψ={b1(ψ),b2(ψ),b3(ψ)}. As in [10], ψ is the oriented polar angle of an infinitesimal dislocation linked to the plane OXY. The field σjkrel is obtained from double integrations: the first one is needed to compute the energy matrix B(ψ) related to each dislocation, while the second one is required to sum the contributions of all the planar set of dislocations.For the particular orientationθ=0, we obtain compact expressions for the image force dF that applies on a TD segment dλ in the plane Y=0 and for the relaxation stress σjkrel. The integral (24) given in [10] can be greatly simplified. If λ is a distance along the TD taken from the emerging point and dλ a length element, we obtain the simple expression(13)dF=bXbZC13C11-C12dλ4πC11λ.This expression shows clearly that dF≠0 for mixed dislocations, while dF=0 for pure screw bX=0 or prismatic edge dislocations bZ=0. This explains the trend, for mixed TDs, to have directions systematically misoriented with respect to thec axis, for example, [2].Whenθ=0 the density b(ψ) becomes simple to express because the calculation shows that the energy coefficient matrix B attached to an infinitesimal dislocation placed in the plane Z=0 is no more depending on its orientation ψ. This constant matrix can be evaluated from expressions (73,74) in [26] or other similar procedures indicated in [13, 27]. In passing, some misprints in [26] are indicated in Appendix A.The Burgers vector of an infinitesimal dislocation has the following simple form:(14)bψ=-1πbZC44B11cos⁡ψbZC44B22sin⁡ψbXcos⁡ψ+bYsin⁡ψC13C11-C12C11B33,in which B11, B22 and B33 are the nonzero elements of the energy coefficient matrix B.The relaxed displacement fieldurel in a point defined by its cylindrical coordinates (R1,ψ1,θ1) is then calculated from the addition of all the contributions of the infinitesimal dislocations with Burgers vector b(ψ)dψ, in the plane Z=0.As an example, for the edge case withb(bX,0,0), the nonzero components of the stress tensor are(15)σ31rel=bXC12-C11C13C442π3C11B33∫Δcos⁡ψ2ρ⌢1ψsin⁡ψ-ψ1dψ,σ32rel=bXC12-C11C13C442π3C11B33∫Δsin⁡ψcos⁡ψρ⌢1ψsin⁡ψ-ψ1dψ,in which [10] (see Appendix B)(16)ρ⌢1ψ=R1cos⁡θ12+sin⁡θ12sin⁡ψ-ψ121/2.For the screw case with b0,0,bZ, the stress tensor is(17)σrel=∫Δ2Hψcos⁡ψ2dψ∫Δ2Hψsin⁡2ψdψ0∫Δ2Hψsin⁡2ψdψ∫Δ2Hψsin⁡ψ2dψ000∫Δ2C13C11+C12Hψdψin which(18)Hψ=bZC12-C11C12+C11C448π3C11B11B22·B11+B22+B22-B11cos⁡2ψρ⌢1ψsin⁡ψ-ψ1.The integration domain Δ is indicated in [10] with a crossed integration symbol. Therefore, it is not performed on π exactly. For numerical applications, to avoid zero denominators in (13)–(16), an angular domain Δ running from (ψ1+ε) to (ψ1+π-ε) should be used with ε tending to zero. Our experience is that convergence is ensured with ε equal to 0.0001, a value sufficiently small to obtain a good precision on the stresses nearby the TD core.


### 2.3. The Stress Potential Approach

The displacement field regarding a mixed dislocation withξ//c(θ=0) was already expressed for a plate-like hexagonal crystal in [14] but unfortunately this author did not consider a semi-infinite crystal. In the present work, the solution is found in searching for the limit of the displacement field urel around the upper emerging point when the plate thickness tends to infinity. The calculation is carried out in the OXYZ frame, which gives the following result:(i) Edge dislocationb//OX: the field urel is derived from two stress potentials ϕ1 and ϕ2 mentioned in [14] but calculated at the limit of an infinite thickness:(19)urel=∂∂Xϕ1+ϕ2,∂∂Yϕ1+ϕ2,∂∂Zk1ϕ1+k2ϕ2with(20)ϕ1=bXC13C11-C12ν22πC11C441+k1ν2-ν1Yr∫λminλmaxJ1λrλ2exp⁡-ν1Zλdλ,in which r is the polar radius, while the parameters ν1, ν2, ν3, k1, k2 are given in Appendix C. The potential ϕ2 is obtained from ϕ1 by interchanging ν1 and ν2 and interchanging k1 and k2. J1(λr) is the Bessel function of the first kind for the integer order 1. In expression (17), since the integration domain cannot be exactly π, λmin should is taken close to 0, while λmax should be large enough. A few numerical applications show that values (λmin,λmax)=(0.0001,100) give similar values to those obtained with the Fourier approach and the integral formalism.(ii) Screw dislocationb//OZ: only the tangential component is nonzero:(21)uψrel=bZ2πν3∫λminλmaxJ1λrλexp⁡-ν3Zλdλ.Expression (21) only depends on C44 and C66 via ν3.


## 2.1. The Fourier Calculation Approach

A description of the elastic field of inclined dislocations lying in an isotropic plate was recently presented from a Fourier calculation method [20, 21]. Below, this method is generalized to treat any anisotropic semi-infinite crystal. For the calculation, we use two Cartesian coordinate systems with the same origin O, Figure 1(a). The system Ox1x2x3 is linked to the free surface and oriented with the axis Ox3 directed outside matter. The direction Ox1 is any direction taken in the surface. The second system denoted Ox1′x2′x3′, is attached to a TD placed and oriented along Ox3′. This latter system is obtained from Ox1x2x3 by a rotation θ around Ox1. The TD terminates in O, is oriented by the unit vector ξ//Ox3′, and has the Burgers vectorb.Figure 1
Cartesian frames and symbols attached to (a) an isolated inclined threading dislocation (TD) oriented byξ//Ox3′ in the plane Ox2x3 and (b), two infinite sets of periodically distributed TDs oriented along ξ. Each of these TDs has a cylindrical core, the section of which is a small square.
(a) (b)The first step of the calculation is to consider a biperiodic series of largely spaced and identical emerging TDs. Its stress field can be written as the sum(j,k=1–3)(1)σjk=σjk∞+σjkrel,in which σjk∞ is the stress field generated in an infinite crystal and σjkrel is the unknown relaxation stress field that cancels the forces applied along the plane chosen for the surface. To fulfil (1) the idea is to use partial results obtained in [20] for σjk∞ and in [22] for σjkrel.Let us first consider the calculation ofσjk∞.In the planex3′=0, Figure 1(b), the theoretical dislocation lines form a square pattern of black points. These lines are distant by the period T along Ox1′ and Ox2′. Each dislocation line is oriented parallel to Ox3′ and is associated with a finite core. According to some high-resolution transmission electron microscopy observations and numerical simulations, for example, [23–25], a TD core can have a variety of atomic structures. These latter have certainly an influence on the elastic properties very near the theoretical dislocation line [25]. However, beyond a few b values and thanks to the St. Venant principle, elasticity can be applied without accounting for true limiting boundary conditions around the core. In the present work, for simplicity, the core section is chosen as a square, limited by an infinite prism with a section e2, limited by portions of planes parallel to planes x1′=x2′=±e/2. Inside the core the dislocation density is taken as b/e2, while outside the prism, the dislocation density is zero. Each prism can therefore be described as an extended core dislocation. Around the axis Ox3′, when T increases and the core size e decreases, the elastic field tends gradually to that of a classical Volterra dislocation. Along the axes Ox1 and Ox2, the periodicity of the defects is T and T/cos⁡θ.According to [20], the displacement uk∞ and stress fields σjk∞ of the infinite dislocation set can be written in the Ox1x2x3 frame in the form of double Fourier series. The σjk∞ stress field can be expressed as (integers m,n≠0)(2)σjk∞=Σm,nFjkm,nexp⁡2iπG1x1+G2x2+G3x3,in which the preexponential terms Fjkm,n can be calculated explicitly with(3)G1=mT,G2=nTcos⁡θ,G3=-nTsin⁡θ.The terms Fjkm,n depend on the elastic constants Cijkl taken in Ox1x2x3, on the Burgers vectorb, on the inclination angle θ, and on the geometrical parameters e and T. However, since their expressions are much too long to be written extensively, they are not given in the present text.Let us now consider the calculation ofσjkrel.Since this stress field should also be periodic, it can be expressed formally as a double Fourier series. As shown in [22], the formal forms of the associated displacement field can be written as (m and n≠0)(4)ukrel=Σm,nΣα=1–3Pαλαkexp⁡2iπpαx3exp⁡2iπG1x1+G2x2,in which the constants Pα are three unknown complex constants attached to each harmonic term. They depend on the order (m,n) of each harmonic term and on the limiting boundary conditions of the problem. The constants pα and λαk are found from the nontrivial solutions of the equations(5)ψjk+pαϕjk+pα2Cj3k3λαk=0.For α=1–3, the constants pα are complex chosen with a negative complex part. The corresponding constants λαk are normalized solutions of the same equation. In addition(6)ψjk=Cj1k1G12+Cj2k2G22+Cj1k2+Cj2k1G1G2,ϕjk=Cj1k3+Cj3k1G1+Cj3k2+Cj2k3G2.For convergence towards x3<0, only the negative complex parts of the constants pα are retained. The stress field is finally calculated from derivation of (4) and application of the Hooke law:(7)σjkrel=2iπΣm,nΣα=1–3PαLαjkexp⁡2iπpαx3exp⁡2iπG1x1+G2x2,The constants Lαjk are(8)Lαjk=λαqG1Cjkq1+G2Cjkq1+pαCjkq3.The unknown constants Pα are found from the following condition: the addition of the stress fields (2) and (7) are such that zero stresses act on the plane x3=0, that is, (k=1–3),(9)σ3k∞+σ3krel=0forx3=0.For each k value, this expression represents in fact a double Fourier series, which must be identical to zero. The only solution is that all its terms (m,n) are zero. This leads, for each term (m,n), to a system of three linear equations, the solution of which gives the values of the three complex unknowns Pα. From the knowledge of the Pα values associated with each term (m,n), the total displacement field is obtained by the addition (u∞+urel), in which u∞ is calculated from expression (5) in [20] and urel from expression (4).The total stress field is finally obtained from the addition of expressions (2) and (7).As a consequence, three parameters are required to evaluate the elastic field around an isolated TD: the periodT should much larger than b, while the size e could be estimated to a few b. Of course, a sufficiently large number of terms (m,n) must be retained to obtain a sufficiently precise solution.For the case of a hexagonal crystal,θ=0 and Ox3//c, simplified expressions can be obtained, briefly listed below. The axis Ox1 can be chosen anywhere in the basal plane.Constantspα are the following:(10)p1=-iC11-C122C44G12+G221/2,p2=-iβ1-β22C33C44G12+G221/2,p3=-i-β2+β12C33C44G12+G221/2.Compact expressions are then obtained for the Pα constants because the terms Lαj1 and Lαj2 are not involved in the calculations. As a result, for brevity, denoting Lαj3(α=1–3,j=1–3) as Lαj and Fj3 as Fj, respectively, the Pα constants become(11)P1=iβ3F3L22L31-L21L32-F2L23L31-L21L33-F1L22L33-L23L32,P2=iβ3F3L11L32-L12L31+F2L13L31-L11L33+F1L12L33-L13L32,P3=iβ3F3L12L21-L11L22+F2L13L21-L11L23+F1L12L23-L13L22with(12)β1=C132-C11C33C13+2C442-C11C331/2,β2=C13C13+2C44-C11C33,β3=2πL13L22L31-L21L32-L12L23L31-L21L33+L11L23L32-L22L32.Similar expressions are found for the particular cases m=0 or n=0.The elastic field described by the above Fourier approach can be compared with that calculated with the integral formalism [10]. After a brief presentation of this latter approach, we present simplified expressions for a hexagonal crystal (parameters a and c) that has a basal free surface.


## 2.2. The Integral Formalism Approach

In the integral formalism, the relaxation stress field of a TD dislocation is calculated for convenience in the frameOXYZ chosen in [10] Figure 1(a). The components of the Burgers vector are denoted b=bX,bY,bZ. It is generated by a plane fan-shaped distribution of straight infinitesimal dislocations placed in the plane x3=0 of the unlimited crystal. This distribution is located in the plane Z=0 and has the Burgers vector density bψ={b1(ψ),b2(ψ),b3(ψ)}. As in [10], ψ is the oriented polar angle of an infinitesimal dislocation linked to the plane OXY. The field σjkrel is obtained from double integrations: the first one is needed to compute the energy matrix B(ψ) related to each dislocation, while the second one is required to sum the contributions of all the planar set of dislocations.For the particular orientationθ=0, we obtain compact expressions for the image force dF that applies on a TD segment dλ in the plane Y=0 and for the relaxation stress σjkrel. The integral (24) given in [10] can be greatly simplified. If λ is a distance along the TD taken from the emerging point and dλ a length element, we obtain the simple expression(13)dF=bXbZC13C11-C12dλ4πC11λ.This expression shows clearly that dF≠0 for mixed dislocations, while dF=0 for pure screw bX=0 or prismatic edge dislocations bZ=0. This explains the trend, for mixed TDs, to have directions systematically misoriented with respect to thec axis, for example, [2].Whenθ=0 the density b(ψ) becomes simple to express because the calculation shows that the energy coefficient matrix B attached to an infinitesimal dislocation placed in the plane Z=0 is no more depending on its orientation ψ. This constant matrix can be evaluated from expressions (73,74) in [26] or other similar procedures indicated in [13, 27]. In passing, some misprints in [26] are indicated in Appendix A.The Burgers vector of an infinitesimal dislocation has the following simple form:(14)bψ=-1πbZC44B11cos⁡ψbZC44B22sin⁡ψbXcos⁡ψ+bYsin⁡ψC13C11-C12C11B33,in which B11, B22 and B33 are the nonzero elements of the energy coefficient matrix B.The relaxed displacement fieldurel in a point defined by its cylindrical coordinates (R1,ψ1,θ1) is then calculated from the addition of all the contributions of the infinitesimal dislocations with Burgers vector b(ψ)dψ, in the plane Z=0.As an example, for the edge case withb(bX,0,0), the nonzero components of the stress tensor are(15)σ31rel=bXC12-C11C13C442π3C11B33∫Δcos⁡ψ2ρ⌢1ψsin⁡ψ-ψ1dψ,σ32rel=bXC12-C11C13C442π3C11B33∫Δsin⁡ψcos⁡ψρ⌢1ψsin⁡ψ-ψ1dψ,in which [10] (see Appendix B)(16)ρ⌢1ψ=R1cos⁡θ12+sin⁡θ12sin⁡ψ-ψ121/2.For the screw case with b0,0,bZ, the stress tensor is(17)σrel=∫Δ2Hψcos⁡ψ2dψ∫Δ2Hψsin⁡2ψdψ0∫Δ2Hψsin⁡2ψdψ∫Δ2Hψsin⁡ψ2dψ000∫Δ2C13C11+C12Hψdψin which(18)Hψ=bZC12-C11C12+C11C448π3C11B11B22·B11+B22+B22-B11cos⁡2ψρ⌢1ψsin⁡ψ-ψ1.The integration domain Δ is indicated in [10] with a crossed integration symbol. Therefore, it is not performed on π exactly. For numerical applications, to avoid zero denominators in (13)–(16), an angular domain Δ running from (ψ1+ε) to (ψ1+π-ε) should be used with ε tending to zero. Our experience is that convergence is ensured with ε equal to 0.0001, a value sufficiently small to obtain a good precision on the stresses nearby the TD core.


## 2.3. The Stress Potential Approach

The displacement field regarding a mixed dislocation withξ//c(θ=0) was already expressed for a plate-like hexagonal crystal in [14] but unfortunately this author did not consider a semi-infinite crystal. In the present work, the solution is found in searching for the limit of the displacement field urel around the upper emerging point when the plate thickness tends to infinity. The calculation is carried out in the OXYZ frame, which gives the following result:(i) Edge dislocationb//OX: the field urel is derived from two stress potentials ϕ1 and ϕ2 mentioned in [14] but calculated at the limit of an infinite thickness:(19)urel=∂∂Xϕ1+ϕ2,∂∂Yϕ1+ϕ2,∂∂Zk1ϕ1+k2ϕ2with(20)ϕ1=bXC13C11-C12ν22πC11C441+k1ν2-ν1Yr∫λminλmaxJ1λrλ2exp⁡-ν1Zλdλ,in which r is the polar radius, while the parameters ν1, ν2, ν3, k1, k2 are given in Appendix C. The potential ϕ2 is obtained from ϕ1 by interchanging ν1 and ν2 and interchanging k1 and k2. J1(λr) is the Bessel function of the first kind for the integer order 1. In expression (17), since the integration domain cannot be exactly π, λmin should is taken close to 0, while λmax should be large enough. A few numerical applications show that values (λmin,λmax)=(0.0001,100) give similar values to those obtained with the Fourier approach and the integral formalism.(ii) Screw dislocationb//OZ: only the tangential component is nonzero:(21)uψrel=bZ2πν3∫λminλmaxJ1λrλexp⁡-ν3Zλdλ.Expression (21) only depends on C44 and C66 via ν3.


## 3. Application to GaN

Gallium nitride is chosen to perform numerical applications because of its technological interest for electronic devices that can work under high power density and at relatively high temperature [4, 5]. The growth of GaN can be stabilized by epitaxial growth on bulky GaN substrates or a foreign substrate, for example, (111)Si, SiC, or sapphire [2, 4, 5, 15]. An investigation on GaN thin films is of particular interest because of its TDs that have anisotropic electrical properties regarding cathodoluminescence and photoluminescence. These properties are sensitive to the Burgers vectors and also to the TD directions [5].The thickness of the GaN layer is assumed to be large enough to neglect the elastic interaction of the TD with the GaN/substrate interface. Different anisotropic elastic constantsCij of GaN are reported in several experimental or theoretical works. These constants are taken from [16–19, 28–30] and are indicated in Table 1. GaN has the structure of wurtzite (space group P63mc) with lattice parameters a=0.319 nm and c=0.519 nm [15].Table 1
Elastic constants of GaN. In anisotropic elasticity, the equilibrium angleθeq of the inclined (a + c) TD is small, weakly depending on the existence of a free surface (columns 8 and 9). It is not the case in isotropic elasticity (columns 10 and 11). Symbols A, S and I means Anisotropy, Surface and Isotropy.
(GPa) C 11 C 12 C 13 C 33 C 44 μ (GPa), ν θ e q θ e q θ e q θ e q Refs.       AS A∞ IS I∞ [13] 367 135 103 405 95 113, 0.27 2.6° 3.8° 4.8° 11.6° [14] 388 154 84 458 85 112, 0.28 2.5° 3.6° 4.9° 11.6° [15] 329 109 80 357 91 107, 0.25 2.4° 3.6° 4.5° 11.6° [16] 390 145 106 398 105 120, 0.26 2.8° 4.4° 4.8° 11.6° [17] 374 106 70 379 101 121, 0.22 2.0° 3.0° 4.2° 11.6° [18] 359 129 92 389 98 113, 0.26 2.4° 3.7° 4.7° 11.6° [19] 334 132 99 372 86 100, 0.28 2.6° 3.8° 4.9° 11.6°Since the isotropic approximation is sometimes used in literature, column 7 indicates the calculated Young modulusμ and the Poisson ratio ν relative to each Cij set. μ and ν are found from an averaging of the elastic properties of GaN according to the method proposed in [31]. In this reference, a way is proposed to define an anisotropy factor coefficient A for any crystal symmetry. In the present work, using the Cij constants of [16–19, 28–30], the A values are found between 1.27 and 1.78. These two values correspond to the Cij constants mentioned in [19] and [17], respectively. Previous transmission electron microscopy investigations, for example, [2, 3, 5, 15], show that a large fraction of TDs cross the free surface with directions ξ parallel or not far from thec axis. They can be of edge-type, screw-type, or mixed-type. The local elastic equilibrium nearby the emerging point of a (a +c)-type TD can be treated from the image force theorem due to Lothe et al. [10]. For instance, taking OX//-a and a dislocation direction ξ in the prismatic plane (a,c), a systematic equilibrium angle θeq can be obtained. In Table 1, this angle is mentioned for each set of Cij constants given in [16–19, 28–30] including isotropic approximation. The comparison of columns 8 and 9 shows a small anisotropic effect since θeq changes from averages equal to 2.5° (free surface, column 8) to 3.7° (no free surface, column 9). On the contrary, with the isotropic approximation, columns 10 and 11 indicate a strong effect of the free surface since θeq increases from averages equal to 4.7° (free surface) to 11.6° (no free surface). These last results demonstrate that the isotropic approximation is not reliable for GaN.Other TDs in GaN are sometimes directed along thec axis, while other TDs have lines inclined noticeably with respect toc. For the screw case (θ=0) andb =c, the dislocation density b(ψ) is given by the expression (14). Taking the Cij constants indicated in [16], Figure 2 illustrates the sinusoidal shape of the two nonzero components of the dislocation density b(ψ), that is, b1(ψ) and b2(ψ). Our experience is that, for θ=0, the stress potential approach is more efficient to calculate the stress relaxation field than the integral formalism or the Fourier series.Figure 2
Case of thec screw TD oriented along the normal OZ. Evolution of the dislocation density b(b1,b2,b3) in the basal plane Z=0 of GaN.Forθ≠0, the Fourier approach and the integral formalism have similar performances. The numerical applications below, which use the elastic constants mentioned in line 2 Table 1, correspond to the geometry depicted in Figure 3. The TD is an inclined (a +c) screw dislocation placed in the Y=0 plane. Figure 4(a) describes the changes of the dislocation density in the plane Z=0. The component b3 is nonzero for this screw dislocation, Figure 4, conversely to the orientation θ=0.Figure 3
Orientationξ of the (a +c) screw TD. The inclination angle θ is 31.5°.Figure 4
Case of the (a +c) screw TD. (a) Change with ψ of the three components of the dislocation density b(b1,b2,b3) in the plane Z=0 and (b) evolution of the σ23 stress component along the straight line (X=0, Y=2 nm) according to different assumptions.
(a) (b)Figure4(b) exhibits, as an example, the evolution of the σ23 stress component along a straight line parallel to OZ and cutting the surface Z=0 at the point (X,Y) = (0, 2 nm). This inclined TD line is chosen because it gives some ideas about the results obtained from different approaches. All obtained curves converge far from the free surface but give different values close to the emerging point of the TD. The lower curve denoted “Aniso-∞” seems a bold line but is in fact the quasi-superposition of two curves, one calculated from the integral formalism and the other calculated from the Fourier series approach. For this latter approach, the calculation parameters were the following: period T=30 nm, core size e=b. The squared geometry of the core has therefore no practical impact on the values of the stress at points placed about 2 or 3 nm from the dislocation axis. The number of harmonic terms (m,n) chosen for summing the Fourier series was tested from visual convergence with the curve obtained from the integral approach. For points at 2 or 3 nm from the core, the number of harmonic terms should include a large number of terms for a good convergence (m and n can reach about a few hundreds). The same kind of difficulty is encountered with the integral formalism when it is applied to points close to the core. In the calculations, we adopted the integration domain π-ε with ε=0.0001 rad. Table 1 shows that, with the other sets of Cij constants, the amplitudes of these curves change hardly.The dashed curve denoted “Iso semi-∞” is calculated from the Yoffe-Shaibani-Hazzledine solution [6, 7] and the upper curve denoted “Aniso-∞” assumes an infinite anisotropic crystal. The part of this latter curve, towards the negative Z, is not represented but is symmetrical with respect to x3=0. The bold line denoted “Aniso-∞” shows that both approaches presented in Sections 2 and 3 give the same result for the σ23 stress and includes in particular the slight negative values nearby x3=0, that is, the region close to the dislocation core (a few nanometers). In this region, the stress is slightly negative, which is not the case in the isotropic approximation. This negative value nearby x3=0 is clearly due to the proximity to the TD core and not a calculation artefact generated by a lack of precision. Anisotropy effect is therefore far to be negligible nearby the free surface and the dislocation core.


## 4. Conclusions

In this work, the elastic field of an inclined TD has been investigated using different approaches. The first one, based on Fourier series and the concept of extended dislocation core, is developed in the present work from an appropriate combination of partial results obtained in [20, 22]. It proves to be a valuable alternative method to a well-known approach, the so-called integral formalism [10]. The computation times required by both approaches are found roughly similar and depend on the wanted precision. There is a strong interest in using these two independent approaches to test the ability of elasticity to give the same values around the theoretical core. Both approaches converge beyond one nanometer from the theoretical core, as verified for some TDs in GaN. A third approach using stress potentials and assuming θ=0 is also numerically tested for the first time (at the author’s knowledge).When the anisotropic constants tend to those of an isotropic crystal, we verified numerically that the three above approaches converge to the same stress field, namely, that obtained from the Yoffe-Shaibani-Hazzledine expressions [6, 7]. Correction of misprints should be introduced in [7]; see [32].For a dislocation normal to the basal free surface of a hexagonal crystal(θ=0), some general expressions of the integral formalism are simplified. They were tested positively in regions far from the free surface and close to the emerging core region of a TD. The numerical results exhibited in the four last columns of Table 1 and the curves in Figure 4(b) show clearly that the neglect of both elastic anisotropy and the presence of a free surface can be a confusion source in the interpretation of the angular equilibrium of a TD.Of course, it is expected that the atomic structure of the dislocation core could influence the elastic field at distance from the core lower than a few Burgers vector moduli. According to our experience the assumption of a square section of the dislocation core has no sensible influence on the stress field beyond a few nanometers from the dislocation line provided that the core sizee is less than the modulus of the Burgers vector.



---

*Source: 3038795-2018-05-10.xml*