Micro-cooler enhancements by barrier interface analysis

A novel gallium arsenide (GaAs) based micro-cooler design, previously analysed both experimentally and by an analytical Heat Transfer (HT) model, has been simulated using a self-consistent Ensemble Monte Carlo (EMC) model for a more in depth analysis of the thermionic cooling in the device. The best fit to the experimental data was found and was used in conjunction with the HT model to estimate the cooler-contact resistance. The cooling results from EMC indicated that the cooling power of the device is highly dependent on the charge distribution across the leading interface. Alteration of this charge distribution via interface extensions on the nanometre scale has shown to produce significant changes in cooler performance.


I. INTRODUCTION
Consumer demand for smaller, more powerful electronic devices is on the rise.The ethos of modern design is to produce devices half the size but double the power of their predecessors.However, each advance carries with it the problem of heat management.The cooling systems, vital to the operational efficiency and extended operational life of these devices, are being left behind.To catch up, these systems must undergo redesign through detailed computer aided analysis to enhance their capabilities.
A possible solution to the heating problem of micron sized devices is the thermionic cooler.These devices allow for integration into the same device chip, giving cooling without the bulky excess inherent to gas and liquid coolers.In thermionic emission, electrons transition over a barrier, transporting heat energy from one side of the device to the other, resulting in a temperature differential.Heat diffusion back across the barrier to the cooled side can potentially ruin this cooling effect hence careful consideration of the barrier structure needs to be taken into account when designing such devices to enhance the electronic transport while inhibiting heat diffusion.
This paper presents a self consistent Ensemble Monte Carlo (EMC) model, capable of simulating the heterobarrier structures such as those found within micro-cooler devices.The model has been previously validated, showing excellent agreement with both theoretical and experimental results from the literature. 1 The earlier analysis of those devices revealed an intricate interplay between the field, charge distribution and band structure which lead to proposed engineering modifications of the interface for enhancements to the thermionic cooling.In this paper, the model, discussed in Section III A, is applied to a unique GaAs based micro-cooler design.outlined in section II previously simulated by an analytical Heat Transfer (HT) model (described in Section III B) 2 matched to experimental results.Those results are built upon and analysed in more detail with the EMC model in Section IV.Furthermore, proposed enhancements to the design are presented which demonstrate, theoretically, that higher cooling powers are possible from the same device with minor modifications to the band structure of the interface.

II. MICRO-COOLER
The device consists of an aluminium gallium arsenide (AlGaAs) based superlattice (SL) sandwiched between contacts constructed from layers of GaAs and graded AlGaAs.A schematic of the cooler is given in Figure 1.From top to bottom, it consists of 300 nm of highly doped GaAs followed by 50 nm of linearly graded Al 0-0.2 Ga 1-0.8As making up the cathode region which is grown on top of a 2000 nm long AlGaAs based SL.The anode region is composed of linearly graded 50 nm Al 0.2-0 Ga 0.8-1 As layer, ending with 300 nm of GaAs.The SL is constructed from alternating layers of Al 0.2 Ga 0.8 As and Al 0.1 Ga 0.9 As, each 10 nm long, for 100 periods ending with Al 0.2 Ga 0.8 As.The resulting band discontinuity is shown in Figure 2. The entire device was grown on a 620 μm semi-insulating GaAs substrate.Contact areas were 5,000 μm 2 with a doping of 8 × 10 24 m −3 while the SL was doped to 2 × 10 24 m −3 .

