# Thermal Modeling of the GaN-based Gunn Diode at Terahertz Frequencies
## 2 Method of Analysis and Physical Models  
Compared with the drift-diffusion (DD) model, the EB model is more accurate for spatiotemporal domain and chaos in the short channel. This is because the EB model includes non-localized carrier transport phenomena, such as the velocity overshoot and the non-local impact ionization which are neglected by the conventional DD model. Therefore, the EB model is more suitable for the transient simulation of terahertz devices and is employed in this paper. Indeed, the non-local transport phenomena can be modeled by using Monte Carlo method or using higher order moments of the Boltzmann transport equation, though huge computation times are required by these approaches. In this paper, higher order solution to the general Boltzman transport equation is solved based on the EB model at the Atlas simulation platform, which consists of an additional coupling of the current density to the carrier' temperature, or energy [18]. MODELS statement is the important part in the device-simulation code, where the physical models of the device is specified. In the simulation code of Gunn diode in Atlas, the energy balance model is selected by specifying the 'HCTE' command in the MODELS statement to enable the energy transport model. The energy relaxation time \(\tau_{\epsilon}\) is defined as 150 fs based on [4; 19]. Meanwhile, the self-heating simulations are incorporated in the GaN Gunn diode structures using GIGA module of Atlas to account for lattice heat flow and general thermal environments. GIGA module accounts for the dependence of material and transport parameters on the lattice temperature. For simulation of self-heating effects, it adds the heat flow equation to the primary equations that are solved by Atlas. The heat flow equation has the following expression:\[\mathrm{C}\frac{\partial T_{L}}{\partial t}=\
abla(\kappa\
abla T_{L})+H \tag{1}\]  
where C is the heat capacitance per unit volume, \(\kappa\) is the thermal conductivity, H is the heat generation, and T\({}_{L}\) is the local lattice temperature. In our simulation, the thermal conductivity is set as a temperature-dependent thermal conductivity: \(\kappa=\kappa_{300}/(\Gamma/300)^{\alpha}\), where \(\kappa_{300}\) is the thermal conductivity at the temperature of 300 K and \(\alpha\) is the temperature dependence coefficient [20]. For GaN, \(\kappa_{300}\) is set equal to 1.6 W/cm\(\cdot\)K and \(\alpha=1.4\)[20]. The heat capacitances C are 3.0135 J/cm\({}^{3}\)/K for GaN. All the interfaces between the diode and external surroundings are set to be adiabatic boundary conditions. In order to simplify the algorithm and solve the non-convergence problem, the heat transport from the metal contact to external surroundings is not included in our simulations, a common practice in the simulation of thermal effect as described in [21; 22; 23]. The LAT.TEMP lattice temperature model is added to the MODELS statement to include heat flow. We have used a non-isothermal energy balance (NEB) model by specifying both HCTE. and LAT.TEMP parameters in the MODELS statement to solve both energy balance equations and heat flow equations. The NEB model is an extension of stratton's energy balance model for the case of nonuniform lattice temperature, which is a set of six partial differential equations for electrostatic potential, electron and hole concentrations, electron and hole carrier tempertures, and lattice temperature [18]. For the impact ionization model, we adopt Toyabe impact ionization model, which is available when the energy balance transport model is applied [18; 24]. A key parameter of the model is the energy balance length L\({}_{\mathrm{REL}}\), which can be calculated by L\({}_{\mathrm{REL}}\) = \(\mathrm{v}_{\mathrm{sat}}\times\tau_{\varepsilon}\), where \(\mathrm{v}_{\mathrm{sat}}\) are the saturation velocities for electrons, and \(\tau_{\varepsilon}\) corresponds to the energy relaxation time.