Towards Monolithic Indium Phosphide (InP)-Based Electronic Photonic Technologies for beyond 5G Communication Systems

: This review paper reports the prerequisites of a monolithic integrated terahertz (THz) technology capable of meeting the network capacity requirements of beyond-5G wireless communications system (WCS). Keeping in mind that the terahertz signal generation for the beyond-5G networks relies on the technology power loss management, we propose a single computationally efﬁcient software design tool featuring cutting-edge optical devices and high speed III–V electronics for the design of optoelectronic integrated circuits (OEICs) monolithically integrated on a single Indium-Phosphide (InP) die. Through the implementation of accurate and SPICE (Simulation Program with Integrated Circuit Emphasis)-compatible compact models of uni-traveling carrier photodiodes (UTC-PDs) and InP double heterojunction bipolar transistors (DHBTs), we demonstrated that the next generation of THz technologies for beyond-5G networks requires (i) a multi-physical understanding of their operation described through electrical, photonic and thermal equations, (ii) dedicated test structures for characterization in the frequency range higher than 110 GHz, (iii) a dedicated parameter extraction procedure, along with (iv) a circuit reliability assessment methodology. Developed on the research and development activities achieved in the past two decades, we detailed each part of the multiphysics design optimization approach while ensuring technology power loss management through a holistic procedure compatible with existing software tools and design ﬂow for the timely and cost-effective achievement of THz OEICs. The compact model for the UTC-PDs reported in [13] is a Verilog-A model for describ-ing the semiconductor physics governing the junction currents and the optical response of the photodiode. The model is further enhanced by incorporating improved descriptions of the dark and photocurrents. The schematic (Figure 3) shows the equivalent circuit of the UTC-PD, including extrinsic pads and coplanar waveguide (CPW) parasitic elements. The intrinsic part is represented by the diode dark current source, I dark , the junction capacitance, C j , the shunt resistance, R sh , the series resistances, R s , and the photocurrent, I ph ( w ) . In order to facilitate model scaling, the diode currents and the junction capacitances are normalized by the device active area, S , of UTC-PDs. of the different characteristics the representative device ) at different forced V BE (high-injection) and I . In all cases, good agreement between the experimental data and model simulation is observed, indicating the accuracy of the temperature-related parameters, high-current and self-heating effects. for 25 and 50 µ m-long photodiode. Good accuracy of the photocurrent model, depicted in (4), has been observed. Moreover, good photocurrent scalability is also achieved. on Figure 10 illustrates the results comparing the magnitudes and phases of the measured low-frequency Y -parameters (symbols), Y 12 and Y 22 , and the scalable thermal model simulation (solid lines) for different device geometries. The results are shown for a bias condition of V BE = 0.85 V and V CE = 1.2 V, chosen speciﬁcally around peak f T of the transistors, where self-heating becomes dominant. Both Y -parameters show good agreement between the measurement and the simulation results. Furthermore, good model scalability has been observed across all geometries.


Introduction
Signals in the frequency range of 30-300 GHz are best suited for ultra-broadband wireless communications. Unfortunately, most of the solid-state sources produce insufficient power (limited to few milliwatts) at the upper end of this spectrum to support high data-rate communication over useful distances (>100 m). To respond to this challenge, novel baseband processing technology, identified as the enabling technology for the future beyond 5G networking, is capable of focusing energy into finer-grained areas and provide significant increase in both transmission energy efficiency and spectrum efficiency. As a consequence, accurate compact models and the associated parameter values for the UTC-PDs, the DHBTs as well as the passive elements should be jointly elaborated in order to guide the development of the OEIC at the early stage of its design. These compact models should be included in a single process design kit (PDK) that could be incorporated into the OEIC design tools. Featuring the ensemble of the building blocks, the PDK could thus allow the design and simulation of UTC-TIA-PA.
The aim of this article was to provide a comprehensive evaluation of technological advances both in terms of sheer technological performances (photonic, electrical and thermal) as well as the software tools available for accomplishing this innovative design. In terms of design tools, in particular, electronic and photonic communities have very different simulation tools available, owing to the current status of their respective state of the arts. Integrated circuit design would thus necessitate multiphysics tools for facilitating the development of this emerging technology. Despite the differences in the design methods of the two communities, electronic design automation (EDA) tools have evolved and matured through years of industrial endeavors, leading to an abundance of SPICE compatible compact models and circuit design environments. While the photonic community is still lacking a key circuit design environment, the best possible solution would be to ensure the compatibility of photonic device simulation with existing EDA environments. Therefore, in the current work, we propose not only an up-to-date DHBT compact model dedicated to devices on InP substrate but we also developed and enhanced a compact model specifically for UTC photodiodes with the goal of making way for a common design framework. Both compact models are fully compatible with the EDA circuit libraries. Moreover, to ensure that common key physical quantities be handled, DHBT and UTC-PD compact models compute either electronic or photonic quantities along with their thermal parameters, which can, in the long run, pave the way for multiphysics tools.
In order to pursue the aforementioned ideas further, in this paper, we present efficient compact models that are compatible with UTC photodiode and DHBT technologies developed by III-V Lab, which may be considered as the first step towards the efficient monolithic co-integration of electronic and photonic parts for this unique OEIC technology of the future. The rest of the paper is organized as follows. Section 2 is dedicated to the detailed technologies under investigation, UTC-PD and InP DHBT technologies, both from III-V Lab. Section 3 presents the electrical compact models associated with the two devices. Section 4 shows the DC characterization of UTC-PDs and InP DHBTs for the purpose of compact model validation. Section 5 provides insights into the thermal characterization and associated compact model sub-circuits that allow the management of selfheating as well as the mutual thermal coupling. Section 6 deals with the high frequency The aim of this article was to provide a comprehensive evaluation of technological advances both in terms of sheer technological performances (photonic, electrical and thermal) as well as the software tools available for accomplishing this innovative design. In terms of design tools, in particular, electronic and photonic communities have very different simulation tools available, owing to the current status of their respective state of the arts. Integrated circuit design would thus necessitate multiphysics tools for facilitating the development of this emerging technology. Despite the differences in the design methods of the two communities, electronic design automation (EDA) tools have evolved and matured through years of industrial endeavors, leading to an abundance of SPICE compatible compact models and circuit design environments. While the photonic community is still lacking a key circuit design environment, the best possible solution would be to ensure the compatibility of photonic device simulation with existing EDA environments. Therefore, in the current work, we propose not only an up-to-date DHBT compact model dedicated to devices on InP substrate but we also developed and enhanced a compact model specifically for UTC photodiodes with the goal of making way for a common design framework. Both compact models are fully compatible with the EDA circuit libraries. Moreover, to ensure that common key physical quantities be handled, DHBT and UTC-PD compact models compute either electronic or photonic quantities along with their thermal parameters, which can, in the long run, pave the way for multiphysics tools.
In order to pursue the aforementioned ideas further, in this paper, we present efficient compact models that are compatible with UTC photodiode and DHBT technologies developed by III-V Lab, which may be considered as the first step towards the efficient monolithic co-integration of electronic and photonic parts for this unique OEIC technology of the future. The rest of the paper is organized as follows. Section 2 is dedicated to the detailed technologies under investigation, UTC-PD and InP DHBT technologies, both from III-V Lab. Section 3 presents the electrical compact models associated with the two devices. Section 4 shows the DC characterization of UTC-PDs and InP DHBTs for the purpose of compact model validation. Section 5 provides insights into the thermal characterization and associated compact model sub-circuits that allow the management of self-heating as well as the mutual thermal coupling. Section 6 deals with the high frequency characterization and modeling of both devices. Section 7 discusses additional effects to be taken into consideration through the technology computer-aided design (TCAD) simulation, the evaluation of key performance metrics though model extrapolation in comparison with the state of the art as well as the methodology for assessing OEIC circuit reliability. Finally, Section 8 provides the conclusion and future perspectives.

Technologies
The III-V Lab develops both UTC-PD and InP DHBT technologies within its clean rooms. UTC-PDs and DHBTs share the same base-collector interface heterostructure and therefore, for this interface, the physical models are similar despite differences in thicknesses and doping levels between optimized UTC-PD structures and optimized DHBT structures. This technological similarity opens up the way to integrate DHBT and UTC photodiodes on the same substrate, in order to suppress the parasitic on the interface between them (Figure 1). This could be done either by finding a shared epitaxial structure with the best compromise between UTC-PD and DHBT structures or by a selective regrowth of UTC-PD layers on a DHBT layer or vice versa.
The InP DHBT structure is grown externally by solid source molecular beam epitaxy on 3 wafers. The vertical structure is optimized in order to achieve balanced frequency performances for 0.7 µm emitter DHBTs with f T and f max around 400 GHz. The structure is similar to the one described in [9] except for the low doped collector region for which Si-doping level is reduced to 1.5 × 10 16 at/cm −3 to improve the breakdown voltage. The 0.5 and 0.7 µm InP DHBTs with emitter length varying from 5 to 10 µm are processed using a wet-etched triple mesa technology. Figure 2a  tion, the evaluation of key performance metrics though model extrapolation in comparison with the state of the art as well as the methodology for assessing OEIC circuit reliability. Finally, Section 8 provides the conclusion and future perspectives.