A. EMC model description
The comparison of the HT model to EMC modelling is akin to macroscopic and microscopic analysis: where HT modelling is dependent on knowledge of the thermal transport parameters EMC modelling allows for those parameters to be innately generated through semi classical analysis of the device, while revealing important details of the operation of the cooler that are difficult to obtain from the analytical model alone.While the EMC model can obtain very detailed information about FIG. 2. Schematic of the device band structure (full line) against the simulated approximation used in the model (dotted) across the cathode region.
the active region of the device, it does not currently include the effects of the contacts.For the fullest analysis of the micro-cooler both EMC and HT models are applied.
The EMC model was written for the self consistent calculation of electron transport, charge distribution, electric field and heat generation in semiconductor devices.The model is semi-classical, using the super-particle approximation 3 to simulate electron movement in the real device.Electronphonon scatter events are accounted for in the model including, acoustic, polar-optical, polar and inter-valley scattering.Impurity scattering was scaled to achieve the correct mobility for the doping. 4he energy transfers between the electrons and phonons during the scatter events are used as the source term for heat generation, which is used as the input of the heat diffusion algorithm for calculate changes in device temperature.Electron position and energy is monitored by the model.Those approaching the barrier with sufficient energy will be allowed to surmount the barrier, undergoing a kinetic energy loss equivalent to the barrier height, while those without sufficient energy will be reflected at the barrier interface.Transitions to the barrier via inter-valley electron transfers are accounted for however these are in the minority as most electrons exist in the central valley prior to reaching the barrier interface.

B. HT model description
The HT model was developed by manipulating a variety of equations that have been previously published and generally accepted, 5,6 below: where T = T h − T c .Here, T j , T a , T, T c and T h , are the junction temperature, ambient temperature, temperature differential and temperatures of the cool and hot sides of the cooler respectively.c , d , and SUB are the thermal impedance of the cooler, device and substrate.R T , R tot and R SL are the top (cathode) contact resistance, total resistance of the cooler and resistance of the superlattice respectively.P d is the heat flux from the device being cooled while P c is the heat flux output from the cooler.S is the bulk Seebeck coefficient and I the current required to power the cooler.The combination of these equations have allowed for some unknown variables to be removed from equations to utilise the measurable parameters.This allows the model to be re-worked to find a large number of parameters based on the available input parameters used, for instance, if temperature and input powers are known, then thermal conductivities can be calculated.The model requires the electrical and thermal resistivities of the included materials be known.The Seebeck coefficient is normally taken from the literature however, the EMC modelling has enabled structure specific Seebeck coefficients to be computed.The HT model of the micro-cooler is set up to plot the absolute temperatures of both the hot and cool sides of the micro-cooler as well as the temperature difference between these two values as a function of the current applied to the micro-cooler.

C. Simulated structure
The quantum well structure of the SL is on the same scale of the electron wavelength and can not be accurately simulated using the semi-classical approach.Therefore an approximation to the structure was used which considers the SL to be a single block of pseudo AlGaAs material with Al molar fraction of 0.15.The resulting change to the conduction band when this approximation is implemented is shown schematically in Figure 2. The material parameters were taken as given in 1 and the SL parameters were assumed to be the same as AlGaAs, but with a lower bulk thermal conductivity.As there is no recorded thermal conductivity value for this particular SL and as previous publications on other SL materials reported values just over half that of the bulk, 7 an assumed thermal conductivity of 0.06 Wcm −1 K −1 was used.

D. Simulation procedure
Each simulation was preformed using an ambient temperature of 401 Kelvin to be consistent with the experimental measurements.The cathode temperature boundary condition was allowed to move freely, as calculated by the model, while the anode temperature was maintained at the ambient temperature.Steps were taken to reduce the electric field noise, inherent in simulations of this type, although due to computational restrictions a balance between accuracy and computational time had to be reached.As such the cathode side GaAs layer was subject to some ohmic heating, an artefact of the model due to the field noise in the region.This heat was found to be generated at the rate equivalent to that expected from a 0.29 resistance contact (i.e.∼12 mW at 0.2 A).The significance of this heat generation is covered in Section IV A. A very large number of super particles were used (in the region of 10 5 particles) alongside rapid calculation of the field changes using a time step of 0.1 fs with a spatial grid spacing of 1 nm.To achieve a high level of accuracy in the determination of the device temperature the scatter arrays were set to a resolution of 0.5 meV for energy and a sensitivity of 0.075 K for device temperature.The ensemble of particles was allowed to settle for a period of 30 ps prior to any recordings being made to prevent any errant heating from the initial state approximations during initialisation.The simulation was performed over a range of voltages from 8 to 25 mV, setting the maximum current to approximately 1 A to reproduce the range used in the experimental measurements.

