Doping Engineering for Optimizing Piezoelectric and Elastic Performance of AlN

The piezoelectric and elastic properties are critical for the performance of AlN-based 5G RF filters. The improvement of the piezoelectric response in AlN is often accompanied by lattice softening, which compromises the elastic modulus and sound velocities. Optimizing both the piezoelectric and elastic properties simultaneously is both challenging and practically desirable. In this work, 117 X0.125Y0.125Al0.75N compounds were studied with the high-throughput first-principles calculation. B0.125Er0.125Al0.75N, Mg0.125Ti0.125Al0.75N, and Be0.125Ce0.125Al0.75N were found to have both high C33 (>249.592 GPa) and high e33 (>1.869 C/m2). The COMSOL Multiphysics simulation showed that most of the quality factor (Qr) values and the effective coupling coefficient (Keff2) of the resonators made with these three materials were higher than those with Sc0.25AlN with the exception of the Keff2 of Be0.125Ce0.125AlN, which was lower due to the higher permittivity. This result demonstrates that double-element doping of AlN is an effective strategy to enhance the piezoelectric strain constant without softening the lattice. A large e33 can be achieved with doping elements having d-/f- electrons and large internal atomic coordinate changes of du/dε. The doping elements–nitrogen bond with a smaller electronegativity difference (ΔEd) leads to a larger elastic constant C33.


Introduction
Piezoelectric materials, which can be applied to Radio Frequency (RF) filters, have drawn much attention with the commercialization of 5G communication technologies [1][2][3][4]. Aluminum nitride with wurtzite structure (w-AlN) is the prevailing piezoelectric material for the body acoustic wave (BAW) filters owing to the advantages of high acoustic velocity, minimal acoustic loss, high thermal stability, and good compatibility with Complementary Metal Oxide Semiconductor (CMOS) technology [5,6]. The critical parameters to evaluate the performance of piezoelectric materials for 5G filters are the mechanical quality factor (Q) and the longitudinal electromechanical coupling constant (k 33 2 ). The higher the Q, the lower the mechanical loss. The higher the k 33 2 , the larger the frequency bandwidth. In general, the Q value of 5G RF filters based on w-AlN thin film (Q = 400) is higher than that based on ZnO thin film (Q = 350), achieving low acoustic loss [7]. However, the k 33 2 (6.1%) [8] of undoped w-AlN is lower than some well-known piezoelectric materials, such as lead zirconate titanate perovskite (PZT) (k 33 2 = 8-15%) [8] and ZnO (k 33 2 = 7.5%) [8]; therefore, undoped w-AlN needs further optimization [9].
As shown in Equations (1) and (2), the characteristic Q and k 33 2 of a BAW RF filter are affected by the piezoelectric strain constant (e 33 ) and elastic constant (C 33 ) of the piezoelectric material [10][11][12] where Q, ε 33 s , ω, and η 33 are the acoustic quality factor, the clamped permittivity, the angular frequency, and the viscosity coefficient (details are shown in the support information) along the c-axis direction, respectively. A high C 33 is favorable to Q, and a high e 33 is favorable to k 33 2 . The piezoelectric material coupling coefficient k 33 2 and resonator effective coupling coefficient K eff 2 are positively related. It is not hard to design a resonator with a high K eff 2 from a material having a high k 33 2 value [13].Consequently, w-AlN should be tailored to have a high C 33 and e 33 , simultaneously, which has been proven to be a difficult task.
For example, first-principles calculations [14] and experiments [15] showed that añ 400% increase in the piezoelectric coefficient (d 33 ≈ e 33/ C 33 ) of w-AlN can be achieved with Sc doping. The increase in the e 33 is caused by the increase in the sensitivity of the internal atomic coordinates in response to the strain (du/dε) [16]. However, there also exists an elastic softening, owing to the elongated energy landscape in the c/a direction [17]. The e 33 of w-X a/2 Y a/2 Al 1−a N (X = Li; Y = V, Nb, Ta; a = 0.125, 0.25, 0.375) is enhanced compared to that of undoped w-AlN [18], while the C 33 decreases simultaneously due to the fact that these dopants can lead to a phase transition to a non-polar hexagonal structure. Hirata et al. [19] used first-principles calculations to investigate the enhancement in piezoelectric properties and the reduction in elastic properties by co-doping w-X a/2 Y a/2 Al 1−a N (X = Mg; Y = Nb, Ti, Zr, Hf; a = 0.125). The bonding analysis of the metal-nitrogen pairs by co-doping Mg + Y into w-AlN was carried out by the crystal orbital Hamilton population (COHP), which showed that weaker bonding energy is one of the reasons for the elastic softening.
The above results showed the need for a new mechanism to achieve a high C 33 and e 33 simultaneously. Manna et al. [20] found that the co-doping of Y and B elements in w-AlN improved the elastic properties while retaining good piezoelectric performance. Subsequently, Jing et al. [21] discovered that the C 33 of B 0.125 Sc x-0.125 Al 1−x N is higher than that of Sc x Al 1−x N with a small enhancement of the e 33 . These results confirm the feasibility of improving the piezoelectric and elastic properties by dual-element co-doping [22,23]. However, there is still a lack of systematic analysis leading to a clear strategy to choose doping elements for the enhancement of both the C 33 and e 33 . Therefore, expanding the map of doping elements and the understanding of the adjustment mechanism is critical to finding new doping schemes with excellent performance.
In this work, a high-throughput workflow is designed to calculate the piezoelectricity and elasticity of 117 X 0.125 Y 0.125 Al 0.75 N compounds. Filtered by the non-magnetic criteria, semiconductor criteria, stability criteria, and performance criteria, three dopants are finally screened out, which are B 0.125 Er 0.125 Al 0.75 N (e 33 = 2.11 C/m 2 , C 33 = 262.2 GPa, d 33  GPa, d 33 = 7.49 pC/N). It is found that the primary factor influencing the C 33 is the electronegativity difference (∆Ed) of the metal-nitrogen bonds, and the primary factor influencing the e 33 is the du/dε of the doping atoms. The bonds with a small ∆Ed in the doped-AlN between the doping elements and nitrogen with stronger strength leads to a larger elastic constant C 33 . The energy competition between the doping atoms and Al mainly affects the internal structural response (du/dε) of the crystal due to the transition elements doping into tetrahedral Al sites, tending to form non-tetrahedral coordinates, and undergoing excursions. The increasing of C 33 from the electronegativity difference and e 33 from the du/dε of the doping atom with d-/f-electrons provides clear ideas to design new piezoelectric materials for 5G filters.