Technologies
The III-V Lab develops both UTC-PD and InP DHBT technologies within its clean rooms. UTC-PDs and DHBTs share the same base-collector interface heterostructure and therefore, for this interface, the physical models are similar despite differences in thicknesses and doping levels between optimized UTC-PD structures and optimized DHBT structures. This technological similarity opens up the way to integrate DHBT and UTC photodiodes on the same substrate, in order to suppress the parasitic on the interface between them (Figure 1). This could be done either by finding a shared epitaxial structure with the best compromise between UTC-PD and DHBT structures or by a selective regrowth of UTC-PD layers on a DHBT layer or vice versa.
The InP DHBT structure is grown externally by solid source molecular beam epitaxy on 3 wafers. The vertical structure is optimized in order to achieve balanced frequency performances for 0.7 µ m emitter DHBTs with fT and fmax around 400 GHz. The structure is similar to the one described in [9] except for the low doped collector region for which Sidoping level is reduced to 1.5 × 10 16 at/cm −3 to improve the breakdown voltage. The 0.5 and 0.7 µ m InP DHBTs with emitter length varying from 5 to 10 µm are processed using a wet-etched triple mesa technology. Figure 2a  Emitter and base contacts as well as the base plug needed for interconnection are defined by e-beam lithography whereas stepper lithography is used for other technological steps. TiPdAu (titanium-palladium-gold), deposited by vacuum evaporation, is used for all contacts. After encapsulation with silicon nitride (SiN) for higher robustness and planarization with polyimide, emitter, base and collector contacts are opened by etching to form device interconnects. Figure 2b presents a scanning electron microphotograph of a 0.7 × 5 µm 2 InP DHBT.
The UTC photodiode structure is grown either by gas-source molecular beam epitaxy (GS-MBE) [10] or by metal-organic vapor-phase epitaxy (MOVPE) as in this paper. A multimode diluted waveguide is used for optimal fiber coupling [11]. The process begins with the etching of the top InP layers, followed by base and collector etching using a mixed of dry and wet etching. Pt is used for a P contact followed by a rapid thermal annealing and TiPtAu are used for the N contacts and electrodes (Figure 2c). Figure 2d shows a 5 × 25 µm 2 InP UTC photodiode with the metal interconnections.
To evaluate the performances of this integrated co-design approach with the UTC-PDs and the DHBTs, accurate and physics-based compact models are required. This is even more important if the same epitaxy is shared, in order to obtain a precise overview of the impact on the performances of individual UTC-PDs and DHBTs as well as the potential gain offered by the co-design on the circuit performances. To that end, in this paper, we present two compact models compatible with UTC-PD and DHBT structures which are first validated on individual and optimized UTC-PD and DHBT chips.

Compact Modelling
Along with the progress in electronic design automation, compact models have continued to serve the electronics community at the core of all industrial integrated circuit design flow. Owing to the analytic forms of its equations, a compact model is significantly efficient in terms of computational burden and practical from a design point of view. Moreover, the correct formulation of the analytic equations can capture device physics with sufficient accuracy and is thus preferable to numerical device simulations by designers in order to speed up the simulation time. For electronic devices, the compact model council has certified a number of such accurate physics-based compact models compatible with SPICE simulation environments, including the High Current Model (HiCuM) model that has been widely accepted by the bipolar transistor community. As justified previously, the photonic community has been missing the equivalent tools in the context of optical devices and therefore, multiphysics formulation is crucial for the successful integration of electronic and photonic devices. With that in mind, as a first step, we present a modeling framework for both devices through an assembly of SPICE compact models for InP bipolar transistor (HiCuM [12]) and UTC photodiodes (in-house [13]). In the following section, we describe the fundamentals of these two models.

DHBT Compact Modeling
The HiCuM model [12] has focused on improving the dynamic description of bipolar transistors by putting the main emphasis on accurately modeling charge storage. This approach/representation can successfully capture both the static and dynamic behavior of the transistor through a set of analytical equations governed by the physics of charge build-up/partitioning within the bipolar transistor. Circuit design and optimization is based on proper device sizing and thus requires a geometry scalable compact transistor model such as HiCuM which is quite scaling-friendly. Except for the thermal and substrate coupling network, a compact analytical description as a function of device configuration exists for each element of the HiCuM equivalent circuit. Additionally, a new scalable model for the thermal impedance of III-V DHBTs, that has been developed based on the physics of heat diffusion within the HBT architecture [14], has been used to replace the existing simplistic thermal network of the HiCuM equivalent circuit for better accuracy (further elaborated in Section 5).

UTC Photodiode Compact Modeling
The compact model for the UTC-PDs reported in [13] is a Verilog-A model for describing the semiconductor physics governing the junction currents and the optical response of the photodiode. The model is further enhanced by incorporating improved descriptions of the dark and photocurrents. The schematic (Figure 3) shows the equivalent circuit of the UTC-PD, including extrinsic pads and coplanar waveguide (CPW) parasitic elements. The intrinsic part is represented by the diode dark current source, I dark , the junction capacitance, C j , the shunt resistance, R sh , the series resistances, R s , and the photocurrent, I ph (w). In order to facilitate model scaling, the diode currents and the junction capacitances are normalized by the device active area, S, of UTC-PDs. the physics of heat diffusion within the HBT architecture [14], has been used t the existing simplistic thermal network of the HiCuM equivalent circuit for bet racy (further elaborated in Section 5).

UTC Photodiode Compact Modeling
The compact model for the UTC-PDs reported in [13] is a Verilog-A mode scribing the semiconductor physics governing the junction currents and the o sponse of the photodiode. The model is further enhanced by incorporating impr scriptions of the dark and photocurrents. The schematic (Figure 3) shows the eq circuit of the UTC-PD, including extrinsic pads and coplanar waveguide (CPW) elements. The intrinsic part is represented by the diode dark current source, junction capacitance, , the shunt resistance, ℎ , the series resistances, , and tocurrent, ℎ ( ). In order to facilitate model scaling, the diode currents and the capacitances are normalized by the device active area, , of UTC-PDs.

Dark Current
In forward bias conditions, the dark current of the UTC-PDs is modeled by sical diode reverse saturation current, IS [13]. In reverse bias, in addition to the s current and the breakdown current modeled in [13], two other contributions h introduced to improve model accuracy: the band-to-band tunneling current, and the trap-assisted tunneling current, (in A). The tunneling effect is requir scribe the dominant current contribution under high reverse bias (V < −1) conditio general case, the tunneling current takes the following form [15]: .
and (in V/m) are parameters to be determined fro optimization, Vd is the diode voltage in V, S is the diode area in m 2 .
is the m electric field in V/m at the collector input. It can be expressed as [15] = 0 ( ) .
where 0 (in F) is the zero bias junction capacitance, is the collector's permittivity, is the junction potential (in V) and M is the grading coefficie junction. The trap-assisted tunnel current is dominant at low reverse bias (V > Typically, it is a Schockley-Read-Hall (SRH) generation current enhanced by t effects. It can be expressed as [16] = .
. E . V . − where (in AV −2 m −1 ) and (V/m) are trap assisted tunneling current par These parameters can be extracted from model optimization at low reverse bias.

Dark Current
In forward bias conditions, the dark current of the UTC-PDs is modeled by the classical diode reverse saturation current, I S [13]. In reverse bias, in addition to the saturation current and the breakdown current modeled in [13], two other contributions have been introduced to improve model accuracy: the band-to-band tunneling current, I BTB (in A), and the trap-assisted tunneling current, I TAT (in A). The tunneling effect is required to describe the dominant current contribution under high reverse bias (V < −1) condition. In the general case, the tunneling current takes the following form [15]: where A BTB (in A/V 3 ) and B BTB (in V/m) are parameters to be determined from model optimization, V d is the diode voltage in V, S is the diode area in m 2 . E max is the maximum electric field in V/m at the collector input. It can be expressed as [15] where C j0 (in F) is the zero bias junction capacitance, ε InGaAsP is the collector's relative permittivity, V j is the junction potential (in V) and M is the grading coefficient of the junction. The trap-assisted tunnel current is dominant at low reverse bias (V > −1) [14]. Typically, it is a Schockley-Read-Hall (SRH) generation current enhanced by tunneling effects. It can be expressed as [16] where A TAT (in AV −2 m −1 ) and B TAT (V/m) are trap assisted tunneling current parameters. These parameters can be extracted from model optimization at low reverse bias. It should Appl. Sci. 2021, 11, 2393 7 of 31 be noted that the temperature dependence of the static behavior of the dark current is captured through J S (T), C j0 (T) and V j (T) and the junction temperature is recalculated dynamically via an external thermal network [13].

Photocurrent
Dynamic photocurrent I ph (ω) mainly depends on transit times across the absorber, τ a (in s), and the collector region, τ c (in s). Its amplitude decreases over the frequency range. The dynamic response of the photocurrent has been modified in the model to take the form [17]: where J ph (ω = 0) is the static photocurrent density (in A/m 2 ) and ω is the angular frequency (in rad/s). Note that due to the low power and current density used in our experimental results, the dependence of the absorption coefficient on temperature has not been taken into account for photocurrent calculation.

DC Characterization and Model Validation
In order to validate the accuracy of the compact models described in the preceding sections, we performed extensive characterization on the two technologies which was followed by model validation and parameter extraction at the device level.