E. Temperature determination
Each temperature profile was allowed to develop until a steady state was reached in which the device temperature no longer changed significantly over a period of several nanoseconds.As electron transport and heat diffusion operate on significantly different time scales, a full simulation of the heat development using EMC would be extremely time consuming, therefore, a procedure in which the model alternated between the coupled EMC simulation mode and the analytical Heat Diffusion (HD) mode was utilised.During this procedure EMC simulation is performed as described in section III A where the changes in lattice energy were self consistently calculated from the electronphonon interactions.Once a sufficient amount of time has passed (in the region of tens of picoseconds) the HD mode is activated, here the net lattice energy recorded over the EMC simulation period is used as the heat generation term during the HD period.The temperature changes in the device are then calculated analytically for a large period of time (1000-1500 ps).At the end of the HD mode the lattice energy is updated to account for the HD run period and the model switches back to EMC simulation.This process is repeated until the steady state device temperature is reached.This 'leap frog' method allows the simulation of temperature growth over a larger time scale than would be possible utilising EMC alone.The results follow 8 where the temperature change recorded at the location where the interface starts is taken as the cooling.

A. Fitting to experimental data
Preliminary simulation of the cooler using the approximated band structure, detailed in section III B, resulted in excessive heating within the cathode side region.Analysis of the data indicated a heightened field level localised across the Al 0-0.15 Ga 1-0.85As interface, indicating the charge distribution was the cause of the heating in the region.Results of the time averaged potential from the simulation shows a clear shift in potential across the interface, compelling the electrons in the region to surmount the barrier through energy excitation due to the enhanced field.The localised excitation of the electrons at the interface region reduces the effective barrier height, lowering the energy requirement for barrier transition.The effect is demonstrated in Figure 3 where the barrier height, determined by EMC simulation (solid line), has been reduced by the field at the interface from the input band gap energy (dotted line).The data from this particular simulation demonstrates approximately 50% decrease in barrier height.The consequence of the lowered barrier height is the diminished removal of heat energy via thermionic emission of electrons passing over the barrier, leaving the contact region subject to Joule heating.These results concur with the previously performed simulations on similar cooler devices 1 which suggested that the initial approximation to the band discontinuity is too abrupt, producing an unnatural redistribution of charge at the interface, leading to the field distortion which reduces the effective barrier height and hence lowers the amount of energy extracted from the region by each electron.AIP Advances 4, 027105 (2014) FIG. 4. Cooling against applied current flow from EMC when 57 nm interface is used (•) against experimental measurements for the device using the original 50 nm interface (+).
The hypothesis is that extending the Al 0-0.15 Ga 1-0.85As interface will produce a more gradual band discontinuity, allowing for a more natural and stable charge distribution to form across the region thus reducing the diminishing effect of the field on the barrier height.In order to accomplish this, a series EMC simulations were performed as before with each successive simulation extending the Al 0-0.15 Ga 1-0.85As interface from the initial 50 nm length, until cooling was observed.
The simulation, being an approximation to the micro-cooler as it contains only the active region with no contacts or substrate losses included, predicts cooling of an idealised device, the results of which will be higher than those obtained experimentally.However, the underlying principles and dynamics of the cooler are present within the simulation.Therefore, the qualitative results will show performance trends.Hence, comparison between experimental and model results should focus on the qualitative similarities rather than the quantitative.
The smallest length of interface which showed cooling was 57 nm, below this the region was subject to heating of a few degrees.Cooling obtained for this length, against applied current, is shown in Figure 4.The cooling is localised on the cathode side of the barrier interface, with the rest of the contact undergoing small amounts of ohmic heating due to the erratic field in this region.This is a very delicate state for the contact region, where the device is sitting on the knife edge between heating and cooling that exists when the interface length increases from 56 to 57 nm.This highlights the sensitivity of the overall cooling effect to the very small changes in interface lengths, which could well be within material growth errors and variation across a wafer.
The cooling results when the 57 nm length was used show a maximum cooling of 0.82 K at an applied current of 0.6 A, roughly coinciding with the peak position in the experimental data at ∼0.5 A. Below 0.3 A the model showed no cooling, however, due to the very small amounts of heating it is hypothesised that the cooling results may have been lost to noise in the field due to the susceptibility of lower currents to the temperature fluctuations.Subsequent EMC simulations were performed which indicated that additional increases in the Al 0-0.15 Ga 1-0.85As interface length could generate further gains in cooling power from the device.To find an improved match to the experimental data, the interface length was increased from 57 to 60 nm. Figure 5  resemblance to the experimental data due to the considerable noise in the cooling result, likely due to the susceptibility of the current to small changes in device temperature.The cooling obtained from using the 60 nm interface shows a close resemblance in pattern to the experimental data, however with higher cooling values.However, this is expected as the idealised cooler has no lossy contact or substrate resistance as already discussed.
The experimentally measured resistance for the real device, including both contacts is approximately 1.3 , with the bottom contact resistance estimated at 1.1 as the calculated values for the top contact resistance (from the specific contact resistance found from TLM analysis) and the bulk of the micro-cooler (from mobility and carrier concentration modelling) came to a combined value of 0.2 .
Using the EMC generated results for conductivity, the resistance of the active region of the micro-cooler was determined to be 0.093 , with 0.01 of that across the cathode side GaAs layer.Taking the experimental measurement of the combined resistance of the top contact and cooler, with the calculated active region resistance above, gives an estimated 0.11 for the top contact.
For a more accurate estimation of the top contact resistance, the analytical HT model was modified to include the 60 nm interface and micro-cooler resistance of 0.093 .Simulations were performed, calculating the resulting temperature changes for a range of different cathode resistances, R c , as seen in Figure 6.The variables stated within the model were adjusted to match the values from the EMC model.Comparing the cooling results from EMC (Figure 5) to those generated by the HT model (Figure 6) for the estimated R c of 0.11 shows that the HT model predicts a higher amount of cooling than seen in the EMC results.The closest match in cooling reported by both models occurs between R c = 200 and R c = 240 m , placing the combined (top-contact and micro-cooler) resistance around 293 to 333 m .
As mentioned previously, the EMC results show some heat generation within the cathode side GaAs layer, caused by the artificial noise in the field.This heating is equivalent to 0.3 of resistance, which is the combination of the nominal cooler resistance of 10 m (in the cathode side GaAs layer) and an effective contact resistance of 0.29 produced by the field noise.This result is higher than the initial estimate of R c (0.11 ), but explains why the best match between both models occurs at the higher range of R c .FIG. 6. Cooling results from analytical HT model for variable cathode contact resistances for the micro-cooler using the 60 nm length interface and active region resistance of 0.093 .
To obtain a value for the effective Seebeck coefficient, the analytical function describing the resulting temperature differential of the device, T, via heat transfer was used.The expression is: 6 where R is the contact resistance.As the model does not incorporate an active device attached to the cooler, the term P d is set to zero.Arranging the expression in the form of y = mx + c yields, Plotting the values of I/ T and I 2 / T using the values obtained from the EMC simulation, the slope of the plot gives the value of ST c /R, allowing S to be calculated for each applied current and resulting T c .Using this method the Seebeck coefficient was calculated to be 271, 281 and 289 μV/K for the 57, 59 and 60 nm interface lengths respectively.For reference, an assumed S of 256 μV/K was used when matching the analytical model to the experimental data. 2 To reach this value, the resistance of the device would need to be reduced to approximately 0.08 , the difference being well within the error of the model.