Computational Details
The 2 × 2 × 2 supercells for w-X 0.125 Y 0.125 Al 0.75 N ( Figure 1b) were built with the special quasi-random structures (SQS) method [24]. The first-principles calculations were performed with the Vienna Ab initio Simulation Package (VASP) [25][26][27]. The Perdew-Burke-Ernzerhof (PBE) type generalized gradient approximation (GGA) as the exchangecorrelation function was implemented [24]. The elastic tensor was determined by performing the finite differences method. Six finite distortions of the lattice were taken, and the corresponding elastic constants could be derived from the strain-stress relationship [28]. The strains for the original structure along each of the Cartesian directions were ±0.5% and ±1%. The piezoelectric tensors were evaluated from the phonon and dielectric response calculations performed from the density functional perturbation theory (DFPT) [29][30][31]. The Monkhorst−Pack method [32] was used to set the k-point mesh. The k-grids used in the calculation of the structural optimization, self-consistent, and C ij /e ij were 30/L+1, 60/L+1, and 30/L+1, respectively, where L is the lattice constant of the systems. The cutoff energy of all calculations was 520 eV. The convergence criteria for the energy and force were set to 10 −4 eV and 10 −2 eV/Å, respectively. The Hubbard U values were from Wang et al. and Dudarev et al. [33,34]. The two-dimensional sandwich structure of the resonator and its geometric parameters is shown in Figure S1. The resonator consists of a piezoelectric material with top and bottom electrodes. COMSOL Multiphysics 6.0 is used to simulate the resonator quality factor(Q r ) and effective electromechanical coupling coefficient (K eff 2 ) of the resonator by using the finite element method [35]. Among them, the 2nd order Taylor approximation was performed to simulate the K eff 2 [36]. The Q r value was calculated using the method proposed by Bode et al. [37]. The physical parameters of the materials utilized in the simulation are shown in Table S1.

