Impact of hydrogenation on the stability and mechanical properties of amorphous boron nitride

Interconnect materials with ultralow dielectric constant, and good thermal and mechanical properties are crucial for the further miniaturization of electronic devices. Recently, it has been demonstrated that ultrathin amorphous boron nitride (aBN) films have a very low dielectric constant, high density (above 2.1 g cm−3), high thermal stability, and mechanical properties. The excellent properties of aBN derive from the nature and degree of disorder, which can be controlled at fabrication, allowing tuning of the physical properties for desired applications. Here, we report an improvement in the stability and mechanical properties of aBN upon hydrogen doping. With the introduction of a Gaussian approximation potential for atomistic simulations, we investigate the changing morphology of aBN with varying H doping concentrations. We found that for 8 at% of H doping, the concentration of sp 3-hybridized atoms reaches to a maximum which leads to an improvement of thermal stability and mechanical properties by 20%. These results will be a guideline for experimentalists and process engineers to tune the growth conditions of aBN films for numerous applications.


Introduction
As the semiconductor industry downscales integrated circuits and power consumption increases, materials and reliability are pushed to their limits.It is thus increasingly important to either improve the performances of existing materials and/or develop new materials to meet these stringent demands of power reduction of electronic devices and circuits.While transistors have gone through several generations of design and the introduction of new materials, the back-end-of-line interconnects have seen fewer changes.Nevertheless, significant efforts have been dedicated to address the reduction of the resistance-capacitance (RC) delay [1,2,3].The resistance-capacitance (RC) reduction could be achieved in several ways, 1) decrease the capacitance density by decreasing the dielectric constant of interlayer dielectric (ILD) and the intermetal dielectric (IMD), 2) decrease the resistivity of the metal interconnect wiring, and 3) increase the cross-section of the interconnects.However, there are several problems associated with decreasing RC time delay as per the above approaches, 1) new materials with lower dielectric constants exist but may not be stable within the required process flow, that is they may not be mechanically and thermally stable as well as being good diffusion barriers, and 2) it is difficult to find metallic systems with electrical conductivity lower than Cu that can also meet the stability requirements of devices [4,5,6].A potential new barrier dielectric has recently emerged, amorphous boron nitride (α-BN).Experimental reports on α-BN indicate that it has a low dielectric constant, kvalues lower than 2, and exhibits higher stability and mechanical properties compared to other low dielectric materials such as organic polymeric materials [7,8,9].In addition, theoretical predictions suggest that a certain density of carbon content improves the structural and thermal properties of α-BN:C [10].More recently there is also a study on the use of α-BN for interconnect capping layers [11].This last work reports on the plasma-enhanced chemical vapor deposition (PECVD) of 3 and 7 nm α-BN as a capping layer to replace PECVD-grown Si 3 N 4 .The study finds that α-BN is an excellent insulator with efficient barrier against Cu diffusion, has good adhesion to copper and SiO 2 , is thermally stable and has a much lower dielectric constant (k=3), than Si 3 N 4 (k∼7) enabling an RC-delay reduction of 10-17%.In Ref. [10], using machine learning techniques and classical molecular dynamics, we explored the effects of C content on the physical properties of α-BN in an attempt to create an even more stable dielectric [10].However, neither these results nor the previous studies on ALD-grown α-BN [12,13], report on the effects of C and H content on the properties of the films.Also, it is well known that PECVD-grown Si 3 N 4 can contain large amounts of hydrogen [14] and it is likely that PECVD-grown α-BN may also contain hydrogen.Recently, Jacquemin et al. [15] showed that in the absence of hydrogen atoms, bromide ions can block the formation of sp 2 -hybridized atoms and produce continuous and stable α-BN films in a PECVD reactor using the BBr 3 and N 2 precursor.Hong et al. [7] measured 5.5% of H contamination to their PECVD grown α-BN when they used borazine as a precursor.Nonetheless impact of H on the stability or mechanical properties of α-BN is not clear.Therefore, it is important to understand the stability of hydrogen in α-BN as it is for Si 3 N 4 , since hydrogen can impact the underlying performance of the Si transistors and affect the dielectric and physical properties of interconnects.In this paper, we clarify the effects of hydrogenation of α-BN:H on its thermal stability and mechanical properties.From a computational perspective, atomistic calculations represent a suitable tool to describe complex structures, giving access to details at the atomic and molecular level to enable a basic understanding of new materials without performing more costly and time-consuming experiments.The extremely disordered nature of amorphous materials requires a computational approach able to capture the interatomic potentials in arbitrary complex local environments, a challenge that can only be tackled with machine learning-based methods [16,17,18].More specifically, classical molecular dynamics with the employment of force fields derived using machine learning and abinitio techniques constitute a powerful methodology to describe disordered material while keeping first-principles accuracy [19,20,21,22].In our theoretical study, we found that hydrogen influences the hybridization of the core structure (sp 2 /sp 3 ) and increases structural disorder as observed through radial distribution function (RDF).Further, in comparison to carbon doped α-BN, whose thermal stability and Young's modulus increase monotonically with C doping, hydrogen doping leads to an α-BN with a nonmonotonic increase in these properties, in fact, they peak at around 8 at% hydrogen.These results are critically important in providing directions to the experimentalist in tuning the deposition processes to meet the electronic device requirements.