B. Further enhancements
Extension of the Al 0-0.15 Ga 1-0.85As interface from 57 to 60 nm resulted in an improvement of cooling of around 1.8 K.It was found, that even greater cooling could be obtained through the use of longer interface regions, without modification to the device volume.
Simulations were performed for the same current range as before, this time using Al 0-0.15 Ga 1-0.85As interface lengths of 70 and 80 nm.The data from these simulations are shown in Figure 7 alongside that of the cooling obtained for the 60 nm interface for comparison.While the 57 nm length produced cooling of under 1 K, extending to 60 nm resulted in cooling of 2.5-4 K being seen.The use of a 70 nm length interface (a total extension of 20 nm) produced an increase in cooling of approximately 3 times that seen for 60 nm, while increasing the length from 70 to 80 nm (total extension of 30 nm) results in the cooling doubling.The drastic changes in cooling occurring for such small changes in interface length, emphasizes the sensitivity of the thermionic effect to the  dynamics of the barrier interface.The high cooling results presented here may not be measurable experimentally, however they can be matched by the HT model under idealised conditions, such as when low contact resistance is used.
With the use of the 70 and 80 nm interfaces, the field energy contribution to the electron barrier transition is greatly reduced, allowing the cooling at the interface to become quite substantial, lowering the temperature across the entire cathode contact region.In doing so, the electron distribution in that region is shifted into the lower energy regime, and with the much reduced field across the interface, the thermionic effect is significantly enhanced.
To obtain the values of S for these simulations, a slightly different approach was used.As mentioned in the description of the model, the energy exchanges between the material lattice and electrons due to electron-phonon interactions are recorded with position, as the net lattice energy.This term allows for the calculation of the heat loading, Q TI , on the cathode side of the SL from the summation of the net energy exchanges occurring within the cathode side region, where the E e-ph (x) term is the net energy recorded at position x over time period, t, in device of cross sectional area A. This data can then be used in conjunction with the heat balance equation, to recover the S value for each simulation. 8Here T is the ambient temperature of the cooler, V the voltage drop across the SL, j the current density, k the thermal conductance and d the length of the SL.Each term in the expression represents a different aspect of the heat movement within the device.From left to right, the STj term represents the heat power removed due to thermionic emission, Vj is the ohmic heating in the barrier where α represents the portion of this heating returning to the cool side and k d t is the back heat flow due to the induced temperature gradient.To obtain the value of S for each interface length, Q TI + k/d T was plotted against j.The slope of this plot gives the value of the "ST-αV" term, where α was assumed to be unity for the worst case scenario where the full ohmic back heating is experienced.Assuming that the majority of the applied voltage lies along the SL barrier, S can be recovered for each voltage.The Averaged S values from these calculations are presented in Table I.The Seebeck coefficient calculated for the 60 nm interface using this method is ∼16% greater than that obtained using the macroscopic heat transfer function, described in Section IV A, possibly due to the recorded Q TI in the model being a more accurate measurement of heat lost.This data indicates that extension of the interface would lead to improvements in the Seebeck value.