HiCuM Parameter Extraction Methodology
InP DHBT DC characterization were performed on several device dimensions (0.5 and 0.7 µm emitter widths and emitter lengths of 5, 7 and 10 µm) and compared with the HiCuM compact model simulation. For geometry scalable parameter extraction, we employed a dedicated macro model containing the scaling laws for all HiCuM Level 2.0 (L2) parameters related to emitter dimensions. Before the extraction of the model parameters using HiCuM/L2, some technological parameter values are declared and/or pre-calculated such as, for instance, sheet resistances, zero-bias hole charge, area and perimeter dimensions as well as given design rules and process layout information [12].
Parameters related to low-injection conditions are first extracted [12] since at highcurrent regime, a strong correlation exists between model parameters. Hence, the next step is the junction capacitances extraction. To extract the junction capacitances, the "cold Sparameters" method has been used, in which the transistor is biased to have negligible DC current flow, rendering forward bias across each junction limited to typically less than 0.5 V. The measured S-parameters under these conditions are then converted into Y parameters from which, after proper de-embedding, the depletion capacitances can be derived. The scalable extraction of collector and base current parameters is performed at low injection to prevent model parameter correlation usually in play between high-injection parameters. Considering the collector current split into its internal region and the periphery of the emitter window dimensions: here A E0 is the effective emitter area (m 2 ) and P E0 is the effective perimeter (in m). I CA (A.m −2 ) and I CP (in A.m −1 ) are the contributions from the intrinsic region and periphery of the device, respectively. Once the low-injection parameters are extracted, model parameters related to the high current effects are considered. The starting point is the extraction of the transit time at V BC = 0 V using the already extracted low-injection parameters. The following step is to take into consideration the high current as well as self-heating effects when extracting the transit time. The first step of the transit time parameter extraction method consists in determining the transit frequency from S-parameter measurements at various collector current densities and various voltages V CE or V BC . The transit frequency f T is determined using the current gain bandwidth product by e.g., choosing a spot frequency in the −20 dB/decade roll-off region of the current gain. The second step is the determination of the low-current transit time τ f0 and the transit time increase ∆τ f at high current densities. Then, in the high-current region, HiCuM/L2 model parameters are extracted by optimizing the region where f T reaches the peak along with the Gummel plot high-V BE region. The extraction of thermal resistance is of particular importance [18] since it is strongly coupled with other model parameters in the high-current regime. To obtain a good estimation of the self-heating effects, we used the intersection method [14] to extract the values of R TH from measurements. Additionally, a more detailed thermal impedance model has been adopted to take into account dynamic self-heating, which is detailed in Section 5. In HiCuM/L2, the temperature dependence of parameters for the collector current and emitter back-injection component of the base current is described through the associated saturation currents [12], while the temperature dependence of the zero bias hole charge is fixed by the previously determined capacitance parameters. In order to extract the model parameters, collector and base current measurements with respect to temperature are employed. Gummel plots at V BC = 0 V are well suited for this purpose. By applying the same extraction method for measurements at other temperatures, temperature coefficients of the remaining parameters (parasitic resistances, thermal resistance transit time, etc.) can be determined. These temperature-related parameters can play a significant role for advanced technologies for which the self-heating effect becomes increasingly important.  Figure 4c shows the scaling of I C as a function of the emitter length (for the 0.5 µm emitter width device) depicting the good accuracy of the scaling laws implemented through Equation (5). Figure 5a shows the temperature dependence of the collector current, at four different temperatures (−40, 0, 25 and 35 • C). Figure 5b,c depict the output characteristics for the representative device (0.7 × 5 µm 2 ) at different forced V BE (high-injection) and I B . In all cases, good agreement between the experimental data and model simulation is observed, indicating the accuracy of the temperature-related parameters, high-current and self-heating effects. transit time. The first step of the transit time parameter extraction method consists in determining the transit frequency from S-parameter measurements at various collector current densities and various voltages VCE or VBC. The transit frequency fT is determined using the current gain bandwidth product by e.g., choosing a spot frequency in the −20 dB/decade roll-off region of the current gain. The second step is the determination of the lowcurrent transit time τf0 and the transit time increase Δτf at high current densities. Then, in the high-current region, HiCuM/L2 model parameters are extracted by optimizing the region where fT reaches the peak along with the Gummel plot high-VBE region. The extraction of thermal resistance is of particular importance [18] since it is strongly coupled with other model parameters in the high-current regime. To obtain a good estimation of the self-heating effects, we used the intersection method [14] to extract the values of RTH from measurements. Additionally, a more detailed thermal impedance model has been adopted to take into account dynamic self-heating, which is detailed in Section 5. In HiCuM/L2, the temperature dependence of parameters for the collector current and emitter back-injection component of the base current is described through the associated saturation currents [12], while the temperature dependence of the zero bias hole charge is fixed by the previously determined capacitance parameters. In order to extract the model parameters, collector and base current measurements with respect to temperature are employed. Gummel plots at VBC = 0 V are well suited for this purpose. By applying the same extraction method for measurements at other temperatures, temperature coefficients of the remaining parameters (parasitic resistances, thermal resistance transit time, etc.) can be determined. These temperature-related parameters can play a significant role for advanced technologies for which the self-heating effect becomes increasingly important.  Figure 4c shows the scaling of IC as a function of the emitter length (for the 0.5 µ m emitter width device) depicting the good accuracy of the scaling laws implemented through Equation (5). Figure 5a shows the temperature dependence of the collector current, at four different temperatures (−40, 0, 25 and 35 °C). Figure 5b,c depict the output characteristics for the representative device (0.7 × 5 µ m 2 ) at different forced VBE (high-injection) and IB. In all cases, good agreement between the experimental data and model simulation is observed, indicating the accuracy of the temperature-related parameters, high-current and self-heating effects.

I-V Characterization and Parameter Extraction
Three UTC-PDs with active areas of 5 × 15 µ m 2 , 5 × 30 µ m 2 and 7 × 25 µ m 2 were characterized under a bias sweep in the range of (−3, 1) V. Current-voltage characterizations were performed on each device geometry at three temperatures: 15, 25 and 35 °C. Figure 6 compares measured the DC characteristics and our developed compact model described in Section 3.2, demonstrating good model accuracy over the entire range of forward and reverse bias at different temperatures. The photodiodes present very low dark current levels of around 4, 20 and 55 nA @ -2V for the geometries 5 × 15 µ m 2 , 5 × 30 µ m 2 and 7 × 25 µ m 2 , respectively, which shows the quality of the fabrication and especially the effective passivation of the photodiode junction. The dark current is very low compared to the photocurrent used for mm-wave generation, which is usually in the mA range to generate sufficient power. Therefore, as the dark current is more than three orders of magnitude below the photocurrent, even if a detailed analysis of the photodiode noise was not realized, we can assume that the impact of a dark current is negligible compared to the global noise of the system. We would like to mention that in order to have a complete view of the impact of the photodiode on the performances of the radio link, one requires the analysis of the nonlinear behavior of the photodiode which is not in the scope of this paper.

Photocurrent Characterization
The photocurrent measurements were performed under an optical power ranging between 0 and 9 dBm by illuminating the UTC-PDs with a laser wavelength of 1.55 µ m. Figure 7 shows the measured and simulated photocurrents for all geometries under test. We estimated a high responsivity of 0.6 A/W for 15 µ m long photodiodes and 0.79 A/W

I-V Characterization and Parameter Extraction
Three UTC-PDs with active areas of 5 × 15 µm 2 , 5 × 30 µm 2 and 7 × 25 µm 2 were characterized under a bias sweep in the range of (−3, 1) V. Current-voltage characterizations were performed on each device geometry at three temperatures: 15, 25 and 35 • C. Figure 6 compares measured the DC characteristics and our developed compact model described in Section 3.2, demonstrating good model accuracy over the entire range of forward and reverse bias at different temperatures. The photodiodes present very low dark current levels of around 4, 20 and 55 nA @ -2V for the geometries 5 × 15 µm 2 , 5 × 30 µm 2 and 7 × 25 µm 2 , respectively, which shows the quality of the fabrication and especially the effective passivation of the photodiode junction. The dark current is very low compared to the photocurrent used for mm-wave generation, which is usually in the mA range to generate sufficient power. Therefore, as the dark current is more than three orders of magnitude below the photocurrent, even if a detailed analysis of the photodiode noise was not realized, we can assume that the impact of a dark current is negligible compared to the global noise of the system. We would like to mention that in order to have a complete view of the impact of the photodiode on the performances of the radio link, one requires the analysis of the nonlinear behavior of the photodiode which is not in the scope of this paper.

I-V Characterization and Parameter Extraction
Three UTC-PDs with active areas of 5 × 15 µ m 2 , 5 × 30 µ m 2 and 7 × 25 µ m 2 were characterized under a bias sweep in the range of (−3, 1) V. Current-voltage characterizations were performed on each device geometry at three temperatures: 15, 25 and 35 °C. Figure 6 compares measured the DC characteristics and our developed compact model described in Section 3.2, demonstrating good model accuracy over the entire range of forward and reverse bias at different temperatures. The photodiodes present very low dark current levels of around 4, 20 and 55 nA @ -2V for the geometries 5 × 15 µ m 2 , 5 × 30 µ m 2 and 7 × 25 µ m 2 , respectively, which shows the quality of the fabrication and especially the effective passivation of the photodiode junction. The dark current is very low compared to the photocurrent used for mm-wave generation, which is usually in the mA range to generate sufficient power. Therefore, as the dark current is more than three orders of magnitude below the photocurrent, even if a detailed analysis of the photodiode noise was not realized, we can assume that the impact of a dark current is negligible compared to the global noise of the system. We would like to mention that in order to have a complete view of the impact of the photodiode on the performances of the radio link, one requires the analysis of the nonlinear behavior of the photodiode which is not in the scope of this paper.

Photocurrent Characterization
The photocurrent measurements were performed under an optical power ranging between 0 and 9 dBm by illuminating the UTC-PDs with a laser wavelength of 1.55 µ m.

Photocurrent Characterization
The photocurrent measurements were performed under an optical power ranging between 0 and 9 dBm by illuminating the UTC-PDs with a laser wavelength of 1.55 µm. Figure 7 shows the measured and simulated photocurrents for all geometries under test. We estimated a high responsivity of 0.6 A/W for 15 µm long photodiodes and 0.79 A/W for 25 and 50 µm-long photodiode. Good accuracy of the photocurrent model, depicted in (4), has been observed. Moreover, good photocurrent scalability is also achieved.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 10 of 3 for 25 and 50 µ m-long photodiode. Good accuracy of the photocurrent model, depicted i (4), has been observed. Moreover, good photocurrent scalability is also achieved.