Gaussian Approximation Potentials (GAP)
Structures for training and validation sets are generated using DFT calculations using the Quantum Espresso package [23,24,25] with Perdew-Burke-Ernzerhof (PBE) [26] exchange-correlation functional and projector-augmented wave (PAW) pseudopotentials.The energy cutoff for the wavefunction and electronic density is 75 Ry and 600 Ry, respectively.Both training and validation sets contain sufficiently large data sets of forces, energies, and stresses.The mentioned data sets involve both crystalline and amorphous BN structures, α-BN:H samples with different levels of H concentration.The density of the structures in the data sets ranges between 1.0 g/cm 3 and 3.0 g/cm 3 .Moreover, they also contain several distinct molecular configurations (H 2 , N 2 , ammonia, ammonium, and borazine) and isolated B, N, and H atoms.The final training and validation sets contain 2500 and 1800 samples, respectively.Such a wide variety of atomistic configurations enables us to model α-BN:H samples with better accuracy.The parameters shown in Table 2.1 has been employed to train the GAP model for α-BN:H using the training database.Smooth Overlap of Atomic Positions (SOAP) descriptor [27] is introduced to model the many-body interactions between atoms, while 2b and 3b descriptors are adopted for two-body and three-body interactions.The model is trained using the total energy function (Eq. 1) and local energy contributions, where 2b, 3b and MB represent the two-and three-body and many-body descriptors, where δ d is a scaling parameter which represents the contribution of each term, ϵ is the local energy contribution, K is the kernel function, q is the training configuration and α t is the all fitting coefficients [19].
Uniform sparsification has been used for 2b and 3b terms, while the CUR method [28] has been chosen for the SOAP kernel.After the training, we compare the energies and forces of structures in both training and validation sets obtained from molecular dynamics simulations with GAP model (GAP-MD) with the results from DFT calculations in order to evaluate the accuracy of the generated GAP model as shown in Fig. 1.A significantly low root mean squared error (RMSE) is calculated for both training and validation sets.

Melt-Quench Protocol for Sample Generation
One of the most common strategies used in molecular dynamics (MD) simulations to generate amorphous samples is the melt-quenching protocol.In this method, the sample is first melted by heating above the melting temperature and then rapidly quenched [29].The α-BN:H samples containing varying amounts of hydrogen are generated following this protocol using GAP-MD simulations with Large-scale Atomic/Molecular Massively Parallel Simulator(LAMMPS) code [30].Each sample has 10000 atoms and an equal number of boron and nitrogen atoms with varying amount of H.Each edge of cubic simulation boxes is ranging between 40 -45 Å depending the H concentration to keep the initial density of samples same.First, all boron, nitrogen, and hydrogen atoms are placed randomly in the simulation cell; then, the melted samples are equilibrated at 5000 K for 50 ps (timestep of 0.1 fs) using a Nosé-Hover thermostat.Later, the tem- perature of the samples is reduced to 2000 K in 100 ps (with a cooling rate of 40 K/ps) and equilibrated at this temperature.Samples are then cooled down to 300 K in 150 ps.After a short relaxation (10 ps) run, we also applied an annealing step where the temperature was increased to 700 K with a heating rate of 20 K/s and decreased to 300 K with a cooling rate of 5 K/s.Finally, annealed samples are relaxed and equilibrated at 300 K for 50 ps in the NPT ensemble.