V. CONCLUSION
The simulation of the GaAs based micro-cooler by the EMC model has shown that the device is very sensitive to small changes in the cathode to superlattice interface length.These small changes alter the charge distribution across the interface, which modifies the conduction band structure resulting in significant changes to the cooling power of the device.
The smallest interface length required to show cooling was found to be 57 nm, 7 nm longer than the initial approximation.Large enhancements of the cooling were possible through extension of the interface by only tens of nanometres.Such differences in layer width are within the wafer growth uncertainty and may account for the variation of cooler performances seen across a wafer.
The minimum interface length required to to reproduce the cooler performance trend was 60 nm.Through comparison of the HT to EMC simulation results, the cathode resistance was estimated in the region of 0.2 to 0.24 , for a contact and cooler combined resistance of 0.293 to 0.333 .Furthermore, it was found that heat was being artificially generated in the cathode side GaAs layer due to field noise, equivalent to that produced by a 0.29 contact.
A bulk Seebeck coefficient of 271 μV/K was found for the smallest non-heating interface length by EMC, in good agreement with the Seebeck coefficient used in the HT model analysis when fitting to the experimental data.Extending the interface length by only 20 nm lead to improvements in the Seebeck coefficient from 346-1227 μV/K.It is improbable that such high values could be obtained experimentally as the simulation does not account for the effect of the contact resistance which can have significant effects on the overall cooling performance of the device.However, the qualitative result still stands, that small modifications of the interface length allows for improvement of S and hence the overall cooling mechanism.
FIG. 1. Schematic of the micro-cooler showing the layer arrangement surrounding the AlGaAs based superlattice.The regions not included in the EMC model are shaded grey.

FIG. 3 .
FIG.3.Comparison of conduction band from simulations of 50 nm AlGaAs interface (solid line) to the nominal band gap energy (dotted line), extracted from data for a bias of 24 mV.The conduction band was calculated from the combined effect of the potential and nominal band energies.

TABLE I .
Calculated average Seebeck values for the AlGaAs interface lengths.