Thermal Characterization
Reliability issues such as self-heating and the rise of junction temperature are amon the crucial concerns for modern electronic technologies and are becoming increasingl important from the scaling point of view. Additionally, heat-dissipation in photodiode is also a prominent source of concern that may be coupled to the thermal instability in th electronic part, and can thus further deteriorate the performance of an integrated hybri monolithic circuit through thermal feedback between the interconnected devices. Therma management is thus absolutely important for both devices, in order to ensure stable ope ation. In order to alleviate the thermal issues, thermal impedance modeling approache have been proposed for InP heterojunction bipolar transistor (HBT) s [14] including th accurate estimation of thermal resistances and thermal R-C network representations. Sim ilar strategies have been adopted for UTC photodiodes, in this work, to provide an est mation of the thermal resistances and ultimately to develop a unified thermal manage ment strategy for monolithic integrated circuits.

Dynamic Self-Heating and Thermal Impedance Modeling in InP DHBTs
The Low-Frequency S-parameter measurement setup consists of a semiconductor pa rameter analyzer, HP 4155, for DC biasing and a vector network analyzer, Agilent E5061 (5 Hz-3 GHz). In order to couple RF and DC bias, bias tees (bandwidth of 30 kHz to GHz) were used. The dynamic self-heating is observed in the 30 kHz-300 MHz rang where thermal impedance is extracted. A standard short-open-load-thru (SOLT) calibra tion was used with an RF input power of −28 dBm, followed by open-short de-embedding DC bias points were chosen with higher VBE and VCE, where self-heating effects are dom nantly visible. Figure 8 shows the magnitude and the phases of the Y12-parameters of one of the HB geometries (0.7 × 5 µ m 2 ) under test, depicting two distinct zones over the entire frequenc range, separated by a threshold thermal frequency (fTH) around 0.1 GHz. This frequenc denotes a phase shift at the minimum amplitude of the Y12-parameter. The first regio (frequencies higher than fTH) reflects only static self-heating and the behavior in this regio is purely electrical. In the other region, at frequencies below fTH, dynamic self-heating dominant and electro-thermal networks can be used to extract the thermal impedanc The extracted Y parameters are used to calculate the normalized thermal impedance usin

Thermal Characterization
Reliability issues such as self-heating and the rise of junction temperature are among the crucial concerns for modern electronic technologies and are becoming increasingly important from the scaling point of view. Additionally, heat-dissipation in photodiodes is also a prominent source of concern that may be coupled to the thermal instability in the electronic part, and can thus further deteriorate the performance of an integrated hybrid monolithic circuit through thermal feedback between the interconnected devices. Thermal management is thus absolutely important for both devices, in order to ensure stable operation. In order to alleviate the thermal issues, thermal impedance modeling approaches have been proposed for InP heterojunction bipolar transistor (HBT) s [14] including the accurate estimation of thermal resistances and thermal R-C network representations. Similar strategies have been adopted for UTC photodiodes, in this work, to provide an estimation of the thermal resistances and ultimately to develop a unified thermal management strategy for monolithic integrated circuits.

Dynamic Self-Heating and Thermal Impedance Modeling in InP DHBTs
The Low-Frequency S-parameter measurement setup consists of a semiconductor parameter analyzer, HP 4155, for DC biasing and a vector network analyzer, Agilent E5061B (5 Hz-3 GHz). In order to couple RF and DC bias, bias tees (bandwidth of 30 kHz to 3 GHz) were used. The dynamic self-heating is observed in the 30 kHz-300 MHz range where thermal impedance is extracted. A standard short-open-load-thru (SOLT) calibration was used with an RF input power of −28 dBm, followed by open-short deembedding. DC bias points were chosen with higher V BE and V CE , where self-heating effects are dominantly visible. Figure 8 shows the magnitude and the phases of the Y 12 -parameters of one of the HBT geometries (0.7 × 5 µm 2 ) under test, depicting two distinct zones over the entire frequency range, separated by a threshold thermal frequency (f TH ) around 0.1 GHz. This frequency denotes a phase shift at the minimum amplitude of the Y 12 -parameter. The first region (frequencies higher than f TH ) reflects only static self-heating and the behavior in this region is purely electrical. In the other region, at frequencies below f TH , dynamic self-heating is dominant and electro-thermal networks can be used to extract the thermal impedance. The extracted Y parameters are used to calculate the normalized thermal impedance using the following [14]:

Thermal Impedance Extraction Method from Low-Frequency S-Parameters
where Y 22 AC corresponds to Y 22 without dynamic self-heating, i.e., when the thermal capacitance, C TH , is zero, Y 12 DC and Y 22 DC are the Y parameters in the DC condition (i.e., ω→0). Here, the units of all Y parameters are in S and I C , V BE and V CE are expressed in A, V and V, respectively. The extraction of Y 22AC has been performed using the HiCuM model simulation, described in the preceding sections, of the Y 22 parameters with C TH set to zero.
where Y22 AC corresponds to Y22 without dynamic self-heating, i.e., when the thermal capacitance, CTH, is zero, Y12 DC and Y22 DC are the Y parameters in the DC condition (i.e., ω→0).
Here, the units of all Y parameters are in S and IC, VBE and VCE are expressed in A, V and V, respectively. The extraction of Y22AC has been performed using the HiCuM model simulation, described in the preceding sections, of the Y22 parameters with CTH set to zero.

Geometry Scalable Thermal Impedance Model
Two scenarios of the pyramidal downward heat flow through the device have been schematically illustrated in Figure 9a,b. The heat diffusion angle, θ, dictates the temperature distribution along the vertical cross section of the device. For an angle θ < θC (where θC is the critical angle at which the heat flow changes from scenario (a) to (b), with θC = tan −1 [(WCOL−WE)/2Dox]), with WCOL and WE being the lateral widths of the collector and the emitter, respectively, and Dox is the height of the passivation layer, the heat flow is not stopped by the isolation surrounding the collector and the heat diffusion is pyramidal along the depth of the transistor (Figure 9a). On the contrary, for an angle θ > θC, the heat flow is first blocked by the mesa edge, thus leading to a uniform heat diffusion in one subsection (see Figure 9b), followed by a pyramidal diffusion pattern toward the substrate. Based on previous works [14,19], the heat diffusion angle is found to be always lower than 65°, leading to the heat-flow scenario depicted in Figure 9a. Hence, for the rest of the thermal modeling approach, scenario (a) was considered. Following our previous modeling approach [14] to obtain an optimal representation of the behavior of dynamic self-heating in the DHBTs, we considered three major vertical sections of the heat diffusion, which included the collector region beneath the heat source located at the B-C junction, the sub collector and the substrate. The upward heat flow through the emitter metal layers is represented only by a thermal resistance, RTHM, in parallel with the remaining electrical part, as shown in Figure 9c. The emitter-metal layer thermal resistance is calculated in a similar manner as the intrinsic device, considering three different thermal conductivities of the emitter layers: InP, InGaAs (indium gallium arsenide) and TiPdAu [20]. A modified Foster-like thermal network was used to represent heat diffusion through these three sections, as illustrated in Figure 9d. Considering the average thermal conductivity for each section, one can write the following expressions for the thermal resistances of each subsections [14]:

Geometry Scalable Thermal Impedance Model
Two scenarios of the pyramidal downward heat flow through the device have been schematically illustrated in Figure 9a,b. The heat diffusion angle, θ, dictates the temperature distribution along the vertical cross section of the device. For an angle θ < θ C (where θ C is the critical angle at which the heat flow changes from scenario (a) to (b), with θ C = tan −1 [(W COL −W E )/2D ox ]), with W COL and W E being the lateral widths of the collector and the emitter, respectively, and D ox is the height of the passivation layer, the heat flow is not stopped by the isolation surrounding the collector and the heat diffusion is pyramidal along the depth of the transistor (Figure 9a). On the contrary, for an angle θ > θ C , the heat flow is first blocked by the mesa edge, thus leading to a uniform heat diffusion in one subsection (see Figure 9b), followed by a pyramidal diffusion pattern toward the substrate. Based on previous works [14,19], the heat diffusion angle is found to be always lower than 65 • , leading to the heat-flow scenario depicted in Figure 9a. Hence, for the rest of the thermal modeling approach, scenario (a) was considered. Following our previous modeling approach [14] to obtain an optimal representation of the behavior of dynamic self-heating in the DHBTs, we considered three major vertical sections of the heat diffusion, which included the collector region beneath the heat source located at the B-C junction, the sub collector and the substrate. The upward heat flow through the emitter metal layers is represented only by a thermal resistance, R THM , in parallel with the remaining electrical part, as shown in Figure 9c. The emitter-metal layer thermal resistance is calculated in a similar manner as the intrinsic device, considering three different thermal conductivities of the emitter layers: InP, InGaAs (indium gallium arsenide) and TiPdAu [20]. A modified Foster-like thermal network was used to represent heat diffusion through these three sections, as illustrated in Figure 9d. Considering the average thermal conductivity for each section, one can write the following expressions for the thermal resistances of each subsections [14]: here h 1 , h 2 and h 3 denote the depths of the three regions calculated depending on the geometry and the heat flow angle. W E and L E are the effective emitter width and length, respectively, and κ avg1 , κ avg2 and κ avg3 are the average thermal conductivities of the three regions calculated for InP following the approach in [21]. Similarly, the expressions of the thermal capacitances of the three regions are as follows: [ ( + 2ℎ 1 )( + 2(ℎ 1 + ℎ 2 ) ) ( + 2ℎ 1 )( + 2(ℎ 1 + ℎ 2 ) ) ] ,3 = 1 2 ( − ) 3 [ ( + 2(ℎ 1 + ℎ 2 ) )( + 2 ) ( + 2(ℎ 1 + ℎ 2 ) )( + 2 ) ] here h1, h2 and h3 denote the depths of the three regions calculated depending on the geometry and the heat flow angle. WE and LE are the effective emitter width and length, respectively, and κavg1, κavg2 and κavg3 are the average thermal conductivities of the three regions calculated for InP following the approach in [21]. Similarly, the expressions of the thermal capacitances of the three regions are as follows: The complete electro-thermal network including both downward and upward heat diffusion, replaces the one-pole thermal network in HiCuM compact model to estimate the junction temperature more accurately. In Equations (7) and (8)  The complete electro-thermal network including both downward and upward heat diffusion, replaces the one-pole thermal network in HiCuM compact model to estimate the junction temperature more accurately. In Equations (7) and (8), the thermal resistances (R THF,1 , R THF,2 , R THF,3 ) are expressed in K/W, thermal capacitances (C THF,1 , C THF,2 , C THF,3 ) in pJ/K, average thermal conductivities (κ avg1 , κ avg2 and κ avg3 ) in W.m −1 , K −1 , θ in degrees, heat diffusion coefficient D th in m 2 /s and all geometrical dimensions are expressed in m.