Morphology Analysis of α-BN:H
We first present the analysis of the morphology of α-BN:H with different H concentrations employing the melt-quench protocol.A subset of the samples generated in this work is shown in Fig. 2. The RDF in Fig. 3 shows the density of surrounding atoms as a function of distance and gives insight into the crystallinity of the material.
Even though the first two peaks (≤ 4 Å) are clearly identified, no peak can be recognizable for longer distances, indicating the lack of long-range order.Hence, the amorphous character of α-BN does not change with the H doping.
The short-range order in α-BN:H is dominated by the first-neighbor distance, which contributes to the first peak in the RDF located at an average distance of ∼ 1.42 Å, as shown in Fig. 3 As α-BN:H is formed by the melt-quenching method from high temperatures, the resulting microstructure frozen in from the quenching process has a large influence in the hybridization of atoms in formed film.Understanding how hybridization can be changed with the fabrication conditions can help us tailor the properties of the material.The coordination number of α-BN:H is calculated to determine the type of hybridization.Fig. 4 shows the ratio of sp 2 (having coordination number 3) and sp 3 (having coordination number 4) hybridized atoms.sp 1 -hybridized atoms (having a coordination number of 2) are also presented.With the introduction of H atoms to α-BN, the ratio of sp 2 /sp 3 drops rapidly with a minimum observed at 8 at% H concentration.A deeper understanding of the chemical composition of the samples can be obtained by investigating the number and nature of chemical species involved in the bonds of the samples as a function of H concentration.As shown in Fig. 4, while B-H and N-H bonds are increasing, others seem to be decreasing monotonically with larger H concentrations.No H-H bonds are observed up to 20 at% of H concentration.While having more H atoms reduces the total mass of the sample, canceling the sp 1 hybridizations and having shorter bonds cause a dramatic shrinkage in the volume of the cell down to an H concentration of 8 at%.Due to this interplay, the density of α-BN:H at low H concentration levels is increased from 2.17 to 2.181 g/cm 3 , however, at larger H concentrations, the density of α-BN:H samples drops rapidly, as low as 2.01 g/cm 3 .In amorphous structures, such as α-BN and α-Si, some dangling bonds, vacancies and voids can be observed.H incorporation to α-Si can reduce the amount of dangling bonds and vacancies at small concentrations.Smets et al. [32] showed that when the H concentration is lower than 14 at%, then H atoms bond with other Si atoms around Si divacancies.However, at higher concentrations, H atoms disrupt the Si network, reduce the density of α-Si:H structure and increase the size and number of voids within the structure [33,34,35].Our simulations show a similar behavior for α-BN.While α-BN amount of sp 1 -hybridized atoms is reduced significantly with H incorporation lower than 8%, some nanovoids occur and amount of sp 1 -hybridized atoms increase.Moreover, a sharp decrease in the amount of B-N bonds can be observed after H concentration of 12%.This indicates that H atoms start to disrupt the B-N network and create these nanovoids, altering the mechanical properties and thermal stability of α-BN:H samples.

Thermal Stability of α-BN:H
Upon investigating the morphology of the generated samples, we also calculate the diffusivity of samples as a function of sample temperature.The diffusivity (as in described in Eq. 2) of samples at any given temperature can be extracted from the mean square displacement (MSD) of atoms where the MSD shows the average mobility of particles, as shown in Eq. 3 where r is position of the atom at any moment.The diffusivity of a sample is zero when MSD approaches a non-zero constant and has a zero slope.However, when the sample under investigation experiences a structural rearrangement and loses its stability, MSD has a non-zero slope.This allows us to evaluate the thermal stability of the samples and understand when they become unstable.Here, thermal sta-bility refers to the material's ability to retain its structural integrity without significant atomic diffusion or rearrangement and loss of short-range order when subjected to high temperatures.In order to assess the thermal conductivity, we calculate the stability of samples between 300 K to 3000 K.We first calculate the MSD and diffusivity at that temperature for 70 ps and then increase the temperature of the samples by 50 K in 30 ps, later we calculate the MSD at the new temperature.The time intervals are determined to be large enough to obtain statistically meaningful data.The NPT (constant number of atoms, pressure, and temperature) ensemble with a Nose-Hoover thermostat has been applied with a timestep of 0.25 fs.
The diffusivity of H atoms in α-BN:H is presented in Fig. 5. Non-zero values at low temperatures are obtained due to the vibration of atoms.Larger MSD values, which indicate structural rearrangement of atoms and unstable structure, are observed starting from 1600 K for 3 at% and 5 at% H-doping.For the case of 8 at% doping, diffusivity is near zero up to 1800 K.For larger H doping this temperature value drops significantly.
For highly H-doped α-BN, diffusivity becomes non-zero between 1000-1200 K.However, B and N atoms become diffusive at higher temperatures compared to H atoms.While H atoms begin to diffuse around 1800 K, B and N atoms still have near-zero diffusivity until 2000 K. Similarly at higher H concentration levels, there is a significant difference between the temperature at which H atoms and B and N atoms start to diffuse.At low-level H doping, this temperature difference is quite low.Regardless of the identity of atoms, a monotonic trend between thermal stability and H concentration levels is observed.α-BN:H samples become more stable until the H concentration reaches 8 at% due to the increase in sp 3 -hybridized atoms and reduction in sp 1 -hybridized atoms.Thereafter, thermal stability decreases rapidly since larger H doping reduces the density and causes some red nanovoids within the sample.At H concentrations larger than 20 at%, samples become unstable.At low hydrogen levels, H doping can lead to a more thermally stable structure since it reduces the amount of sp 1 -hybridized B and N atoms, reduces the number of dangling bonds, and increases the density.However, at larger concentrations, H atoms disrupt the stability of the α-BN samples.Kaya et al. [10] showed that the thermal stability of C-doped α-BN samples can be improved with larger amounts of sp 3 -hybridized atoms.Similarly, Cloutier et al. [36] and Liu et al. [37] showed that increasing the density of α-C films and number of sp 3 -hybridized atoms can improve the stability of films.