Results
To explore the theoretical feasibility of doping engineering to obtain materials with a high performance of large e 33 and C 33 , 117 dopants of X 0.125 Y 0.125 Al 0.75 N without toxic elements were tested. As shown in Figure 1a, the orange, green, blue, and gray spheres indicate X, Y, Al, and N, respectively. Considering the charge conservation law, the reasonable elements X and Y are substituted to the Al sites by 1:1. Moreover, Sc, Y, La Er, B, Ga, and In elements can be doped into either the X site or Y site due to the valence of +3. To effectively screen the piezoelectric and elastic performance of X 0.125 Y 0.125 Al 0.75 N materials, a high-throughput workflow was designed ( Figure 1c). First, the entries with complex magnetism were removed due to the difficulties to accurately calculate the properties of the magnetic materials for the high-throughput method. Second, the non-semiconductor systems were removed. If the band gap of X 0.125 Y 0.125 Al 0.75 N is less than 0, it indicates that the system is metallic and is not suitable for making piezoelectric layers for 5G filters. Then, the mechanical criterion was tested by the Born-Huang criteria of hexagonal structures [38]: It is clear that all of the models we considered were mechanically stable, and the detailed results are listed in Table S2. (Details can be found in Figure S2).
The detailed results of the 67 mechanically stable dopants are shown in Table 1. The modulation ranges of the e 33 and C 33 are 0.064~2.408 C/m 2 and 165.556~396.671 GPa, respectively. Table 1 shows that the e 33 of the dopants with small atomic radii elements and transition elements is high. The e 33 of the dopants with large atomic radii, such as K, Rb, Ca, Sr, Ba, and La, is smaller than that of those with small atomic radii, such as Li and Mg. Furthermore, the dopants that have one small radii element and one transition element (e.g., Mg co-doped with Ce, Ti, Hf, and Zr) show a higher e 33 than Mg co-doping with carbon group elements (i.e., C, Si, Ge, Sn, and Pb). For the C 33 , when the difference between the electronegativity of the doping atom and the N element is small, the C 33 is always high.

Analysis of Elastic Properties
As shown in Figure 2, we explored in detail the mechanism of co-doping to enhance the characteristics of the C 33 and e 33 , respectively. The hardness of the crystal is positively related to the bond density and negatively related to the ionicity indicator f i [42][43][44]. Figure 2a is the relationship of the C 33 and the electronegativity difference ∆Ed, where E X , E Y , and E N are the electronegativity of elements X, Y, and N, respectively. The electronegativity difference indicates the ionicity indicator (f i ) of the chemical bonds according to the Pauling for AB-type compounds [45], where f i indicates the degree of ionization of the hybrid bonds with a larger f i indicating that the chemical bond is closer to an ionic bond. Figure 2a shows that the C 33 is negatively related to the ∆Ed (i.e., the smaller the difference of electronegativity, the smaller the f i and the larger the C 33 ). Moreover, other factors, such as the bond density induced by lattice distortion, also slightly influence the C 33 . A specific mechanistic explanation of the effect of lattice distortion on the C 33 can be found in the supporting information. Generally, the smaller the electronegativity difference, the smaller the degree of ionization of the metal-N in X 0.125 Y 0.125 Al 0.75 N and the larger the hardness of the crystal. Thus, the electronegativity difference could be a criterion for the selected doped-AlN with a high C 33 .  Figure 2b shows the distribution of the e 33 , which comprises an electronic-response part and ion-polarization part [46].