Thermal Impedance (Z TH ) Extraction and Model Validation
Low-frequency S-parameters have been measured on six device geometries. Figure 10 illustrates the results comparing the magnitudes and phases of the measured low-frequency Y-parameters (symbols), Y 12 and Y 22 , and the scalable thermal model simulation (solid lines) for different device geometries. The results are shown for a bias condition of V BE = 0.85 V and V CE = 1.2 V, chosen specifically around peak f T of the transistors, where self-heating becomes dominant. Both Y-parameters show good agreement between the measurement and the simulation results. Furthermore, good model scalability has been observed across all geometries.
1, x FOR PEER REVIEW 13 of 31

Thermal Impedance (ZTH) Extraction and Model Validation
Low-frequency S-parameters have been measured on six device geometries. Figure 10 illustrates the results comparing the magnitudes and phases of the measured low-frequency Y-parameters (symbols), Y12 and Y22, and the scalable thermal model simulation (solid lines) for different device geometries. The results are shown for a bias condition of VBE = 0.85 V and VCE = 1.2 V, chosen specifically around peak fT of the transistors, where self-heating becomes dominant. Both Y-parameters show good agreement between the measurement and the simulation results. Furthermore, good model scalability has been observed across all geometries.
(a) (b) Figure 10. Comparison of the magnitudes (a) and phases (b) of the measured low-frequency Yparameters (solid symbols), Y12 and Y22, and the scalable thermal model simulation (solid lines) for different geometries. Figure 11 presents the magnitudes of the thermal impedance, ZTH, extracted using (6), comparing the measurements (solid symbols) and the scalable thermal impedance model simulation (solid lines), depicting good model accuracy and scalability. Figure 11a shows the magnitude of the ZTH for emitter length scaling for the fixed emitter width of (a) 0.7 µ m and (b) 0.5 µ m. In both cases, the model describes the behavior of the dynamic self-heating (below fTH) quite well.
The thermal resistances and capacitances for the three regions of the intrinsic device are shown in Table 1, extracted from compact model simulations, for six geometries with  Figure 11 presents the magnitudes of the thermal impedance, Z TH , extracted using (6), comparing the measurements (solid symbols) and the scalable thermal impedance model simulation (solid lines), depicting good model accuracy and scalability. Figure 11a shows the magnitude of the Z TH for emitter length scaling for the fixed emitter width of (a) 0.7 µm and (b) 0.5 µm. In both cases, the model describes the behavior of the dynamic self-heating (below f TH ) quite well. different emitter lengths and two chosen emitter widths of 0.5 and 0.7 µ m. The total thermal resistance is calculated taking into account both intrinsic devices and the contribution from the emitter metal layers, RTHM. RTH is extracted from the measurement using the intersection technique [14]. The extraction results are shown in Table 2 for all geometries, depicting a very good estimation of the thermal resistances by the developed model.

Thermal Characterization and Modeling of UTC Photodiodes
To take into account the self-heating effect in the UTC photodiodes, a classical thermal network, consisting of a thermal resistance and capacitances, has already been implemented in the compact model, as described in [13]. From the extraction point of view, the The thermal resistances and capacitances for the three regions of the intrinsic device are shown in Table 1, extracted from compact model simulations, for six geometries with different emitter lengths and two chosen emitter widths of 0.5 and 0.7 µm. The total thermal resistance is calculated taking into account both intrinsic devices and the contribution from the emitter metal layers, R THM . R TH is extracted from the measurement using the intersection technique [14]. The extraction results are shown in Table 2 for all geometries, depicting a very good estimation of the thermal resistances by the developed model.

Thermal Characterization and Modeling of UTC Photodiodes
To take into account the self-heating effect in the UTC photodiodes, a classical thermal network, consisting of a thermal resistance and capacitances, has already been implemented in the compact model, as described in [13]. From the extraction point of view, the thermal resistance, R th , of this network can be estimated using the following equations [22]: with where I d and V d are the current (in A) and voltage (in V) across the UTC-PD, N is diode current ideality factor, k B is the Boltzmann constant (in eVK −1 ), T j and T amb are the junction and ambient temperatures (both in K), respectively, q is the electronic charge (in C), I s is the reverse saturation current (in A) and R is the total resistance (in Ω) in series with the photodiode. The parameters in Equation (9) have been extracted under high-current injection (at 1 V) conditions, where self-heating occurs. In order to converge on a solution for T j from both sides of Equation (9), an iteration of the values of R th was performed to solve the expression. In order to facilitate optimal convergence, the model parameters I S , N or R were slightly varied, within an acceptable tolerance that is summarized in Table 3. Subsequently, the extracted values of the thermal resistances of the UTC-PDs are also summarized in Table 3. Interestingly, the extracted values of the thermal resistances of the UTC-PDs are not significantly different (slightly higher) from what can be predicted for these geometries using the R TH scaling law implemented for the DHBTs. Since the two technologies have been developed by the same foundry with similar epitaxial compositions and similar substrate thermal conductivity (InP), the results are coherent and promising for achieving homogeneous monolithic integration with the two devices. Despite the reasonable R TH values obtained for the UTC-PDs, further model improvement and refinement of the extraction methodology lie within the scope of future work. As emphasized earlier, unlike the HBT modeling framework, which has been refined and standardized for a long time, the UTC compact modeling approach is under development and still has several aspects to improve. While the current framework is reasonably accurate, further improvements can include the addition of more material parameters, the description of photocurrent saturation and the more in-depth analysis of self-heating.

High Frequency Characterization and Modelling
Since monolithic integrated optoelectronic circuits are targeted to operate at terahertz frequencies in beyond 5G communication systems, it is crucial to analyze the performances of the electronic and photonic devices at a very high frequency of operation. The difficulties of this approach lie with the conventional test-structures and on-wafer characterization setup that often suffer from characterization technique, calibration and de-embedding techniques adapted to (sub) millimeter-wave frequency bands. The extrapolation of measurements to higher frequencies can often lead to errors and even the models might lack additional effects that can only be understood from actual measurements beyond 200 GHz. Hence, this section sheds light on this aspect through test structure design strategies, in order to correctly "de-embed" the real intrinsic device characteristics, as well as to demonstrate model validity in comparison with the on-wafer measurements performed at very high frequencies.

Dedicated Test Structure Design
An on-wafer thru-reflect-line (TRL) calibration kit was developed and implemented in the InP DHBT process. Its purpose is to allow on-wafer RF transistor measurements beyond 110 GHz. In fact, at such high frequencies, conventional calibration methods fail to provide sufficiently accurate measurements as demonstrated in [23], by comparing short-open-load-thru (SOLT) on commercial calibration substrate on alumina and on-wafer TRL calibrations up to 500 GHz. Thus, on-wafer calibrations are of particular interest for device characterization in the sub-millimeter-wave frequency ranges [24].
In [25], coplanar waveguide (CPW) transmission lines in InP technology proved to be suitable for on-wafer TRL calibration for InP DHBT characterization up to 500 GHz. Hence, the same methodology for designing calibration structures was applied to the III-V Lab InP DHBT process ( Figure 12). As shown in Figure 12, the calibration kit contains one 50 Ω Thru and two 50 Ω Lines of length 180 and 120 µm, covering, respectively, the bandwidths of 35-315 GHz and 70-630 GHz for the TRL calibration. The layout of the RF pads allowed a probing pitch of 50-100 µm. After the on-wafer TRL calibration, the reference plane is set at the inner edge of the pads (see Figure 12b), which is also the beginning of the 50 Ω Thru standard. Hence, the RF pads are included in the calibration error terms.
in order to correctly "de-embed" the real intrinsic device characteristics, as well as to demonstrate model validity in comparison with the on-wafer measurements performed at very high frequencies.

Dedicated Test Structure Design
An on-wafer thru-reflect-line (TRL) calibration kit was developed and implemented in the InP DHBT process. Its purpose is to allow on-wafer RF transistor measurements beyond 110 GHz. In fact, at such high frequencies, conventional calibration methods fail to provide sufficiently accurate measurements as demonstrated in [23], by comparing short-open-load-thru (SOLT) on commercial calibration substrate on alumina and on-wafer TRL calibrations up to 500 GHz. Thus, on-wafer calibrations are of particular interest for device characterization in the sub-millimeter-wave frequency ranges [24].
In [25], coplanar waveguide (CPW) transmission lines in InP technology proved to be suitable for on-wafer TRL calibration for InP DHBT characterization up to 500 GHz. Hence, the same methodology for designing calibration structures was applied to the III-V Lab InP DHBT process ( Figure 12). As shown in Figure 12, the calibration kit contains one 50 Ω Thru and two 50 Ω Lines of length 180 and 120 µ m, covering, respectively, the bandwidths of 35-315 GHz and 70-630 GHz for the TRL calibration. The layout of the RF pads allowed a probing pitch of 50-100 µ m. After the on-wafer TRL calibration, the reference plane is set at the inner edge of the pads (see Figure 12b), which is also the beginning of the 50 Ω Thru standard. Hence, the RF pads are included in the calibration error terms.

