Effect of Si doping on the thermal conductivity of bulk GaN at elevated temperatures – theory and experiment

The effect of Si doping on the thermal conductivity of bulk GaN was studied both theoretically and experimentally. The thermal conductivity of samples grown by Hydride Phase Vapor Epitaxy (HVPE) with Si concentration ranging from 1.6×1016 to 7×1018 cm-3 was measured at room temperature and above using the 3ω method. The room temperature thermal conductivity was found to decrease with increasing Si concentration. The highest value of 245±5 W/m.K measured for the undoped sample was consistent with the previously reported data for free-standing HVPE grown GaN. In all samples, the thermal conductivity decreased with increasing temperature. In our previous study, we found that the slope of the temperature dependence of the thermal conductivity gradually decreased with increasing Si doping. Additionally, at temperatures above 350 K the thermal conductivity in the highest doped sample (7×1018 cm-3) was higher than that of lower doped samples. In this work, a modified Callaway model adopted for n-type GaN at high temperatures was developed in order to explain such unusual behavior. The experimental data was analyzed with examination of the contributions of all relevant phonon scattering processes. A reasonable match between the measured and theoretically predicted thermal conductivity was obtained. It was found that in n-type GaN with low dislocation densities the phonon-free-electron scattering becomes an important resistive process at higher temperatures. At the highest free electron concentrations, the electronic thermal conductivity was suggested to play a role in addition to the lattice thermal conductivity and compete with the effect of the phonon-point-defect and phonon-free-electron scattering.The effect of Si doping on the thermal conductivity of bulk GaN was studied both theoretically and experimentally. The thermal conductivity of samples grown by Hydride Phase Vapor Epitaxy (HVPE) with Si concentration ranging from 1.6×1016 to 7×1018 cm-3 was measured at room temperature and above using the 3ω method. The room temperature thermal conductivity was found to decrease with increasing Si concentration. The highest value of 245±5 W/m.K measured for the undoped sample was consistent with the previously reported data for free-standing HVPE grown GaN. In all samples, the thermal conductivity decreased with increasing temperature. In our previous study, we found that the slope of the temperature dependence of the thermal conductivity gradually decreased with increasing Si doping. Additionally, at temperatures above 350 K the thermal conductivity in the highest doped sample (7×1018 cm-3) was higher than that of lower doped samples. In this work, a modified Callaway model adopted for n-type GaN at high...


I. INTRODUCTION
In recent years the electronic and optoelectronic devices based on group-III nitrides, such as high-electron-mobility transistors (HEMTs), Schottky barrier diodes (SBDs), light-emitting diodes (LEDs), and laser diodes (LDs) have proved their potential as core elements for power electronics systems and solid state lighting technology.The performance and reliability of these devices, especially when operating at high temperatures and output powers, critically depends on the efficiency of heat removal from the active region of the device.Thus, the thermal management is an essential part of the device design, and the accurate knowledge of the materials thermal properties (thermal conductivity, thermal expansion, heat capacity) becomes of crucial importance.
The thermal conductivity of group-III nitride materials has been extensively investigated in the past but the reported experimental data varies greatly due to the differences in the quality of the samples examined and the specific challenges in the measurement methods applied.2][3] For bulk GaN, the advances in high-pressure high-temperature (HPHT), ammonothermal (AT) and hydride phase vapor epitaxy (HVPE) growth, has led to higher quality material with lower dislocation densities 4 and consequently larger values of the thermal conductivity have been reported.For example, room temperature thermal conductivities of 220 W/m.K 5 and 230 W/m.K 6 have been reported in a HPHT grown n-type GaN and in an AT grown semi-insulating GaN sample (co-doped with O and Mg), respectively.In free-standing undoped GaN samples grown by HVPE, thermal conductivities of 227 W/m.K 7 and 253 W/m.K 8 have been measured by the steady-state heat flow method and laser flash method, respectively.For semi-insulating (Fe-doped) HVPE grown GaN with a dislocation density as low as 5×10 5 cm -2 , Mion et al. 9 has demonstrated a thermal conductivity of 230±15 W/m.K measured by 3ω-method.More recently, a thermal conductivity of ∼294±44 W/m.K was reported using the laser flash method for unintentionally doped HVPE GaN. 10 Theoretical studies of the thermal transport in bulk GaN by first-principles calculations 11,12 and molecular dynamics simulations 13,14 have predicted an upper limit for the room temperature thermal conductivity in perfect isotopically pure GaN to be between 185 W/m.K 13 and 404 W/m.K. 12 These values represent the intrinsic lattice thermal conductivity which is only governed by the phonon-phonon scattering.In crystals with a naturally occurring Ga and N isotope concentrations the thermal conductivity is expected to be by 35-40% lower due to the phonon-isotope scattering. 13,14part from the "isotope" effect the reduction in the ideal intrinsic thermal conductivity in GaN is mainly determined by the presence of point defects (intrinsic or extrinsic) and extended defects (dislocations and stacking faults).Indeed, the large variation of the experimental data reported so far is a result of the fact that samples of different defect density were measured.While the role of dislocations on the thermal conductivity of bulk and heteroepitaxial GaN have been well investigated both experimentally 9,15 and theoretically, 16,17 the effect of point defects is not completely understood.Previous studies of doped GaN have been done on material with a high dislocation density, 16,17 on samples grown by different methods, 18 or on samples doped with O and Mg, or with Fe. 6,9 There has been no systematic study of the effect of a single impurity on the GaN thermal conductivity, especially in bulk material.
Currently, there is an expanded use of native GaN substrates for commercial group-III nitridebased electronic and optoelectronic devices.The reason is not only the better crystal quality of the epitaxial device structures, but also the reduced thermal boundary resistance at the homoepitaxial interface in high-power devices working at higher temperatures. 19For different device applications, different GaN substrates are needed -undoped, n-type doped, or semi-insulating.Thus, for proper thermal management, an accurate knowledge of the effect of doping on the thermal conductivity of GaN is required.A detailed study of the role of Si doping is highly needed since Si is the main dopant for achieving n-type conductivity in GaN.
In our previous study, we have examined the thermal conductivity of Si-doped bulk GaN grown by HVPE using 3ω method measurements. 20The slope of the temperature dependence of the thermal conductivity above the room temperature was found to gradually decrease with increasing Si doping.Additionally, at temperatures above 350 K, the measured thermal conductivity in the samples with higher doping was higher than that in lower doped samples.In this work, we extend our study by measuring the thermal conductivity in low-dislocation-density bulk GaN samples with controllable Si doping in a wider temperature range.Also, we develop a model for the thermal conductivity of GaN material with high free electron concentrations at higher temperatures.The model is based on the Callaway formalism 21 with taking into account recently updated and more correct material parameters of GaN, as well as with a refined treatment of both the phonon-point-defect scattering and the phonon-free-electron scattering as compared with the previously studies.A reasonable match between the measured and the simulated thermal conductivity in the examined GaN samples was obtained.

