Modelling anisotropic lateral oxidation from circular mesas

In this paper, an iterative method to model the anisotropic lateral oxidation of circular structures is proposed and validated by confrontation to experimental data. The described model enables the efficient calculation of the temporal bi-dimensional evolution of the oxidation front shape, starting from a circular mesa, and progressing inward as a result of an anisotropic process combining an isotropic diffusion with an anisotropic reaction. The result of the developed model shows that the oxide aperture smoothly deforms from a circle to become more diamond-like, mimicking the experimental situation encountered when fabricating Vertical-Cavity Surface-Emitting Lasers (VCSELs) on (100) wafers or, more generally, when oxidizing circular mesas of aluminum-containing III-V semiconductor on similarly oriented substrates. © 2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement OCIS codes: (160.2100) Electro-optical materials (160.6000) Semiconductor materials (250.7260) Vertical cavity surface emitting lasers (310.1860)


Introduction
The wet oxidation of aluminum-containing III-V semiconductors is a process where water vapor induces a compositionally-selective transformation of high-refractive-index (n Al-III- V ~2.9) high-Al-content semiconductor alloys into low-refractive-index (n AlOx ~1.6) insulating aluminum oxides [1].This change in material properties is typically used to define the electrical injection path and the optical field distribution in a variety of devices [2] including Vertical-Cavity Surface-Emitting Lasers (VCSELs) [3] [4], optical waveguide devices [5] [6] [7] [8] [9], transistors [10] and transistor lasers [11].In the former two cases, the device internal and external (fiber-coupling) efficiencies critically depend on the size and shape of the oxide aperture i.e. depend on the geometrical characteristics of the interface between the (unoxidized) semiconductor and the oxidized areas resulting from the partial lateral oxidation of mesas whose etched sidewalls have been exposed to water.To address the above-mentioned requirement and potentially enhance the device manufacturability, several models have been developed to reproduce the process kinetics [12] [13] [14] [15] [16] [17] [18] [15] [19] [20] [21] [22].These models are, in essence, all based on the empirical model established by Deal and Grove which describes and quantifies the temporal evolution of the oxidation depth of surface oxidized silicon wafers in terms of the interplay between the diffusive and reactive nature of the process [23].As the Al-III-V semiconductor oxidation progresses laterally from the sidewalls of etched mesas rather than from the wafer surface, the reported models have primarily been concerned with the extension from one dimension to 2-to-3-dimensions as well as with the analysis of the associated consequences.It has thus been shown that the process kinetics depends on the etched mesa geometry [12] [13] [14] [15] [16] [17], on the thickness of the to-be-oxidized layer [18] [15] [19] [20] and of its (oxidizing) surroundings [21] [22].In all these models, the process is assumed to be isotropic and, therefore, the shape of the oxide aperture would then readily be inferred from the etched mesa geometry by a homothetic transformation.It has, however, been experimentally established that the Al-III-V semiconductor oxidation is anisotropic [24] [25], which, in turns, means that the shape oxide aperture can significantly differ from the scaled etched mesa contour [26] [27] [28].
In this paper, we introduce a perturbative refinement to include the process anisotropy into an existing bi-dimensional model capable of simulating the temporal characteristics of the oxidation of circular mesas [15].We subsequently validate this approach by reproducing the gradual deformation of a thin-layer oxide aperture starting from a circular mesa and progressively morphing into a diamond-like shape resulting from the oxidation of Al(Ga)As disk patterns etched into (100)-oriented wafers [26] [27] [28] as commonly observed when fabricating VCSELs (see Fig. 1).