Analysis of Piezoelectric Properties
clamped + e 33 non_clamped (5) e 33 clamped represents the electronic response under strain, which is evaluated by fixing the internal atomic coordinates at their equilibrium positions. e 33 non_clamped represents the ion polarization under strain, which is derived from the internal atomic coordinate changes. The mean and standard deviation of the e 33 non_clamped are 2.001 and 0.689, respectively. However, the mean and standard deviation of the e 33 clamped are −0.435 and 0.077, respectively. Obviously, the e 33 non_clamped mainly contributes the e 33 of w-AlN, owing to wider adjustable values and larger weights. Here, we focus on the derivation of the ion-polarization part, where n runs on all atoms in the supercell, e is the elementary charge, and a is the equilibrium lattice constant. Z 33 is the c-axis component of the dynamic Born charge tensor, and du/dε is the strain sensitivity. u is the ratio of the length of the metal-N along the c-axis (uc) to the lattice constant c in w-AlN (Figure 2c), which can be changed by the strain in the c direction. du/dε is the factor about the c-structure change, and Z 33 is the factor about the piezoelectric polarization variation on the structure change. Based on the first-principles calculation, the average Z 33 is 2.77 and can be adjusted from −6.48% to 8.53%; the average du/dε is 0.17 and can be adjusted from −90.89% to 27.10%. The variation of the du/dε is particularly large, which may significantly affect the e 33 non_clamped [16,47]. Figure 3a shows that there is a linear correlation between the du/dε along the c-axis and the e 33 . The du/dε of w-AlN is calculated by varying the doping atoms with an adjustment of the internal structure parameter, especially the structure parameter along the c-axis.  Figure 3b shows that the variation range of |du/dε| of the doping elements X and Y is much larger than that of Al and N. The average |du/dε| of the doping elements X and Y is 0.195 and 0.184, respectively, while that of the elements Al and N is only 0.0597 and 0.0836, respectively. Thus, the doping elements affect the e 33 dominantly compared to Al and N. The systems with large lattice distortion are doped by Sc, Y, and other transition elements with d-electrons and f-electrons.  Figure 4f). For the transition elements X or Y, such as Ti, Zr, Hf, Er, and Ta, they tend to format other non-tetrahedral coordination (e.g., octahedral coordination of Figure 4e). Octahedral coordination will compete against the tetrahedral coordination of the substituted Al and is more unstable than the tetrahedral coordination of Al. Figure 4c Figure 4d, for a non-transition element, the electron cloud of the regular tetrahedron geometry to bond to the nitrogen atom does not aggregate in the c-axis. For a transition element, it might bond to the nitrogen atom along the c-axis (Figure 4c). When a strain is performed on Mg 0.125 Ti 0.125 Al 0.75 N with unstable coordination, atoms move away from their regular tetrahedral positions and induce a larger du, which is due to the bond along the c-axis. As a result, the non-tetrahedral coordination of transition elements X or Y is easier to increase |du/dε| than the main group doping atoms with tetrahedral coordination under sp3 hybridization. It should be noted that the atomic radius also affects the e 33 . While the atomic radius of the doping atom is excessively large, it will produce a large local distortion in the lattice leaving a small space for an atom to move under the strain. For example, in Ba 0.125 Ti 0.125 Al 0.75 N, the atomic radius of Ba is 2.78 Å, and the du/dε is only 0.139. In a word, a small atomic radius and d/f-electrons are two parameters for finding doped-AlN with a large e 33 . As shown in Table 2, the Q r value of all three selected systems is higher than that of Sc 0. 25

Conclusions
Based on the high-throughput workflow, more than 117 X 0.125 Y 0.125 Al 0.75 N compounds were examined. In addition, B 0.125 Er 0.125 Al 0.75 N, Mg 0.125 Ti 0.125 Al 0.75 N, and Be 0.125 Ce 0.125 Al 0.75 N were screened out as having a higher e 33 , C 33 , and d 33 than Sc 0.25 Al 0.75 N. The Q r of the resonators made with these three systems was higher than that of Sc 0.25 AlN. The effective coupling coefficient (K eff 2 ) of B 0.125 Er 0.125 AlN and Mg 0.125 Ti 0.125 AlN was also higher than that of Sc 0.25 AlN, except for Be 0.125 Ce 0.125 AlN due to the high permittivity. The C 33 is affected by the electronegativity difference. There is a negative correlation between the ∆Ed and C 33 . The doping elements-nitrogen bond with a small ∆Ed leads to a larger elastic constant C 33 of the doped-AlN because the strength of the bond is stronger. The e 33 is affected by the du/dε of the doping atoms. The large du/dε comes from the competition between the tetrahedra coordinates [AlN4] of w-AlN and the non-tetrahedra coordinates of the doping elements with d-/f-electrons. This work provides a new way to find promising doped-AlN materials for 5G filters.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/ma16051778/s1, Figure S1: The two-dimensional sandwich structure of the resonator; Figure S2: The calculated and experimented e 33 , C 33 , and d 33 of Sc x Al 1−x N (x = 0~0.5); Figure Table S1: Physical parameters of the materials utilized in the simulation; Table S2: Dopants considered in this study and the result of C 11 -C 12 , 2C 13 2 -C 33 (C 11 + C 12 ), C 66 . Reference [48] are cited in the supplementary materials.