A. Sample growth and characterization
The samples studied were cut from bulk-like GaN grown by HVPE in a vertical reactor by Kyma Technologies Inc.The growth was done along the [0001] GaN crystallographic direction on sapphire substrates.The typical thickness of GaN layers, self-separated from the sapphire during cooling down, was about 1 mm.The lateral dimensions of the samples used in thermal conductivity measurements ranged from 5×5 mm 2 to 10×10 mm 2 .
The silicon doping was achieved using silane (SiH 4 ) introduced into the reactor with nitrogen (N 2 ) carrier flow.The Si concentration was controlled by varying the silane flow in a range of 0-250 sccm.All other growth conditions, such as the growth temperature, N 2 carrier flow and III/V ratio, were purposely kept constant.Secondary ion mass spectroscopy (SIMS) measurements were used to determine the total Si concentrations, as well as some of the typical residual impurities in HVPE grown GaN -oxygen, hydrogen and carbon.The Si concentration was found to range from 1×10 16 cm -3 in the undoped sample to 7×10 18 cm -3 in the highest doped sample.Also, the SIMS data revealed that the concentrations of [O] ≈ 2.3×10 16 cm -3 , [H] ≈ 1.5×10 18 cm -3 , and [C] ≈ 7×10 15 cm -3 were unaltered by the varying the Si doping.The free electron concentration, n e , as determined by conventional Hall effect measurements at room temperature showed a nearly complete donor ionization, i. e. the n e scales closely with the Si concentration in an agreement with previous studies. 22The Ga vacancy density was studied by positron annihilation spectroscopy (PAS).In all samples, no vacancy signal was detected which leads to the conclusion that the Ga vacancy density was below 1×10 15 cm -3 (detection limit of the PAS).All samples exhibited a threading dislocation density of < 5×10 6 cm -2 as determined by etch pit density and cathodoluminescence.The data from SIMS, PAS and Hall measurements for all studied samples are summarized in Table I.