Experimental data
The sample under study was grown onto a (100) GaAs wafer by molecular beam epitaxy and consisted of a 70-nm-thick AlAs-layer capped by a 50-nm GaAs layer.Circular mesas were subsequently lithographically defined and transferred by dry etching to expose the edges of the AlAs layer.The oxidation was carried out in a furnace with an in-situ infrared reflectometry imaging system allowing the time-dependent evolution of the oxidation front to be monitored with a spatial resolution of 0.8 µm.The sample was oxidized in conditions leading to anisotropic oxidation [28], more precisely, using a reduced-pressure (∼0.5 atm.) environment, a substrate temperature of 400 °C and a mixed H 2 /N 2 /H 2 O gas steam generated by an evaporator-mixer system operating at 95 °C.Chemically, the oxidation can be described using the following set of thermodynamically-favorable reactions [29]: t=51'13'' t=0'0'' t=10'57'' t=19'22'' t=39'25'' t=45'40'' Fig. 1.Left: Timed microscope images taken during the process showing of the evolution of the oxidation front for a 30µm-radius mesa and, right, extracted evolution along the <100> and <110> directions together with the fitted isotropic model simulation.
Fig. 1 (left) shows a set of images illustrating the gradual evolution of the oxidation front starting from a 30-µm-radius circle to a close-to-diamond shape as the oxidation progresses.Because of the anisotropy, the oxidation extent along the <100> and <110> directions (corresponding to the oxidation fast and slow axis respectively) differs and this difference increases with time until it reaches 1.65 µm after 33 min 11 s of oxidation (corresponding to an average oxidation length of 15.5 µm or equivalently to an aperture radius of 14.5 µm).Beyond that time, the reduction in the aperture size becomes dominant and the anisotropic difference in oxide extent between the above-mentioned axes diminishes.Contributing factors to the slower progression along the <110> directions (compared to the <100> directions) may include the higher densities of Al-atoms and surface-accessible bonds (including in-plane bonds) in the <110> directions mitigated by their smaller inter-atomic-plane distance.

General presentation
As stated above, the proposed model is a variant on an established model which can reproduce the oxidation dynamics of circular mesas [15].A generalized formalism of the model used in [15] is presented hereafter allowing both isotropic and anisotropic oxidations of a thin layer to be simulated.
The lateral oxidation of Al(Ga)As compounds is relatively complex process which combines several steps: the transfer of the reactant (H 2 O) from the gas supply into the (oxidized) peripheral area of the etched mesa, its diffusion from the edge of the mesa to the oxide/semiconductor interface where a set of cascaded chemical reactions occur with the Alcontaining semiconductor material to form the oxide [29].These stages are followed by the diffusion of the by-products (H 2 , AsH 3 ) of the chemical reactions through the oxidized area towards the edge of etched mesa where they transfer to the gas supply.To date, all the models rely on a simplified analysis of the above-described process where the system is supposed to be controlled by the diffusion or the reaction of a limiting species (the reactant or one of the by-products) which does not change with time/oxidation extent.Mathematically, the process can thus be analyzed by solving Fick's second law of diffusion in the oxidized area: where C( ⃗, ) is the considered species concentration at position  ⃗ and at the oxidation time, t, and  � �⃗ stands for the diffusion tensor.The v( ⃗, ) coefficient represents the convection term (including bulk convection and other diffusing-limiting phenomena such as strain-induced effects [20] or/and pore blocking [30]).The associated boundary conditions are set by the gas transfer at the edge of the mesa (external periphery,  ⃗  , with gas transport coefficient, h) and by the consumption/production of the element at the oxide/semiconductor interface ( ⃗  ) with an orientation-dependent distribution of reaction rates, { �⃗ }: )  �⃗  and  �⃗  being the unitary outward-pointing surface-normal vector at the external and internal interfaces respectively.Eq. ( 3) states that the oxidation front locally progresses with an extent which is proportional to the local concentration in the considered species and also proportional to the largest projection of the distribution of reaction rates in the direction orthogonal to the oxide/semiconductor surface.As such, it is worth noting that, even if the reaction rate were to be null in a specific direction, the oxidation of a continuous front would still proceed at a nonzero speed in that particular direction.The speed at which the semiconductor/interface progresses is then given by: where N is the number of reactant/by-product molecules needed to create a unit volume of oxide.
Assuming that the diffusion is isotropic (  � �⃗ becomes a constant, D, since the oxide is amorphous and considered spatially homogeneous) and that the convection coefficient is also spatially uniform and stationary (no delamination), equations ( 1) to ( 3) can be re-written into cylindrical coordinates as follows: ρ ����⃗ ,  θ ����⃗ and   ����⃗ being respectively the unit vectors in the radial, azimuthal and vertical directions of the cylindrical coordinate system.Following the approach described in [15], ( ⃗  ) can be calculated by solving Eq. ( 5) in the stationary regime with the concentration distribution, C, and the oxidation speed, σ, separated into their vertical and lateral profiles (i.e.C(ρ,θ,z)=C(ρ,θ)C(z)).As such, the oxidation simulation essentially becomes an in-plane (bi-dimensional) problem.The in-plane spatial dependence of σ is then related to the oxidation time through: (8) taking  ⃗  and  ⃗  along the shortest oxidation path.