Mechanical Properties of α-BN:H
To calculate the mechanical properties of α-BN:H samples as a function of H concentration, we compute the elastic constant of the samples by using the stress fluctuations and Born matrix (i.e., the second derivatives of energy with respect to the strain) in a NVT ensemble at 300 K [38].For each sample generated in this study, we calculate the full elastic stiffness tensor C ij using the LAMMPS.Later, we calculate Young's modulus, shear modulus, bulk modulus, and Poisson's ratio from the stiffness matrix.  of α-BN films changes with the doping, we compare the temperature at samples lose their thermal stability of B and N atoms, and Young's modulus of α-BN:H and α-BN:C reported by earlier [10].The temperature values presented in the Fig. 7 show the approximate temperature where the diffusivity of B and N atoms becomes a non-zero value.Since C atoms can bond with four atoms instead of one and are heavier than H atoms, they lead to more sp 3 bonds and denser samples than α-BN:H.This leads to a higher Young's modulus and more stable structures, even at higher doping values (20 at%).Variation in stability and mechanical properties due to doping shows that the fundamental properties of α-BN films can be altered in different fabrication conditions.

Conclusion
Excellent mechanical properties and ultralow dielectric constant of α-BN opens a new pathway for microelectronics and neuromorphic computing technologies.This study reveals how H doping tunes the morphology, stability, and mechanical properties of α-BN.We first develop a machine-learning interatomic potential for H dopant in α-BN and performed GAP-driven MD simulations to generate realistic structures.Thanks to the accurate machine learning approaches, we show that thermal stability and mechanical properties of α-BN are improved by small H doping levels.Despite α-BN:H's extraordinary properties, it is crucial to perform a thorough benchmark analysis on growth conditions and corresponding properties.The results obtained in our study will provide a guide to process engineers to optimize the growth conditions to achieve the optimum materials performance in the context of microelectronics.

Figure 1 .
Figure 1.Scatter plots of the energies and forces obtained using newly generated GAP model of the corresponding DFT values computed over both training and validation datasets.
. A closer look at the first peak reveals an increased broadening and a shift to the left side with increased doping concentration, suggesting that it is induced by the presence of hydrogen atoms.Such change occurs due to the formation of chemical bonds different from B-N (B-H, N-H, B-B, and N-N), whose average bond length is different from B-N.The average bond lengths of these bonds are approximately ∼ 1.21 Å, ∼ 1.05 Å, ∼ 1.81 Å, and ∼ 1.44 Å, respectively.

Figure 2 .
Figure 2. A subset of generated α-BN:H samples with varying H concentration using VESTA[31] software where B, N, and H atoms are represented as black, grey, and black, respectively.

Figure 3 .
Figure 3. Average RDF of melt-quenched α-BN:H samples.All lines are averaged over five samples.

Figure 4 .
Figure 4.The ratio of sp 2 and sp 3 -hybridized atoms of α-BN:H as a function of H concentration (left).Number of the chemical bonds with respect to the H concentration (right).

Figure 5 .
Figure 5. Diffusivity of H atoms as a function of temperature for samples of α-BN:H with varying H concentration. Inset: Diffusivity of H atoms at 2200 K as a function of level of H doping.

Fig. 6 showsFigure 6 .
Figure 6.Mechanical properties of α-BN:H with respect to the H concentration.

Figure 7 .
Figure 7.Comparison of the temperature that the α-BN:H and α-BN:C samples lose their stability and their Young's modulus with respect to the doping concentration.

Table 1 .
Parameters used to train the GAP potential for H-doped α-BN.