B. Thermal conductivity measurements by the 3ω method
The thermal conductivity was measured using the 3ω method. 23,24In the 3ω method, a thin metal wire with four contact pads are photo-lithographically deposited on the sample surface with a standard design. 24The wire serves as both a heater and a sensor, where the entire pattern is referred to as a 3ω device.For our measurements, several 3ω devices with wires of length 3 mm and variable width (4, 6, and 10 µm) were fabricated on the surface of each sample in accordance with the constraints of the boundary conditions required for the sample and frequencies used. 23Gold was chosen as a material for the metal wire due to its high temperature coefficient of resistivity (α T = 0.0031 K -1 ) and its immunity to oxidation.Due to the electrical conductivity of the samples measured, a 100-nm-thick insulating layer of SiO 2 was deposited on the surface of the GaN before the deposition of the metal wire.The measurements were performed in a temperature range of 295-470 K by mounting the samples on a temperature controlled plate.An additional thermocouple was fixed at the top of the plate in the vicinity of the sample in order to verify the temperature of the samples during the measurements.The temperature stability was controlled within ± 2 K by surrounding the sample setup in temperature isolating box and using long dwell times at each temperature.
The 3ω method is based on applying an alternating current at angular frequency ω along the wire and measuring the voltage drop at 3ω as a function of ω.The Joule heating of the wire generates a heat TABLE I. Impurity concentrations, Ga vacancy (V Ga ) concentrations, free electron (n e ) concentrations, and the resistivity of samples as determined by SIMS, PAS and Hall measurements.Sample SiH4 flow (norm) Si (cm -3 ) O (cm -3 ) H (cm -3 ) C (cm -3 ) V Ga (cm -3 ) n e (cm were α T is the temperature coefficient of resistivity, ∆T is the temperature oscillations, and V 1ω is the amplitude of the voltage drop at fundamental frequency ω.The precise value of α T for each heater/sensor device is experimentally determined by measuring the change in resistance over the entire temperature range used in the experiment.Two differential amplifiers were used to subtract the fundamental frequency component from the signal and to amplify the third harmonic, which is typically 1000 times smaller in magnitude than the first harmonic.
Provided that the thermal penetration depth in GaN is much smaller than the sample thickness and much larger than the wire width, the boundary conditions are satisfied, and the power normalized temperature variation can be approximated by 23

∆T
Here P is the power applied to the wire, l is the length of the wire, k is the thermal conductivity, and C is a constant independent of the frequency and the wire length.Since ∆T/P is directly proportional to V 3ω the thermal conductivity can be extracted from the slope of the V 3ω versus ln(ω) dependence.
Note that the linearity of this dependence observed for the entire frequency range used (50-4000 Hz) at different temperatures (see Fig. 1) is an experimental verification that the boundary conditions are fulfilled in the measurements.We note that the insulating SiO 2 layer on the surface of the samples does provide non-negligible thermal boundary resistance.However, tests have shown that this thermal boundary resistance simply gives a shift in magnitude of the temperature fluctuations of ∆T /P and does not change the slope of ∆T /P vs. ln(ω) used to determine thermal conductivity. 24The presence of the thin thermally resistive SiO 2 layer does not provide a frequency dependent shift in ∆T /P and therefore does not change the measured thermal conductivity.

III. THEORETICAL MODEL
The common treatment of the lattice thermal conductivity in solids is based on Callaway's phenomenological formalism, 21 where a Debye-like phonon spectrum is assumed and the various phonon scattering processes are accounted for by introducing frequency-and temperature-dependent relaxation times.During the years, several modifications and extensions of the original Callaway model have been proposed [26][27][28][29] and verified through comparison to experimental data for a variety of semiconductor materials, including GaN. 16,30,31 In our analysis, we are using the modified Callaway model where the contributions of longitudinal and transverse acoustic (LA and TA) phonons are treated separately and the thermal conductivity is given by 29,30 Here, k B is the Boltzmann constant, T is the temperature, is the reduced Planck constant, ν s is the acoustic phonon velocity, θ s D is the Debye temperature (determined from the maximum phonon frequency at the zone boundary), and

The variable of integration
x is related to the phonon frequency ω and is defined by x = ω/k B T .The summation is over all three acoustic modes in GaN -one longitudinal and two transverse, i. e. s = L, T1, T2.In Eq. ( 3) τ s N denotes the relaxation time for normal (N) phonon-phonon scattering, τ s R is the relaxation time for all resistive (R) scattering processes, which are taken to be additive, and τ s C is the combined relaxation time given by ( The resistive scattering processes considered here are the Umklapp (U) phonon-phonon scattering, the phonon-isotope (I) scattering, the phonon-pointdefect (PD) scattering, the phonon-dislocation (D) scattering, the phonon-boundary (B) scattering, and the phonon-free-electron (FE) scattering.Then, the resistive scattering rate is expressed as The various scattering mechanisms have different temperature and frequency dependencies which leads to different weights in the thermal conductivity in different temperature regions.This in turn results in the application of different approximations and simplified models when comparison with the experimental data is made, depending on the particular materials studied or the temperature region of interest.Here, we are following the modified Callaway model 21,29,30 without further approximations and the relaxation times are taken in their general forms but accounting for the specifics of the crystal structure, phonon spectrum, and defects in GaN.We concentrate our analysis in the temperature region of T ≥ 250 K where the measurements are done.The explicit scattering rate expressions for different phonon scattering processes are presented and discussed below.
Normal phonon-phonon scattering (N-scattering) is not resistive process, i. e. does not directly degrade the heat flow, however it has an essential role in the phonon distribution and thus influences all resistive processes.In many practical cases, e. g. non-perfect crystals at high temperatures, the N-scattering rate appears to be much smaller than the scattering rate of the resistive processes, (τ s N ) −1 (τ s R ) −1 , so that τ s C ≈ τ s R and therefore in the calculations of the thermal conductivity the second term in Eq. ( 3) is neglected. 16Several forms of the normal phonon scattering rate in solids are suggested. 32Generally, the N-scattering rate is given by 30,32 ( or in terms of the dimensionless variable x where M is the average mass of an atom in the crystal, V o is the volume per atom, and γ s is the phonon-mode Grüneisen parameter, which characterizes the anharmonicity of the crystal lattice.The values of a and b depend on the crystal symmetry, the type of phonon mode and the temperature range. 32In the case of a hexagonal crystals for TA phonons. 32,33mklapp phonon-phonon scattering (U-scattering) is a three-phonon process in which phonon momentum is not conserved.It is the main resistive process in defect-free crystals and leads to an exponential decay of the thermal conductivity at high temperatures.Based on the analysis of Slack and Galginaitis, 34 the U-scattering rate is expressed (in terms of the dimensionless variable x) as 30 ( The form of the U-scattering rate is assumed to be identical for all phonon modes.The difference in the contributions of longitudinal and transverse phonons comes only from the mode-dependent material parameters.
Phonon-isotope scattering (I-scattering) is essentially scattering due to the mass fluctuation in crystals where various isotopes of the constituent chemical elements are present.This scattering is shown to reduce the maximum thermal conductivity (typically occurring at T ≈ 0.05θ D ) by 30-50% and to affect even room-temperature values in various undoped and low-doped semiconductors. 29,30s shown by Klemens 35 the scattering of phonons by a defect with mass different from that of the host atom follows a Rayleigh scattering law, therefore the I-scattering rate is given by 30 ( Explicitly for GaN the strength of the phonon-isotope scattering Γ I can be written as 30 Here Γ(p) is the scattering parameter for a single chemical element (p = Ga, N): with where m p i is the atomic mass of ith isotope and c p i is the fractional atomic natural abundance.There are two Ga isotopes -69 Ga and 71 Ga with 60.108% and 39.892% abundance, respectively, and two N isotopes -14 N (99.632%) and 15 N (0.368%). 36Using Eqs.(8-10) and the exact atomic mass of isotopes 37 we calculate Γ I = 2.744×10 -4 .
Phonon-point-defect scattering (PD-scattering) originates from the disturbance of the lattice periodicity by the point defects (substitutional impurities or vacancies).The disturbance is generally due to (a) the mass difference between the impurity and host atoms, (b) the change in the force constant of the bonds in a crystal with defects, (c) the strain field around imperfection produced by the difference in the size of the defect and the host atom.The phonon-point-defect scattering rate also follows a Rayleigh scattering law and is expressed as 35 ( When several substitutional impurities are present in the crystal the strength of the PD-scattering Γ PD is given by 27,35 where f i and M i are the fractional concentration and the atomic mass of the ith impurity, G i is the average stiffness constant of the nearest-neighbor bonds between the ith impurity and host atoms, G is the average stiffness constant of the bonds in the host crystal, ∆R i is the average local displacement in the host lattice due to the ith impurity, and R is the average nearest-neighbor atomic radius in the host crystal.The displacement ∆R i is generally defined as the difference between the Pauling ionic radii of the impurity and the host atom.The constant Q equals 4.2 if the nearest-neighbor bonds of the impurity have the same anharmonicity as all other bonds, or Q = 3.2 if the anharmonicity of the nearest-neighbor bonds is excluded. 35The second term in Eq. 12 is usually neglected because of the lack of data about the change in the stiffness constants in the presence of a substitutional impurity. On the other hand, the Pauling ionic radii are known to depend on the charge state, the electronic spin state and the coordination number.Moreover, the ionic radii of atoms in a semiconductor crystal differ from the ionic radii of "free" atoms due to the partial or complete covalent character of the bonds. 38So, the third term in Eq. ( 12), which account for the phonon scattering by the elastic strain field, is in some way ambiguously defined.To refine the analysis, we set which reflects the fact that in most semiconductors the elastic constants depend on the atomic volume only 27 and Q = 4.2 (i.e. it is assumed that the impurities do not change the anharmonicity of the bonds 35 ).Then, to be consistent with the above approximation the displacement ∆R i is calculated from the difference between the average bond length in the host crystal and the average bond length in the presence of an impurity for a relaxed lattice.The vacancies are treated in the same way but with M i = 0, G i = 0, and Phonon-dislocation scattering (D-scattering) has two distinctive mechanisms: scattering on the core of dislocation lines and scattering by the elastic strain field on the dislocation lines.Taking into account the contributions of edge and screw dislocations, the total D-scattering rate is 16,35 ( Here The parameter P E depends on both the ratio between the phonon velocities for two polarizations, ν L /ν T , and the Poisson ratio. 35The first and the second terms in Eq. ( 13) represent the scattering on the core and the scattering by the elastic field, respectively.Such a treatment of the phonondislocation scattering predicts that the thermal conductivity of GaN is independent of the dislocation density for N D < 10 10 cm -2 . 16However, the analysis of experimental data for a variety of materials has showed that the Klemens formulation underestimates the phonon-dislocation scattering rates by a factor of 10 2 -10 3 . 39In fact, previous data for thick HVPE grown GaN has revealed that the dislocations start to affect the thermal conductivity for N D > 10 6 cm -2 . 9Remarkably, the thermal conductivity dependence on the dislocation density experimentally observed for N D > 10 6 cm -2 9 is very similar to that calculated through Eq. ( 13) for N D > 10 10 cm -2 . 16Therefore, a dimensionless pre-factor of 1000 is introduced in the right part of Eq. ( 13) in order to match the threshold value of N D . 17This does not affect our analysis because the dislocations density in all samples investigated is < 5×10 6 cm -2 .
Phonon-free-electron scattering (FE-scattering) is usually neglected in the analysis of the thermal conductivity because it is believed that its contribution is much smaller than other resistive scattering processes, especially in undoped semiconductors at low temperatures.However, for doped semiconductors at elevated temperatures (above the room temperature), where in the wide-bandgap materials most of donors are ionized, the role of the FE-scattering becomes essential.Within the effective mass approximation, the FE-scattering rate is expressed by 40

(τ
where Here, E 1 is the deformation potential for the electron-acoustic phonon interaction, m e is the electron effective mass, ρ is the mass density, and E F is the Fermi-level energy measured from the bottom of the conduction band.Phonon-boundary scattering (B-scattering) dominates at very low temperatures, where the phonon mean free-path becomes comparable with the dimensions of the sample under investigation.The phonon-boundary-scattering rate is independent of the phonon frequency and leads to T 3 dependence of the thermal conductivity below its maximum at T ≈ 0.05θ D .The B-scattering rate is usually expressed in terms of effective phonon mean path, L ph : 29 where ab is the cross section of the samples and d is the sample length in the direction of the heat flow.
Finally, the lattice thermal conductivity can be calculated by inserting the expressions for relaxation times for all scattering processes (Eqs.(5-16)) into Eq.(3).

IV. RESULTS AND DISCUSSION
All GaN samples were measured at identical conditions, namely the same dimensions of the 3ω devices and the same frequency range.At room temperature (T = 295 K), a thermal conductivity of k = 245±5 W/m.K was measured for the undoped sample (sample A).8][9] With increasing Si concentration, the thermal conductivity gradually decreased and for the sample with the highest doping (sample F), k = 210±6 W/m.K was found.Such a behavior is quite reasonable and can be simply explained by an increased contribution of phonon-point-defect scattering.At elevated temperatures (T > 295 K) the thermal conductivity for all samples decreased with increasing temperature (Fig. 2).However,  the rate of the decrease (i.e. the slope of the temperature dependence of k) was different for different Si concentrations, which will be discussed later.Usually, to analyze the temperature dependence of the thermal conductivity of solids and to distinguish between the contributions of different phonon scattering processes the experimental data is fitted by the Callaway 21 or Holland 26 models.In doing so, the coefficients in front of the frequency and temperature dependence in the scattering rate expressions for different scattering processes are used as fitting parameters. 29For undoped and isotopically pure materials with high crystalline quality, such an approach works very well because only the phonon-phonon scattering contributes to the thermal conductivity and few adjustable parameters are needed (e. g. the Nand U-scattering rate coefficients for longitudinal and transverse phonon modes).However, the presence of dopants, free carriers, and dislocations in most of the technologically important semiconductor materials requires all scattering mechanisms to be included in the analysis.Thus, the number of the fitting parameters increases considerably and the fitting procedure becomes more complicated.The Callaway model, with some approximations, has been successfully used to fit thermal conductivity data for various GaN samples. 6,19,30,31,41However, the extracted values for the scattering rate coefficients, which in principle only depend on the material parameters, are quite different.
To analyze our experimental data a slightly different approach was employed.The modified Callaway model (as described in Sec.III) was applied to fit the temperature dependence of the thermal conductivity, but instead of the numerous scattering rate coefficients, only the Grüneisen parameters for longitudinal and transverse phonons, γ L and γ T , were used as adjustable parameters.This approach is justified by the fact that all other material parameters in Eqs. ( 3)-( 16) are fairly well known for GaN and are summarized in Table II.The mode-dependent Grüneisen parameters represent the anharmonicity of the atomic interactions in the crystal structure and are defined as 42 where ω s is the phonon frequency, q s is the phonon wavevector, and V is the volume.Experimentally, the Grüneisen parameters can be obtained from the pressure dependence of the phonon frequencies.
Unfortunately, there is no such data for the acoustic phonon modes in GaN.For zinc-blende GaN,  ab initio lattice-dynamic calculations revealed γ s values in the range 0.4 -1.3 for different phonon modes and for some modes negative values were found. 42Mode-averaged Grüneisen parameters for wurtzite GaN of γ = 0.74 43 and γ = 0.85 44 have been extracted from the specific heat and the thermal expansion measurements, respectively.Due to the uncertainty in the suggested values we prefer to use the Grüneisen parameters as fitting parameters.The values for γ L and γ T listed in Table II are those obtained from the best fit of thermal conductivity for the undoped sample (sample A).Note that the same Grüneisen parameter is assumed for the T1 and T2 phonon modes.Another parameter which should be carefully chosen for a proper fitting of the thermal conductivity data is the Debye temperature.As pointed out by Morelli et al. 30 the calculation of the Debye temperature from the specific heat or from the acoustic phonon velocities results in a significant overestimation.Since the heat is carried by acoustic phonons with frequency up to the maximum frequency at the zone boundary, ω max s , the correct Debye temperature to be used in Eqs. ( 3) and ( 6) is given by θ s D = ω max s /k B .The Debye temperatures for the longitudinal and transverse phonons, θ L D , θ T 1 D , and θ T 2 D were calculated from the zone boundary frequencies along Γ-M and Γ-K directions of the Brillouin zone as measured by inelastic x-ray scattering. 45The phonon velocities were determined from the experimental elastic constants in bulk HVPE grown GaN. 46or the phonon-point-defect scattering, substitutional Si and C on Ga sites, substitutional O on N sites, and Ga vacancies were considered.The atomic mass of the impurities (listed in atomic mass units in Table II) were estimated with by considering the natural abundance of the corresponding isotopes. 36,378][49] In the treatment of the phonon-dislocation scattering we assumed a random orientation of the dislocation lines with respect to the temperature gradient, i. e. η = 0.55 35 and an equal density of screw and edge dislocations, i. e. N S = N E = N D /2.From the elastic constants of GaN 46 the parameter P E was estimated to be 1.1.For the phonon-free-electron scattering the Fermilevel energy E F at each temperature was calculated from the measured free electron concentration by an analytical approximation of the Fermi-Dirac integral. 50An isotropic electron effective mass of m e = 0.23m o was used, where m o is the free electron mass. 51The electron-acoustic phonon deformation potential was taken as E 1 = 9.1 eV. 52o test the model and the fitting procedure we first analyze the thermal conductivity data reported by Slack et al. 7 using the steady-state heat flow method for an undoped free-standing GaN sample grown by HVPE.The residual impurity concentration (Si and O) and the threading dislocation density in this sample were ∼10 16 cm -3 and ∼10 6 cm -2 , respectively.In other words, the sample is similar to our undoped sample (sample A).The experimental data and the fitted curves are shown in Fig. 3 with Grüneisen parameters γ L = 0.9 and γ T = 0.65.We note that the same experimental data analyzed by the Callaway model including only N-, U-, Iand B-scattering contributions and γ L = 1.0 and γ T = 0.7 were extracted from the fit. 30The small difference in the Grüneisen parameters can be attributed to the slightly different phonon velocities used in Ref. 30.As seen in Fig. 3 the LT model underestimates the thermal conductivity at T > 250 K.In this temperature region, a more satisfactory fit (using the same Grüneisen parameters) is achieved by the HT model.At lower temperatures, however, the HT model is obviously inappropriate.These results indicate that for the analysis of thermal conductivity data in our samples, which were measured at temperatures above 295 K, the HT form for the N-scattering rate should be used.
The best fit of the experimental data for undoped sample (Fig. 2 (a)) yields Grüneisen parameters γ L = 0.8 and γ T = 0.61.Using these values the temperature dependence of the thermal conductivity for the other four samples is calculated by the model described in Sec.III.Only the Si concentration (i.e. the fractional concentration f Si in Eq. ( 12)) and the free electron concentration (i.e. the Fermilevel energy in Eq. ( 15)) are varied.Figure 2 (a) shows that the calculated thermal conductivity of the low doped sample (sample B) matches the experimental data very well.The difference between the two curves is quite small which is not surprising given that a small amount of impurities (and free carriers) are not expected to significantly disturb the phonon propagation.In fact, the best fit of the data for sample B, gives the same Grüneisen parameters.For the moderately doped sample (sample C), the correspondence between the calculated and measured thermal conductivity is also rather good (Fig. 2(b)).With increasing doping, however, the simulated curves deviate from the experimental data especially at higher temperatures.For the highest doped sample (sample E), the calculations underestimate the measured thermal conductivity in the entire temperature range considered.In order to provide a better match, we attempted to fit the data for samples D and E in order to find new Grüneisen parameters.For example, the k values measured near the room temperature (290-310 K) for sample E can be reproduced by a slight reduction of γ L or γ T .Nevertheless, for both samples it was impossible to reach a satisfactory fit of the experimental temperature dependence.The experimental results clearly show that the higher the doping, the weaker the temperature dependence.In other words, the decrease in k with increasing temperature becomes smaller.As seen in Fig. 2 (b) the thermal conductivity of sample E is higher than that in lower doped samples at temperatures above 330 K.The calculations also indicate a change in the temperature dependence with Si concentration, but at each temperature the thermal conductivity is found to decrease monotonically with increasing doping.
Following the Callaway model, the change in the strength of a particular phonon scattering process cannot entirely explain the change in the behavior of the thermal conductivity with temperature.It is likely that the peculiarity of the temperature dependence is a result of the interplay between different scattering processes.In our analysis of the effect of Si doping on the thermal conductivity of GaN, only the rates for PDand FE-scattering are varied.We have further studied the role of the FE-scattering in more details.The thermal conductivity as a function of the Si doping at different temperatures is presented in Fig. 5.It is again evident that the inclusion of the FE-scattering in the model is indispensable for a proper explanation of the experimental data at T = 300-350K.Above 350 K, however, the contribution of the FE-scattering seems to decrease as seen by the proximity of the two simulated curves.At higher temperatures the thermal conductivity is predicted to be independent of Si concentration up to [Si] ∼ 1×10 18 cm -3 , while at room temperature k is only independent on Si concentration up to [Si] ∼ 5×10 16 cm -3 .The trend of the room temperature thermal conductivity show in Fig. 5 is generally consistent with previously reported data for thin Si doped GaN layers residing on sapphire. 16However, all values measured here were significantly higher due to the fact that the bulk GaN is of much better quality in terms of the structural defect density as compared to the thin heteroepitaxial GaN layers.
Although we were able to successfully model both the temperature dependence of k at low and to moderate doping levels and the Si concentration dependence of k at temperatures up to 350 K, the discrepancy at higher temperatures and higher doping remains puzzling.The discrepancy could be attributed to some shortcomings of the Callaway model, such as the simple Debye-like phonon dispersion or the independence of the phonon modes.In that case, however, one would also expect a discrepancy in the model for undoped materials, which is not seen here, nor was in previous studies. 30,31Since the strength of the PD-scattering increases with increasing Si concentration the thermal conductivity should gradually decrease independent of the temperature.Therefore, it is most likely the FE-scattering and its handling in the model are responsible for the deviation between the calculated and measured thermal conductivity values well above room temperature.Ziman 40 has showed that at higher temperatures the number of electrons available for scattering with phonons decreases due to the requirement for energy and momentum conservation, which means that the contribution of the FE-scattering to the thermal conductivity will diminish with increasing temperatures.This effect is accounted for in Eqs. ( 14) and ( 15) but likely the temperature dependence is not covered precisely because of the simplified treatment of the problem. 40Another plausible explanation for the weaker temperature dependence of thermal conductivity at higher temperatures and doping levers is the role of the free electrons in the heat transport which is not considered in the model.The electronic thermal conductivity was estimated to contribute to the total room temperature thermal conductivity of GaN at free carrier densities above 1×10 18 cm -3 . 53These estimates have been done for a material with dislocation density > 5×10 8 cm -2 , where a strong electron-dislocation scattering contribution is present.However, in the bulk GaN samples measured here the dislocation density is two orders of magnitude lower and one can expect a larger contribution of the electronic thermal conductivity especially above room temperature.Additional insight to the dominant scattering processes above room temperature can be obtained from a simple power law fit of the measured temperature dependence of the thermal conductivity, namely Here k o is the thermal conductivity at room temperature and T o = 295 K.This empirical fit could be useful for thermal management and device design since it predicts the thermal conductivity behavior both at different temperatures and different doping levels.An example for such a fit for all measured samples is shown in Fig. 6.The variation of the slope, α, with the Si concentration is depicted in Fig. 7.It is found that the slope gradually decreases with increasing doping from α = 1.3 for the undoped sample to α = 0.55 for the highest doped sample.The slope found for the undoped sample is consistent with previously reported numbers for high quality HVPE grown GaN -1.43 9 and 1.22. 30It is generally assumed that at high temperatures (above the Debye temperature) the thermal conductivity is determined by the Umklapp phonon-phonon scattering and the k is expected to vary as ∼ 1/T. 54t these temperatures, the contribution of all resistive scattering processes related to imperfection (I-, PD-, and D-scattering) is independent of the temperature. 54However, the Debye temperatures of GaN for longitudinal (439 K) and transverse (310 K) acoustic phonons are within the temperature range considered here.Thus, the above simple approximation does not hold.Moreover, as it is evident from Fig. 3  as ∼ 1/T and ∼ T 2 , respectively. 54A similar behavior of the temperature dependence slope above room temperature has been observed for the thermal conductivity of nand p-type 6H SiC. 55Figure 7 also shows the temperature dependence slope extracted from the calculated thermal conductivity.The deviation from the values determined from the experimental temperature dependence at high Si concentration can be attributed to the overestimation of the FE-scattering and the neglecting of the electronic thermal conductivity as discussed above.

V. CONCLUSIONS
The thermal conductivity of HVPE grown bulk GaN with varying Si doping concentration was measured by the 3ω method in the temperature interval of 295-470 K.It was found that the thermal conductivity decreases with both increasing temperature and increasing Si concentration.The experimental data was analyzed by use of a modified Callaway model which included all relevant phonon scattering mechanisms.From the fit of the experimental data for an undoped sample, Grüneisen parameters for acoustic phonons were extracted and then used for calculation of the thermal conductivity in the doped samples.A very good correlation between the simulations and the experiment was found for Si concentration up to 1×10 18 cm -3 .It was shown that for a proper explanation of the measured data the phonon-free-electron scattering was necessary to include in the model.At higher doping levels the model underestimated the thermal conductivity, especially at higher temperatures.The discrepancy was attributed to the electronic thermal conductivity which appears to contribute essentially to the total thermal conductivity in addition to the lattice thermal conductivity especially at higher doping and higher temperatures.The slope of the temperature dependence of the thermal conductivity gradually decreases with the increasing Si doping.It was found that at T > 350 K the thermal conductivity in the highest doped sample (7×10 18 cm -3 ) was higher than that in lower doped samples.The results from this study could be of particular importance for thermal management of devices grown on n-type GaN substrates and operating at high power and elevated temperatures.

FIG. 1 . 5 Paskov
FIG. 1. Measured drop voltage V 3ω as a function of frequency (symbols) together with a linear fit in a semi-logarithmic plot (solid lines) at four different temperatures for undoped GaN (sample A).

095302- 6 Paskov
et al.AIP Advances 7, 095302 (2017) N D is the total dislocation density, N S and N E are the densities of screw and edge dislocations, respectively, η is the weight factor accounting for the mutual orientation of the dislocation line and the temperature gradient, and b S = c o and b E = a o √ 2/3 are the magnitudes of Burgers vectors for the screw and edge dislocations (a o and c o are the lattice constants in a wurtzite crystal).

FIG. 2 .
FIG. 2. Thermal conductivity of undoped and Si doped HVPE grown GaN.The solid lines represent the best fit of the experimental data (a) and the calculated thermal conductivity with the Grunstain parameters extracted from the fit (b).
FIG. 4. Comparison of the different models for the thermal conductivity of the moderately doped sample (sample C).The best match with the experimental data is obtained for the model with the high temperature (HT) form of the N-scattering and the inclusion of the FE-scattering.

FIG. 5 .FIG. 6 .
FIG. 5. Thermal conductivity of bulk GaN as a function of Si concentration at different temperatures.The experimental values (symbols) and the theoretical dependences calculated with (solid lines) and without (dot-dashed lines) including the FE-scattering are shown.
FIG. 7. The slope of the temperature dependence of the thermal conductivity as a function of Si concentration.
25wer oscillating with frequency 2ω which is dissipated in into the sample beneath.The heat waves from the sample modulate the resistance of the metal wire due to the non-zero temperature coefficient of resistivity.As a result, a third harmonic component appears in the voltage drop across the wire.The amplitude of this component, V 3ω , is given by25

TABLE II .
Material parameters of GaN used in the thermal conductivity calculations.