Measurement Results and HiCuM Compact Model Validation Up to 330 GHz
The HiCuM model parameter extraction was performed with S-parameter measurements up to 40 GHz and validated against the measurements up to 330 GHz. Among these two separate measurement campaigns, the former was dedicated to parameter extraction while the latter targets high-frequency characterization and model validation. The HiCuM compact model scalability was first verified with RF measurements up to 40 GHz. Five

Measurement Results and HiCuM Compact Model Validation Up to 330 GHz
The HiCuM model parameter extraction was performed with S-parameter measurements up to 40 GHz and validated against the measurements up to 330 GHz. Among these two separate measurement campaigns, the former was dedicated to parameter extraction while the latter targets high-frequency characterization and model validation. The HiCuM compact model scalability was first verified with RF measurements up to 40 GHz. Five different geometries of transistors were measured and the transit frequency, f T , was extracted at 15 GHz for different bias points. The complexity of the HiCuM model requires a systematic extraction of device parameters, especially when different parameters are strongly correlated at high-current conditions. Hence, the interlinked parameter extraction flow described in Section 3 has been used to extract the parameters related to transit time along with the capacitance and current parameters. Moreover, the temperature-related parameters are also extracted for each quantity at all steps. The comparison between the measurement and simulation results is depicted in Figure 13a. A good agreement is observed between the measurements and the geometry-scalable compact model simulation with HiCuM, where the 0.7 × 5 µm 2 InP DHBT presented the maximum RF performance. The corresponding Gummel plot and static gain are shown in Figure 13b. a systematic extraction of device parameters, especially when different parameters are strongly correlated at high-current conditions. Hence, the interlinked parameter extraction flow described in Section 3 has been used to extract the parameters related to transit time along with the capacitance and current parameters. Moreover, the temperature-related parameters are also extracted for each quantity at all steps. The comparison between the measurement and simulation results is depicted in Figure 13a. A good agreement is observed between the measurements and the geometry-scalable compact model simulation with HiCuM, where the 0.7 × 5 µ m 2 InP DHBT presented the maximum RF performance. The corresponding Gummel plot and static gain are shown in Figure 13b. Then, the S-parameter measurements were carried out up to 330 GHz using three measurement benches: (i) up to 110 GHz using an Agilent ® E8361 PNA, with frequency extenders for the 67-110 GHz frequency range; (ii) in the G-band (140-220 GHz) and Jband (220-330 GHz) using a Rohde & Schwarz ® ZVA24, coupled with Rohde & Schwarz ® frequency extenders. Hence, Picoprobe ® probes were used in the 1-110 GHz and J-band, and Cascade Infinity ® probes in the G-band.
As mentioned earlier, the transistor under test has an emitter length of 5 µ m and an emitter width of 0.7 µ m, since it showed the maximum current gain among all the geometries measured. The measurements were made up to 330 GHz using on-wafer TRL calibration. The transistor was biased under a fixed 1.6 V collector-to-emitter voltage and a sweep of the base-to-emitter voltage from 830 to 910 mV with a 10 mV-step. The maximum gain is obtained for VBE = 910 mV.
The gain curves, cut-off frequencies and S-parameters are depicted in Figure 14, where a comparison is drawn between the measured quantities and the ones obtained after simulation of the HiCuM compact model up to 330 GHz. A pretty good agreement can be observed between measurement and simulation up to 330 GHz for the different bias points. Interestingly, a slight discrepancy can be observed between the measurement and the HiCuM simulation for the Mason's gain (Figure 14b). This is due to the limitation of the TRL calibration below 35 GHz. The design of the on-wafer TRL calibration kit focuses on the very high-frequency measurements beyond 110 GHz, and due to the constraints of wafer area, a compromise was made in selecting the line dimensions to be fabricated. Hence, a maximum length of 180 μm (equivalent to a minimum frequency of 35 Then, the S-parameter measurements were carried out up to 330 GHz using three measurement benches: (i) up to 110 GHz using an Agilent ® E8361 PNA, with frequency extenders for the 67-110 GHz frequency range; (ii) in the G-band (140-220 GHz) and J-band (220-330 GHz) using a Rohde & Schwarz ® ZVA24, coupled with Rohde & Schwarz ® frequency extenders. Hence, Picoprobe ® probes were used in the 1-110 GHz and J-band, and Cascade Infinity ® probes in the G-band.
As mentioned earlier, the transistor under test has an emitter length of 5 µm and an emitter width of 0.7 µm, since it showed the maximum current gain among all the geometries measured. The measurements were made up to 330 GHz using on-wafer TRL calibration. The transistor was biased under a fixed 1.6 V collector-to-emitter voltage and a sweep of the base-to-emitter voltage from 830 to 910 mV with a 10 mV-step. The maximum gain is obtained for V BE = 910 mV.
The gain curves, cut-off frequencies and S-parameters are depicted in Figure 14, where a comparison is drawn between the measured quantities and the ones obtained after simulation of the HiCuM compact model up to 330 GHz. A pretty good agreement can be observed between measurement and simulation up to 330 GHz for the different bias points. Interestingly, a slight discrepancy can be observed between the measurement and the HiCuM simulation for the Mason's gain (Figure 14b). This is due to the limitation of the TRL calibration below 35 GHz. The design of the on-wafer TRL calibration kit focuses on the very high-frequency measurements beyond 110 GHz, and due to the constraints of wafer area, a compromise was made in selecting the line dimensions to be fabricated. Hence, a maximum length of 180 µm (equivalent to a minimum frequency of 35 GHz) was chosen for the lines. In addition, due to the high sensitivity of the measurements of the Mason's gain, even the slightest variation in the measurement environment may lead to deviation from what is predicted by the model. Despite the use of three different on-wafer probing setups, an excellent band continuity could be achieved among the measurement frequency bands.
GHz) was chosen for the lines. In addition, due to the high sensitivity of the measurements of the Mason's gain, even the slightest variation in the measurement environment may lead to deviation from what is predicted by the model. Despite the use of three different on-wafer probing setups, an excellent band continuity could be achieved among the measurement frequency bands. Overall, this demonstrates again the accuracy of the HiCuM compact model to accurately reproduce the behavior of InP DHBT transistors in the sub-millimeter-wave frequency range, provided accurate and reliable measurements are performed with on-wafer TRL calibration. Overall, this demonstrates again the accuracy of the HiCuM compact model to accurately reproduce the behavior of InP DHBT transistors in the sub-millimeter-wave frequency range, provided accurate and reliable measurements are performed with on-wafer TRL calibration.

Extrinsic Elements
The pads and the coplanar waveguide (CPW) transmission lines have been modeled using a network of lumped elements located in the extrinsic part of the model. The extrinsic circuit includes the series resistance R p , the inductance, L p and the capacitance, C p .
In order to extract the values of these parasitic lumped elements, UTC-PD test structures in open ( Figure 15a) and short (Figure 15b) configurations are characterized. One-port S-parameter measurements were performed up to 40 GHz with an RF power of -25 dBm using a vector network analyzer (VNA). Figure 16 shows the Smith chart depicting the measured reflection coefficients for the open and short test structures. From the measured real (Re) and imaginary (Im) parts of the Y (in S) and Z 11 (in Ω) parameters of the open and short structures, the values of C p (in F), L p (in H) and R p (in Ω) can be extracted using the following expressions: The average values of resistance R p , C p and L p extracted at high frequency are summarized in Table 4.    The average values of resistance Rp, Cp and Lp extracted at high fr marized in Table 4. The process flow of the on-wafer RF measurements and the param summarized in Figure 17. First, the S-parameters of the open and short measured up to the targeted frequency. Then, one-port S-parameters o measured under various bias conditions. This is followed by an open-sh process to extract the intrinsic frequency-dependent behavior of the UT

De-Embedding (Extraction of Parasitic Elements)
The process flow of the on-wafer RF measurements and the parameter extraction is summarized in Figure 17. First, the S-parameters of the open and short test structures are measured up to the targeted frequency. Then, one-port S-parameters of the UTC-PDs are measured under various bias conditions. This is followed by an open-short de-embedding process to extract the intrinsic frequency-dependent behavior of the UTC-PD.
Appl. Sci. 2021, 11, x FOR PEER REVIEW Figure 17. The flow chart of measurements and parameters extraction process.
The intrinsic admittance of a device (Yint) can be extracted from the measu admittance of the device (Ytot), by subtracting the admittance of the measured op and short (Yshort) test structures. The following extraction step is the model optimi the magnitude and the phase of the S11 under different bias to extract the frequ well as the bias dependence of the junction-related parameters of the UTC-PD. O S-parameters, measured up to 40 GHz under a bias range of 0 to −2 V, are shown Figure 18 for the three geometries under test, comparing the de-embedded expe data and the compact model simulation, affirming a good model accuracy over t The intrinsic admittance of a device (Y int ) can be extracted from the measured total admittance of the device (Y tot ), by subtracting the admittance of the measured open (Y open ) and short (Y short ) test structures. The following extraction step is the model optimization of the magnitude and the phase of the S 11 under different bias to extract the frequency as well as the bias dependence of the junction-related parameters of the UTC-PD. One-port S-parameters, measured up to 40 GHz under a bias range of 0 to −2 V, are shown in Figure 18 for the three geometries under test, comparing the de-embedded experimental data and the compact model simulation, affirming a good model accuracy over the entire frequency and bias range.
The intrinsic admittance of a device (Yint) can be extracted from the measure admittance of the device (Ytot), by subtracting the admittance of the measured open and short (Yshort) test structures. The following extraction step is the model optimiza the magnitude and the phase of the S11 under different bias to extract the freque well as the bias dependence of the junction-related parameters of the UTC-PD. On S-parameters, measured up to 40 GHz under a bias range of 0 to −2 V, are shown in Figure 18 for the three geometries under test, comparing the de-embedded experim data and the compact model simulation, affirming a good model accuracy over the frequency and bias range.