Isotropic scenario
To begin with, we consider the simplest scenario possible: the isotropic case where, in addition to the diffusion and convection, the reaction rate is taken to be independent of the direction i.e. { �⃗ } is a sphere and the oxide progression is equal to kC (right term of Eq. ( 7)).As a result, the in-plane dependence of C (and σ) becomes one-dimensional (purely radial) and the inward-progressing oxidation of a circular mesa of etched radius, R mesa , will lead to circular oxide apertures.Assuming that the oxidation only occurs in one AlGaAs layer and that its thickness is sufficient to avoid impeding the diffusion (~50 to ~100 nm), then, the oxidation time, t, needed for the semiconductor/oxide interface to reach a radius, R, is given by [15] [14]: where A=2D(k+h)/(kh), B=2DC * /N and β=k/(k+h).
The spatial dependence of the concentration is given by [15]: 6

Anisotropic(-reaction) scenario
In the remainder of this article, we investigate the situation where the process is anisotropic.We still consider that the oxidation lateral extent is much larger than the thickness of the oxidizing layer and, therefore, that the simulation remains a bi-dimensional (ρ,θ) problem.In this context, the chemical reaction anisotropy is introduced by setting the in-plane distribution of reaction rates to be directionally-dependent such that: where ϕ is the in-plane angular direction (referenced to   ����⃗) and (  ����⃗,   ����⃗,   ����⃗) is a fixed unitary set of orthogonal axes oriented with respect to the semiconductor crystallographic axes.The inherent 2-π rotational symmetry of the reaction rate distribution means that the realvalued k(ϕ) function can be defined for ϕ∈[0;2π[ and be expressed as a Fourier series: (12) k s,n and k c,n being the Fourier decomposition coefficients.With a N-fold symmetry imposed by the semiconductor crystallography k(ϕ) becomes: The procedure to implement the model with an anisotropic reaction (and an isotropic diffusion and convection) is based on the iterative use of Eq. ( 5) to (7) with { �⃗ } defined in (11) and (12).It relies on understanding that the progression of the oxidation front with time is the result of a succession of an ensemble of elementary point-like oxidations of duration, dt, whose local extent around the point,  ⃗, in the direction, ϕ, is given by ( ⃗, ) = ()( ⃗)/.Assuming that the anisotropy leads to a perturbative modification of the oxidation front compared to the isotropic (circular aperture) case, the concentration at the oxide/semiconductor interface ( ⃗) is then taken to only depend on radius, and thus ( ⃗) = ().The oxidation front at t+dt can therefore be deduced by performing a (morphological) erosion of the oxidation front at time, t, with a fixed shape kernel defined by ()()/.In practice, the etched mesa contour (corresponding to the initial oxidation front) is first discretized regularly (with an angular step, Δθ) into N s samples (R(θ j ,t=0) with θ j =jΔθ ; j∈[1;N s ]).The iterative procedure involves performing the following three consecutive steps repeatedly: 1. Estimating of the (initial) mean oxidation speed, σ mean,i , using Eq. ( 10) using k mean , the mean of k(ϕ), as the k parameter and R mean,i , the mean of R(θ j ,t i ), for R. In doing so, the additional assumption made is that, irrespective of the oxidation front shape, the concentration at the oxide/semiconductor interface is uniform (i.e. the effective diffusion length is large compared to reaction-induced inhomogeneities or, by extension, compared to the size of the oxide aperture).2. Defining a time step, Δt i, longer than the chosen time resolution limit, Δt lim , and such that the point oxidations originating from the R(θ j ,t i ) will overlap (Δt i =3πR mean,i /σ mean,i ). 3. Calculating the (sampled) shape of the oxide aperture R(θ j ,t i+1 ) at t i+1 = t i +Δt i .
According to Eq. ( 7) and ( 8), this is achieved by defining R(θ j ,t i+1 ) as: R(θ j ,t i+1 )=R(θ j ,t i )-max� �⃗ •  �⃗  � /  •  , • ∆  .(15) Once again, eq. ( 15) can be physically interpreted as carrying out a set of fixedshape anisotropic oxidations starting from each of the individual R(θ j ,t i ) points followed by the selection of the point the furthest away from the previous oxidation front in each angular sector.Practically, k(ϕ) is thus discretized regularly (with an angular step, Δϕ) into N ϕ samples.Using k(ϕ m )= k(mΔϕ; m∈[1; N ϕ ]), a set of N ϕ point-oxidation boundary points, R point (ϕ m ,θ j ,t i ), is defined from each R(θ j ,t i ) as the result of the (Minkowski) vectorial summation of the (θ j ,R(θ j ,t i )) position vectors with the (ϕ m , k(ϕ m )/  •  , • ∆  ) position vectors.The resulting oxidation front, R(θ j ,t i+1 ), is then taken to be the R point (ϕ m ,θ j ,t i ) with the minimum radius in the θ j ±Δθ/2 angular sector.This 3 rd -step calculation process illustrated in Fig. 2. Fig. 2. Illustration of the calculation of the sampled "new" oxidation front (R(θj,ti+1), red points) by point-like anisotropic oxidation (black crosses) starting from the sampled "old" oxidation front (R(θj,ti), blue crosses).The dashed grey lines show the limits of the considered angular sectors.

Parametric study
The first part of the numerical study aims to get acquainted with the influence of the various parameters used.The above-described model was thus implemented using, as established in [14], N/C * =6 and β=0.95 i.e. k=19h.The oxidation fronts were sampled to provide a calculation spatial resolution (=2πR mesa /N s ) smaller than ~0.5 µm i.e. with N s =361 points for a 30-µm-radius mesa.Finally, the point-like oxidations were evaluated using N ϕ =5(N s -1)+1 points and, in light of the experimental data (see Fig. 1), k(ϕ) included a 4-fold symmetry oriented with ϕ=0 corresponding to the slow oxidation axis (<110> direction).To limit the number of simulation parameters, the Fourier decomposition of k(ϕ) was restricted to its first order and re-written as: with k mean being the mean reaction rate of the oxidation and b∈[0;2] representing the degree of anisotropy.
Unless specified differently, the parameters are taken to be D=78 µm 2 /min, k mean =55.1 µm/min (corresponding to h=2.9 µm/min) and b=1.98.Fig. 3 shows the timedependent evolution of the oxidation extent in the fast and slow directions when varying one parameter at a time.As it can be observed, changing the diffusion coefficient, D, mostly impacts the dynamics when the oxidation extent becomes large since the diffusion across the oxide becomes the limiting process factor.In contrast, the reaction rate, k mean , (and associated gas transport coefficient, h=19k mean ) controls the early stages of the oxidation since, then, the contributing species can be readily exchanged with the supply gas.Finally, as seen in Fig. 4, the effect of anisotropy coefficient can mostly be detected for aperture sizes around R mesa /2, as a result of the trade-off between the divergence induced by the anisotropy and convergence imposed by the closure of the aperture for large oxidation extents.As could be anticipated, the degree of anisotropy leaves mostly unaffected the general oxidation dynamics but increases the difference in oxidation extent between the slow and fast axis as it increases as explicitly shown on the inset of Fig. 4.

Experimental validation
Having established that the results of the proposed model are consistent with the lateral oxidation underlying physical mechanisms, the remaining task was to ensure that it could be successfully used to fit experimental data.
As suggested in the previous paragraphs, the procedure is performed in two stages, starting by fitting the time dependence of the averaged oxidation extent neglecting the anisotropy and, subsequently, adjusting the degree of anisotropy to obtain a match with the spatial degeneracy introduced between the oxidation slow and fast axes.Fig. 5 displays the result of this fitting procedure which demonstrates that the dataset can be reproduced with greater accuracy using the anisotropic model using as parameters D=78 µm 2 /min, k mean =55.1 µm/min and b=1.98.Indeed, the root-mean-square error between experimental and modelled extensions along the <100> and <110> directions is 0.62 µm for the anisotropic model and 0.84 µm for the isotropic case.In Fig. 6, the calculated oxidation profiles are superimposed to the corresponding microscopic image (of Fig. 1).This set of images not only shows the gradual deformation from a circle to a diamond-like shape but also demonstrates a visual agreement between the simulation and experimental data throughout the oxidation duration.Typical quantitative analysis of the error between the extracted experimental oxidation profile and the modelled ones for isotropic and anisotropic scenarios is shown in Fig. 7 in a case where the anisotropic deformation can be clearly identified.It demonstrates that the anisotropic model exhibits a significantly reduced error compared to the isotropic case, more specifically along the slow axes (0°-180° ; 90°-270°).The error extrema and sign can also be traced back to the difference between the simulated and experimental curves shown on Fig. 5. Furthermore, the slight asymmetry in the angular distribution of the error between the [0 1 0] (45°-225°) and the [0 -1 0] (135°-315°) directions suggests that the assumed 4-fold symmetry might not be entirely respected, a feature rather difficult to observe from the images presented in Fig. 1 and which may actually result from an slight irregularity in the mesa circularity rather than from the semiconductor crystalline structure (as it was not seen either in [25]).Finally, the comparison of the root-mean-square errors resulting from the use of the isotropic and anisotropic models is presented on Fig. 8.It shows that the anisotropic model is more precise as soon as the aperture becomes diamond-like (i.e., here, for oxidation times greater than 20 min) and can lead to a reduction of the absolute distance error by a factor of up to 1.5.
At this point, it should be emphasized that increasing the spatial resolution of the calculation does not affect the observed error as the radial progression is calculated by proxy and not discretized on a fixed grid.The fact that the error is maximal around the <100> directions (45°-225° and 135°-315°) suggests that the curvature of k(φ) is not steep enough around these directions.To further improve the match between the simulated contours and the experimental data, the only solution would be to use a higher-order Fourier decomposition of k(φ) and therefore increase the number of fitting parameters.

Conclusions
In this paper, a model permitting the quantitative bi-dimensional calculation of the timedependent evolution of the shape of the oxide apertures resulting from the anisotropic lateral oxidation of circular mesas has been established.The approach is based on an iterative procedure which implements the anisotropic progression as the superposition of (anisotropic) point-like oxidations whose mean extent is ruled by the isotropic dynamic behavior, given here by an analytical formula.The proposed model was proven to lead to results that are not only in quantitative agreement with experimental data but also allow a more accurate determination of the oxide aperture shape than the previously-established isotropic models.
The match between the model and the experiment may be further improved with refined (and thereby more complex) shapes for the point-like oxidations.Nevertheless, the presented version of the model is sufficiently accurate and simple to be used in three-dimensional VCSEL models to enable a more precise fitting and prediction of their opto-electronic characteristics or to be introduced in oxidation monitoring softwares to facilitate the control of the produced oxide apertures.Being a refinement over an isotropic model, it retains the ability to analyze the process dependence to temperature and to the composition and thickness of the to-be-oxidized alloy with the added benefit of being able to track the anisotropy evolution.Furthermore, extension of this model to predict the anisotropic lateral oxidation of other mesa shapes as well as work on the inverse problem (i.e.determining the shape of the mesa to be etched to obtain a specific aperture) would ideally complete this study to provide an extensive predictive tool for oxide-based device simulation and manufacture.

Fig. 3 .
Fig. 3. Impact of the diffusion coefficient and reaction rate on the dynamics of the oxidation of a 30-µm-radius mesa.The dashed curves represents the dimensions along the <110> direction while the solid lines are taken along the <100> direction.

Fig. 4 .
Fig. 4. Impact of the degree of anisotropy on the dynamics of the oxidation of a 30-µm-radius mesa.The dashed curves represents the dimensions along the <110> direction while the solid lines are taken along the <100> direction.Inset: Calculated difference between the <100> and <110> oxidation lengths at t=32 min (L~Rmesa/2) as a function of the anisotropy coefficient, b.

Fig. 5 .
Fig. 5. Temporal dynamics of the experimental anisotropic oxidation of a 30-µm-radius mesa fitted using anisotropic model.

Fig. 7 .
Fig. 7. Angular distribution of the error between the modelled and extracted experimental oxidation fronts after 39 min 25 s of oxidation.

Fig. 8 .
Fig.8.The root-mean-square error between the modelled and extracted experimental oxidation fronts as function of oxidation time.