Optical Frequency Response
The optical frequency response of the photodiodes was measured using a heterody bench measurement setup. All UTC-PDs were biased at −2 V under a photocurrent o mA. In addition, 3-dB bandwidths of above 35 GHz for the largest geometry, 7 × 25 µ m and about 50 GHz for the smallest UTC-PD geometry, 5 × 15 µ m 2 , can be observed Figure 19. The model simulation shows good agreement with the measured data valid ing the accuracy of Equation (4) implemented in the model.

Optical Frequency Response
The optical frequency response of the photodiodes was measured using a heterodyne bench measurement setup. All UTC-PDs were biased at −2 V under a photocurrent of 1 mA. In addition, 3-dB bandwidths of above 35 GHz for the largest geometry, 7 × 25 µm 2 , and about 50 GHz for the smallest UTC-PD geometry, 5 × 15 µm 2 , can be observed in Figure 19. The model simulation shows good agreement with the measured data validating the accuracy of Equation (4) implemented in the model.

Optical Frequency Response
The optical frequency response of the photodiodes was measured using bench measurement setup. All UTC-PDs were biased at −2 V under a phot mA. In addition, 3-dB bandwidths of above 35 GHz for the largest geometry and about 50 GHz for the smallest UTC-PD geometry, 5 × 15 µ m 2 , can be Figure 19. The model simulation shows good agreement with the measured ing the accuracy of Equation (4) implemented in the model.

Co-Simulation for Design of Future Electronic-Photonic Integrated Circ Performance Prediction
While the previous sections presented detailed characterization and com ing approaches of the individual electronic and photonic modules, in order predictive, multiphysics simulation framework for the OEICs, in-depth und underlying transport mechanisms is crucial. With that in mind, this section

Co-Simulation for Design of Future Electronic-Photonic Integrated Circuits and Performance Prediction
While the previous sections presented detailed characterization and compact modeling approaches of the individual electronic and photonic modules, in order to develop a predictive, multiphysics simulation framework for the OEICs, in-depth understanding of underlying transport mechanisms is crucial. With that in mind, this section provides insights into these aspects through a complementary TCAD simulation framework, performance evaluation using compact model extrapolation to scaled dimensions and ultimately a unified reliability-aware design methodology for the conception of the OEICs.

TCAD and Monte Carlo Simulations
To strengthen the physical basis of the developed compact models, numerical simulations can be used to unravel the principles of additional effects. For example, the high-frequency S-parameters show a deviation from the expected behavior especially at very high V CE and beyond 200 GHz, in comparison with the measurements obtained at low to medium V CE (around 1V). The observed discrepancies between the experimental measurements and the HiCuM model at high V CE have been investigated using a calibrated multi-scale TCAD model [26]. By combining a full-band, atomistic, and ballistic quantum transport (QT) solver [27] with the hydrodynamic (HD) model of the Sentaurus Device (S-Device) tool [28], qualitatively and quantitatively reliable agreements between the simulated and the experimental data can be produced. Utilizing the configured simulation framework, the influence of a high electric field in the collector region is studied. Figure 20a shows the bulk bandstructure of In 0.53 Ga 0.47 As, as calculated with tight-binding after calibrating the model with the recommended band gap energies of Ref. [29]. The position of the different valleys of In x Ga 1-x As are plotted as a function of the Indium composition x in Figure 20b. Note that to achieve a balanced agreement among all targeted parameters, slight discrepancies exist in the calculated and recommended bandgaps, especially in E X . However, the bowing parameters obtained from our tight-binding approach remains in the recommended range. After calibrating the bandstructure model, the conduction band diagrams of the considered DHBT have been simulated at different V CE with S-Device, as shown in Figure 21a. As V CE goes from 1 to 2 V, the conduction bands in the collector region are pushed downwards and the electric field at the base-collector interface increases. These potential profiles have then been inserted into a QT solver to compute the electron distribution inside the active E-B-C domain. Results are summarized in Figure 21b. For V CE > 1.6 V, a significant increase in the electron concentration inside the collector region is observed. The origin of this phenomenon has been investigated utilizing the energy and spatially resolved electronic distribution inside the simulated structure, as depicted in Figure 21c,d. At low V CE (Figure 21c), most of the electrons are confined in the emitter region, while the base and collector electron concentrations are almost constant and relatively low. However, as the electric field increases with V CE , electrons situated in the Γ valley start to be transferred to two of the six X valleys, those that can be reached without scattering, i.e., those whose principle axis is aligned with the electron transport direction. This can be directly observed in Figure 21d. In contrast to Figure 21c, the region with the highest electron density in Figure 21d is located in both the emitter region and the X valley of the collector. The transfer of electrons from the Γ to the X valleys is commonly observed in III-V materials, where the electronic occupancy probability of the X valleys increases with electric field [29]. It should be noticed that the transition of electrons from the Γ to the L valleys could also contribute to the observed charge increase, but this would require the absorption or emission of a phonon to compensate for the momentum mismatch. The probability for such events to occur is much lower than the elastic transfer of electrons from Γ to X.
A discrepancy between the experimental and theoretical onset of this effect can nevertheless be observed. It can be attributed to the energy level of the X valley, E X , which, on one hand, is not exactly known, and whose behavior, on the other hand, is not accurately reproduced by our tight-binding model.

Extrapolation of Figures of Merits to State-of-the-Art (0.13× 2 um 2 ) DHBTs Using Compact Model
Another important aspect of a co-simulation framework is the predictive capabilities of the compact models. In order to achieve the desired functionality of an integrated circuit, it is imperative that the mathematical model formulation ensures a strong correlation Another important aspect of a co-simulation framework is the predictive capabilities of the compact models. In order to achieve the desired functionality of an integrated circuit, it is imperative that the mathematical model formulation ensures a strong correlation between geometry scaling and model parameters to be able to correctly predict performances for an extrapolated dimension. In order to continually enhance figures of merits along the technology scaling roadmap, various scaling approaches have been adopted or envisioned including process optimization to reduce the capacitances and transit delays [31], or even a hybrid approach for silicon-germanium (SiGe) HBTs consisting doping profile optimization through TCAD simulations and predictive compact modeling through HiCuM/L2 model simulation [32]. With the same fundamental approach for III-V DHBTs, we used geometry scaling laws for HiCuM model parameters to extrapolate the major figures of merit ( f T , f MAX , β and R TH ), without modifying the epitaxial structure or any other process features, to obtain reasonable prediction results. This requires an accurate scalable model card extracted for the geometries under test, which was further used to extrapolate the figures of merit of the current technology generation to DHBTs with the emitter dimensions of 0.13 × 2 µm 2 , in order to compare them with those of the state-of-the-art InP HBT [1,33]. Figure 22 shows the scaling prediction results for these figures of merit using our scalable HiCuM model card for the current technology generation (G1). The extrapolated f T and f MAX are compared with that of the 0.13 × 2 µm 2 state-of-the-art for a V CE of 1 V.
The predicted peak f T [ Figure 22a] shows that a higher value can be expected compared to the reported value for the experimental devices [33]. The new technology generation shows a similar f T compared to the previous generation [34]. The reduction in the transit time with emitter width scaling has been observed at high current densities, which can be attributed to a pronounced collector current spreading in smaller emitter dimensions [34], leading to larger critical current and smaller transit time values. Small emitter dimension particularly poses a challenge to maintain an acceptable current gain (β), of which rather low (15)(16)(17)(18)(19)(20) values have been reported for the 0.13 µm InP HBT process [1,33]. On the other hand, the prediction for the current generation shows an extrapolated current gain of 25-30 for the 0.13 × 2 µm 2 DHBT (Figure 22b). The predicted current gain is comparatively lower than the prediction of the previous generation [33], owing to the higher effective base doping in the new technology generation, which in turn lowers the base sheet resistance and results in a lower static current gain as well as a higher f MAX . Subsequently, in Figure 22c, a lower peak f MAX was observed (even though the new generation exhibits a slightly higher f MAX compared to what had been predicted in [34]), in comparison with the values reported in [33] for the 0.13 × 2 µm 2 state of the art InP HBTs. In fact, the technologies under study have a balanced f T − f MAX performance, and extrapolation to smaller geometries maintains this feature. Despite the improvement in overall f MAX for the new generations, due to higher significance of peak f MAX for circuit applications, further efforts are required to achieve enhanced performances, in order to attain the targeted terahertz operations. These efforts mainly require the modification of the epitaxial base, either at the epitaxial level to induce a higher quasi-electric field [33] or base emitter spacer optimization at the process level [1]. In terms of thermal performances of the DHBT technology, we obtain very similar thermal resistances compared to the preceding generation [34], as depicted in Figure 22d. Similar to our predictions for the previous technology generation, the 0.13 × 2 µm 2 HBT shows significantly elevated values of R TH . Despite appearing as a crucial drawback at such scaled dimensions, we previously observed that the increase in the internal junction temperature (∆T j = R TH × V CE I C | f T, peak ) at peakf T (Figure 22a) (plotted with solid symbols in Figure 22d, right y axis) exhibited a slightly downward trend with scaling. The prediction for the current technology generation (G1) demonstrates that there is an even stronger reduction in the internal device temperatures in the scaled version of this technology. This is due to the fact that the scaling of the DHBTs reduces the V CE I C product predominantly, which in turn results in the drastic reduction in the junction temperature, even though a high thermal resistance was predicted. This is promising from the scaling point-of-view indicating an already improved thermal response has been achieved between subsequent process runs. Approaches such as transferring epitaxial structures on a substrate of higher thermal conductivity [35,36] without modifying the epitaxial structure [36], substrate thinning or even using a quaternary graded base alloy gallium-indium-arsenic-antimonite (GaInAsSb) [37] could be followed in future process runs that can further enhance the figures of merit. Guidelines from the predictive scalable compact model simulations will continue to accompany the future technology generations through further process optimization at the epitaxial level along the InP HBT scaling roadmap.

Toward a Reliability-Aware Design Framework
To emphasize the aforementioned discussions, one of the prerequisites of the monolithic integration of electronic and photonic technologies, that are capable of functioning at terahertz/sub terahertz frequencies, is a systematic assembly of existing design software that will not only take into account multi-physics transport, electro-magnetic simulation as well as interconnect modeling, but will also address emerging reliability concerns. epitaxial structure [36], substrate thinning or even using a quaternary graded base alloy gallium-indium-arsenic-antimonite (GaInAsSb) [37] could be followed in future process runs that can further enhance the figures of merit. Guidelines from the predictive scalable compact model simulations will continue to accompany the future technology generations through further process optimization at the epitaxial level along the InP HBT scaling roadmap.

Toward a Reliability-Aware Design Framework
To emphasize the aforementioned discussions, one of the prerequisites of the monolithic integration of electronic and photonic technologies, that are capable of functioning at terahertz/sub terahertz frequencies, is a systematic assembly of existing design software that will not only take into account multi-physics transport, electro-magnetic simulation as well as interconnect modeling, but will also address emerging reliability concerns.
Optical devices, in particular photodiodes, are prone to significant heat dissipation [38], which may consequently impact the electronic part and can eventually become a serious reliability issue. Catastrophic failure (burn-out) can occur at constant power dissipations of 240 mW for the UTC-PDs [39]. On the other hand, owing to continuous scaling of electronic devices, self-heating is already a severe problem in highly scaled electronic components such as in advanced bipolar transistor technologies including SiGe or InP HBT [40]. As made evident in Section 7.2, thermal resistances (RTH) of DHBTs increase rapidly as device geometries are scaled down. Self-heating management has become the Optical devices, in particular photodiodes, are prone to significant heat dissipation [38], which may consequently impact the electronic part and can eventually become a serious reliability issue. Catastrophic failure (burn-out) can occur at constant power dissipations of 240 mW for the UTC-PDs [39]. On the other hand, owing to continuous scaling of electronic devices, self-heating is already a severe problem in highly scaled electronic components such as in advanced bipolar transistor technologies including SiGe or InP HBT [40]. As made evident in Section 7.2, thermal resistances (R TH ) of DHBTs increase rapidly as device geometries are scaled down. Self-heating management has become the major issue of further HBT process development. The first problem due to increased junction temperature is the impact on its electrical performance. The second problem is reliability related. In fact, increasing junction temperature has a negative impact on the device lifetime [21]. Most of the degradation modes of InP HBTs are strongly accelerated by temperature. Rise in device temperature will ultimately lead to pronounced self-heating at the circuit level and the failure mechanisms are likely to activate leading to degradation in the circuit performances.
To address the emerging reliability concerns in InP-based monolithic integrated circuits, [41] has proposed a reliability-aware design methodology to "pre-diagnose" failure mechanisms in advanced circuits and to design optimized "robust" circuit architectures in order to maximize circuit operation lifetime. The overall goal of this methodology is to provide a holistic reliability analysis solution for next-generation electronic system design that will meet the pre-requisite reliability-constraints prior to manufacturing. This unique methodology is achievable with already available design tools and can be applied regardless of the technology type, thus allowing a straightforward implementation into commercial computer-aided design (CAD) flows. In terms of resources, the methodology aims to work with a single simulator-based design approach, i.e., to offer unified simulation capabilities for all parts of the integrated circuit (IC). In the context of monolithic integrated circuits, multiphysics simulation is thus in high-demand, as emphasized in Section 1. Moreover, this methodology can help drastically shorten design time and costs, while ensuring the same accuracy of results compared to conventional reliability analysis methodologies, by avoiding rigorous circuit aging tests and IC redesign phases for technology qualification. Figure 23 shows the different phases of the IC qualification process which are typically repeated several times in order to meet the qualification criteria. Apart from the obvious time constraints owing to redesign cycles, process optimization and especially mask redesigning costs tens of thousands of dollars for InP-based technologies [41]. The reliability-aware design and optimization at the circuit simulation stage thus allows to carefully choose most "robust" circuit architecture leveraging dynamic degradation laws directly implemented in the compact models compatible with all Berkeley-SPICE derivative simulators. A "virtual degradation acceleration mechanism" is pivotal for ensuring the simulation efficiency of the aging model in order to emulate years of circuit aging by picoseconds of simulation time [41]. While each phase of the reliably aware design methodology can be performed on both electronic and optical parts of the Monolithic Microwave Integrated Circuit (MMIC) until the circuit simulation stage, the reliability improvement short loop is applicable for the entire monolithic integrated circuit. The feasibility and advantages of this methodology has already been demonstrated for an InP/InGaAs DHBT process [41]. The performance of complex circuits such as a transimpedance amplifier, consisting 24 transistors based on this technology, have been analyzed over long-term aging. Design optimization to minimize the circuit failure rate has been demonstrated, making use of the reliability improvement loop. The following steps can be to demonstrate this methodology for photonic components, such as UTC photodiodes, and ultimately to leverage this methodology for the development of monolithic integrated circuits. Owing to the emerging nature of this technology, reliability-aware design will be even more important in order to ensure the desired level of maturity and demonstrate the desired technology readiness level.

Conclusions
In this review, we addressed several of the most important aspects of monolithic integrated circuit design for future applications in beyond-5G communication systems. For the design of such novel technology, InP-based bipolar transistors and UTC photodiodes integrated on the same substrate are viewed as the most promising solution, owing to the respective merits of these two technologies. However, aspects such as the absence of a multiphysics simulation tool, requirement for interconnect electromagnetic modeling at high-frequency operation and tedious thermal management are some of the roadblocks for the current state-of-the art. Most importantly, reliability and thermal management framework at the system level is a prerequisite in the electronic-photonic co-simulation infrastructure. In this paper, design methodologies and modeling approaches as well as

Conclusions
In this review, we addressed several of the most important aspects of monolithic integrated circuit design for future applications in beyond-5G communication systems. For the design of such novel technology, InP-based bipolar transistors and UTC photodiodes integrated on the same substrate are viewed as the most promising solution, owing to the respective merits of these two technologies. However, aspects such as the absence of a multiphysics simulation tool, requirement for interconnect electromagnetic modeling at high-frequency operation and tedious thermal management are some of the roadblocks for the current state-of-the art. Most importantly, reliability and thermal management framework at the system level is a prerequisite in the electronic-photonic co-simulation infrastructure. In this paper, design methodologies and modeling approaches as well as integrated reliability analysis modules are discussed in detail in order to converge towards a multiphysics design approach equipped to tackle these emerging concerns. Extensive characterization, including DC, RF (up to 300 GHz) and optical measurements, have been realized on InP DHBTs and UTC photodiodes, based on the same InP substrate, both technologies developed at III-V Lab. Compared to our previous study [34] where the previous generation of the DHBT technology was investigated up to 40 GHz, the characterization and parameter extraction have been extended to 300 GHz for the new generation of the DHBTs. To accommodate the needs of both time and cost-effective design as well as accuracy of the simulations, computationally efficient as well as physics-based compact models have been developed. While the industry-standard HiCuM model was used to extract the parameters of this new generation of DHBT technology, the UTC-PD compact model developed in [13] has been further improved, including more accurate description of the dark and photo currents, and utilized on the UTC-PD technology from the III-V lab. From an industrial design viewpoint, the models are made compatible with existing electronic design framework. In addition, thermal management and reliability analysis modules are integrated within the models. Particularly, inferences are drawn on the similar values of thermal resistances in the two technologies, owing to their similar epitaxial composition. We have noted, however, a higher rise in junction temperature in the DHBTs compared to that of the UTC-PDs owing to their larger device geometry and consequently more efficient heat evacuation. Nevertheless, further investigations are required in order to quantify the heat dissipation and mutual thermal coupling at more severe operating conditions of the UTC-PDs, which were not addressed in this work. Furthermore, in order to study the expected performances of this technology, model projections are compared with state-of-the-art technologies and their predictive capabilities are demonstrated. Coupled TCAD-QT simulations are leveraged to better understand the underlying physics in order to improve the formulation of analytic models. Finally, an elaborate design methodology was detailed that assembles all the preceding modules under one roof, taking into account compact models, thermal management and reliability as well as failure mechanisms, in order to address the system level design of the monolithic optoelectronic integrated circuits.
Among the aspects discussed earlier, electromagnetic (EM) simulations have proven to be quite essential as the frequency of operation of the active devices fast approaches the sub-terahertz range. As illustrated in [28], as the frequency increases beyond 110 GHz, the electromagnetic coupling between the measurement probes and the substrate can have significant as well as detrimental impacts on the correct assessment of device performance. From the metrology viewpoint, the co-simulation of the characterization setup using electromagnetic and standard SPICE simulations [42] can be used to validate the accuracy of the measurements. Moreover, it can also provide insights into the decorrelation of the parasitics, stemming from the environment of the measurement setup, and the intrinsic device architecture in itself. This method has proven to be particularly essential at the (sub)millimeter wave range as the effects of the test-fixtures and passive components on the intrinsic device correctly "de-embedded" provided on-wafer calibration techniques are utilized. Moreover, along with the thermal management between electronic and photonic modules, intra-chip electromagnetic coupling can be foreseen as another factor to be accounted for in OEIC design. Thus, the SPICE-EM co-simulation can be used to further explore the performances of this key-enabling technology and can be implemented within the proposed co-design tool. Thus, this unique, single simulator-based methodology offers the groundwork for a multiphysics design environment that can be further refined along with the progress of the state-of-the art.