A Large-Signal Graphene FET Model

Chalmers Publication Library Copyright Notice ©2012 IEEE. Personal use of this material is permitted. However, permission to reprint/republish this material for advertising or promotional purposes or for creating new collective works for resale or redistribution to servers or lists, or to reuse any copyrighted component of this work in other works must be obtained from the IEEE. This document was downloaded from Chalmers Publication Library (http://publications.lib.chalmers.se/), where it is available in accordance with the IEEE PSPB Operations Manual, amended 19 Nov. 2010, Sec. 8.1.9 (http://www.ieee.org/documents/opsmanual.pdf) (Article begins on next page) 1 A Large Signal Graphene FET Model Omid Habibpour, Student Member, IEEE, Josip Vukusic, Jan Stake, Senior Member, IEEE Abstract—We propose a semiempirical graphene field effect transistor (G-FET) model for analysis and design of G-FET based circuits. The model describes the current-voltage characteristic for a G-FET over a wide range of operating conditions. The gate bias dependence of the output power spectrum is studied and compared with simulated values. A good agreement between the simulated and the experimental power spectrum up to the 3rd harmonic is demonstrated which confirms the model validity. Moreover, S-parameter measurements essentially coincide with the results obtained from the simulation. The model contains a small set of fitting parameters which can straightforwardly be extracted from S-parameters and DC measurements. The developed extraction method gives a more accurate estimation of the drain and source contact resistances compared to other approaches. As a design example, we use a harmonic-balance load-pull approach to extract optimum embedding impedances for a subharmonic G-FET mixer. Index Terms—Graphene, Microwave FETs, semiconductor device modeling, subharmonic mixer, harmonic balance analysis. I. I NTRODUCTION RAPHENE field effect transistors (G-FETs) are now approaching millimeter wave operation. G-FETs with intrinsic cutoff frequency (fT ) of more than 100 GHz have been demonstrated [1], [2]. Several circuits including frequency doublers [3], [4], fundamental frequency mixers [5], [6] and a subharmonic frequency mixer [7] based on G-FETs have been presented at microwaves frequencies. Moreover, a wafer-scale graphene MMIC mixer has been demonstrated [8]. Hence, in the near future, system-level graphene integrated circuits are a possible reality. Consequently, there is a demand for a device model that can predict a G-FET’s behavior based on its basic transport parameters and can be implemented in standard circuit simulator programs. Device models can be divided into two major groups: physical models and empirical models. Physical models are important in the early stages of device development and they provide better understanding of the device behavior based on its carrier transport parameters. Several physical models for GFETs have been demonstrated [9], [10], [11], [12]. However, since the physical models are based on device physics, they are usually too intricate for circuit level modeling. They are neither fast nor easily implemented for use in circuit design tools. Empirical (semiempirical) models can provide acceptable accuracy with less calculation time and they can be implemented in standard Electronic Design Automation (EDA) tools. Recently several models for G-FETs have been G Manuscript received October 25 and revised December 28, 2011. This research was supported in part by Swedish Foundation of Strategic Research (SSF) and in part by the Knut and Alice Wallenberg (KAW) Foundation. The authors are with the Terahertz and Millimetre Wave Laboratory, Department of Microtechnology and Nanoscience, Chalmers University of Technology, SE-412 96 Göteborg, Sweden (e-mail: [email protected]). proposed [13], [14], [15]. These models can describe the current-voltage characteristics. However, they utilize the same carrier mobility for electrons and holes, and the same contact resistance independent of the carrier type in the channel. Due to the substrate effect, the carrier mobility differs for electrons and holes [7], [16]. Moreover, because of the charge transfer between graphene and metal contacts, G-FETs generally experience different contact resistances depending on the carrier type in the channel [18]. Consequently, there is an asymmetry in the transfer characteristic of G-FETs, especially for short gate-length devices. The above models cannot predict this phenomenon. Lately, an empirical model for a short-channel Si MOSFET has been used to model G-FETs [17]. It allows the calculations of I-V characteristics. Nevertheless, the model uses the same carrier mobility for both electrons and holes. Up to now, all the proposed models are only validated by DC measurements. RF characterization methods, especially power spectrum analysis, reveals more details about the model accuracy and this validation method is missing in the all previous proposed models. In this paper we utilize the semiempirical square-root charge-voltage relation and present an analytical model for G-FETs [19]. This model is derived for a single layer zerobandgap graphene and it can accept different electron and hole mobilities. The model reflects the inherent symmetry of the intrinsic G-FETs and, is similar to the model developed for HEMTs and MOSFETs [20]. Moreover, it can take into account the asymmetric contact resistance for the n-type and p-type channel. The model has only a few fitting parameters and the parameter extraction is performed by means of both S-parameters and DC data. It is experimentally verified for a G-FET under both DC and RF operation. The RF verification includes S-parameters and power spectrum measurements. The results agree well with the model. Furthermore, as a design example, load-pull harmonic balance simulation is performed to extract optimum embedding impedances for a subharmonic G-FET mixer [7]. II. D EVICE M ODEL A. Intrinsic Device The object of this work is to develop a closed-form large signal analytical model for the G-FET, which can be utilized in EDA tools for designing and analyzing G-FET circuits. Since the intrinsic device of the equivalent circuit (Fig. 1) has three terminals, gate, drain, and source, it is possible to express the drain current, Ids as a function of any two voltages between these terminals. We select Vgs and Vgd as independent variables because of the symmetric structure of the G-FETs and since the drain and source are interchangeable. In order to distinguish between intrinsic and extrinsic voltages, capital letters are used for the latter case, i.e VGS . 2 Fig. 1. Large signal model of a G-FET. Cpg , Cpd , Lg , Ld and Ls are pad parasitic capacitances and inductances. Rg is the gate resistances and, Rs and Rd are the source and drain resistances including contact and access resistances. Fig. 2. Majority carrier type in a G-FET channel at 4 different Vgs - Vgd plane quadrants. The type of majority carriers in the G-FET channel can be determined depending on which quadrant of the Vgs − Vgd plane the device bias is : the carriers are electrons for (Vgs > 0, Vgd > 0), both electrons (near the source) and holes (near the drain) for (Vgs > 0, Vgd < 0), both electrons (near the drain) and holes for (near the source) (Vgs < 0, Vgd > 0), and holes for (Vgs < 0, Vgd < 0) (see Fig.2 for illustration). The current in the channel can be expressed as L W (1) ∫ n(x)vdrif t (x) dx L 0 where n(x), vdrif t (x), L and W are the carrier density, the carrier velocity, the channel length and the channel width respectively. The carrier density n(x), should be equal to the thermally generated carriers, nth (8 × 1010 cm−2 at T = 300 K) for disorder-free graphene layer at zero gate voltage. However, the measured data shows that the carrier concentration at zero gate voltage is higher than nth [19]. This is due to the charge impurities located at the graphene/dielectric interface or inside the dielectric, generating extra carriers (electrons or holes) in graphene layer [21]. This effect is modeled in [19] by the following semi-empirical charge-voltage relation, √ (2) n(x) = n20 + (C × V (x)/q)2 Ids = q where n0 is the residual carrier density due to disorder and thermal excitation. This is the minimum charge concentration of the graphene layer and C = (Cgs + Cgd )/(LW ) is the gate capacitance per area. The total gate capacitance consists of the gate oxide capacitance, Cox in series with the graphene quantum capacitance, Cq . In this model, since Cq ≫ Cox , the effect of Cq is ignored. This is valid except for ultra thin gate dielectrics [24]. We note that above equation reduces to the familiar n(x) = C ×V (x)/q at high gate voltage, and to n = n0 at zero gate voltage. Moreover, the carrier drift velocity is approximated by the following velocity saturation model [22], [23]. µE(x) vdrif t (x) = √ m )m 1 + ( µ∣E(x)∣ vsat (3) where vsat is the saturation velocity of the carrier, µ is the carrier mobility and m is a fitting parameter. The vsat depends √ on carrier concentration in this model with vsat (n) = vF β/ n where vF = 108 cm/s and β relates to the optical phonon wavelength of the dominant scattering phonon and for graphene on SiO2 β = 4×105 cm−1 [13]. The above expression is valid for n > 1.1×1011 cm−2 [23] where the minimum value of n is n0 . To the best of the authors knowledge the minimum reported value of n0 is 2.2 × 1011 cm−2 [19] in supported and top gated graphene. With this range of n0 , vsat < vF , although generally the quoted values of n0 are significantly higher. We also assume a constant carrier mobility for electrons and holes (µe and µh ). In the first quadrant( Vgs > 0, Vgd > 0) Ids = q L W µE(x) dx ∫ n(x) √ L 0 µ∣E(x)∣ m m 1 + ( vsat ) (n(x)) (4) knowing that dV = E(x)dx Ids = q Vgs µn(V ) W dV √ ∫ L Vgd m )∣ m 1 + ( µ∣E(V ) vsat (V ) (5) For simplicity, the velocity saturation value at the average gate √ voltage is used in the above expression v sat = vF β/ 4 n20 + (C(Vgs + Vgd )/2q)2 and also ∣E(V )∣ is approximated by ∣Vgs − Vgd ∣/L and √ µe W Ids1 = √ (Vgs Q20 + (CVgs )2 − µ ∣V −V ∣ m 1 + ( e Lvgssat gd )m L √ Vgd Q20 + (CVgd )2 + √ Q2 + (CVgs )2 + CVgs 2 +q /C ln √ 0 ) Q20 + (CVgd )2 + CVgd (6) where Q0 (q × n0 ) is the residual charge density. In order to simplify the equation notation, the following unitless function 3 is defined: √ √ f (x, y) = x 1 + x2 − y 1 + y 2 √ 1 + x2 + x + ln √ 1 + y2 + y (7) and defining V̄gs = Vgs /V0 and V̄gd = Vgd /V0 where V0 = Q0 /C and therefore, W µ e V0 Q 0 f (V̄gs , V̄gd ) Ids1 = √ µe ∣Vgs −Vgd ∣ m L m 1 + ( Lvsat ) (8) In the second quadrant ( Vgs > 0, Vgd < 0) there is a point in the channel ( x = L1 ) where the voltage is the same as the gate voltage. At this point, the type of majority carriers changes from electrons to holes. Ids = W (∫ L 0 L1 e × n(x)vdrif t,e dx + ∫ L L1 e × n(x)vdrif t,h dx) (9) after integration we have µ e V0 Q 0 Ids2 = √ µ ∣V −V ∣ m 1 + ( e Lvgssat gd )m µ h V0 Q 0 + √ µ ∣V −V ∣ m 1 + ( h Lvgssat gd )m W f (V̄gs , 0) L W f (0, V̄gd ) L (10) In the third quadrant ( Vgs < 0, Vgd > 0) using the same procedure as before µ h V0 Q 0 Ids3 = √ µ ∣V −V ∣ m 1 + ( h Lvgssat gd )m µ e V0 Q 0 + √ µ ∣V −V ∣ m 1 + ( e Lvgssat gd )m W f (V̄gs , 0) L W f (0, V̄gd ) L (11) and finally for the forth quadrant ( Vgs < 0, Vgd < 0) we have W µ h V0 Q 0 Ids4 = √ f (V̄gs , V̄gd ) µh ∣Vgs −Vgd ∣ m L m 1 + ( Lvsat ) (12) The above equations describe Ids based on the Vgd and Vgs in different quadrants. We combine the above equations and write the Ids for all bias conditions as Ids = Ids1 Θ(Vgs )Θ(Vgd ) + Ids2 Θ(Vgs )Θ(−Vgd ) Ids3 Θ(−Vgs )Θ(Vgd ) + Ids4 Θ(−Vgs )Θ(−Vgd ) Fig. 3. Contour plot of drain current (Ids (mA/µm)) of a G-FET versus the intrinsic Vgs and Vgd . µe = 2000 cm2 /Vs, µh = 2400 cm2 /Vs, n0 = 6 × 1011 cm−2 , L = 1µm In the above expressions we assume that there is no unintentional charging (doping) in the channel, i.e. VDirac = 0 at low drain voltages (Vds ≪ 1). Hence, in order to include this effect above equation should be modified by putting Vgs −VDirac and Vgd − VDirac instead of Vgs and Vgd respectively. The convergence of some simulation techniques such as harmonic-balance and transient analysis, requires continuous first- and second-order derivatives [25]. Consequently, the model should have continuous high-order derivatives. Since n f (x,y) ∂ n f (x,y) , ∂yn ), f (x, y) has high-order partial derivatives ( ∂ ∂x n by replacing Θ(x) with a smooth analytic function in equation 13, Ids will have high order continuous derivatives. For that reason, Θ(x) ≃ U (x) = (1 + tanh(x/V1 ))/2 where V1 is a fitting parameter. Finally, from above equations the transconductance gm = dId /dVgs ∣Vds =const of the intrinsic device can be approximated in the different quadrants as, 2µe,h W √ 2 2 − gm1,4 = √ C( V0 + Vgs µe,h ∣Vds ∣ m L m 1 + ( Lvsat ) √ V02 + (Vgs − Vds )2 ) (14) 2µe,h W √ 2 2 − gm2,3 = √ C V0 + Vgs µe,h ∣Vds ∣ m L m 1 + ( Lvsat ) 2µh,e W √ 2 C V0 + (Vgs − Vds )2 √ µh,e ∣Vds ∣ m L m 1+( ) (15) Lv sat (13) where Θ(x) is the step function. Fig. 3 shows a contour plot of the drain current of a G-FET based on intrinsic Vgs and Vgd with µe = 2000 cm2 /Vs, µh = 2400 cm2 /Vs, n0 = 6 × 1011 cm−2 , L = 1µm and m = 1. From above equations, gm,M ax = 2µe,h W √ µ ∣Vds ∣ m L 1+( e,h ) Lv m sat C∣Vds ∣ and for high carrier mobility, gm,M ax ≃ 2v sat W C. It should be noted that this value is obtained from the intrinsic part (without including parasitic drain and source resistances) and therefore it is an upper limit of the G-FET transconductance. 4 B. Drain and Source Resistance Model Measurement data [18], [26] shows that there is a difference in the conductivity of G-FETs depending on whether the majority carriers in the G-FET channel are electrons ( Vgs > 0, Vgd > 0) or holes ( Vgs < 0, Vgd < 0). This effect mainly comes from the charge transfer between the graphene and the metal contact. The charge type can be electrons or holes depending on the work function of the deposited metal [27]. When the carrier type of the G-FET channel is the opposite of the charge in the metal-graphene contact, a p-n barrier is formed and adds an extra resistance [28]. For submicron gatelength G-FETs, this phenomenon dominates the whole channel resistance [2], [26]. Therefore including this term is necessary for G-FET modeling. In order to model this effect we add a carrier dependent term to the drain and source resistors, Rs = Rs0 + Rext (Vgs , Vgd ) Rd = Rd0 + Rext (Vgs , Vgd ) (16) (17) 1 + tanh(Vgs /V2 ) 1 + tanh(Vgd /V2 ) ×Rexto 2 2 In the above expression when the majority carriers in the GFET channel are electrons, an extra resistance, Rexto appears in the drain and source resistances and V2 is a fitting parameter. Since the metal-graphene contact scales with the contact width rather than the contact area [28] and the drain and source have the same contact width, we can assume Rs0 =Rd0 =R0 . Since in the large area graphene there is no bandgap, the contact resistance becomes more critical as the gate length decreases. The main effect than can be seen in very short gatelength (Lg < 100-150 nm) G-FETs, is that by changing the gate voltage, the current variation for n- (p-) channel becomes much smaller than the current variation for p- (n-) channel [29], [2]. In other words G-FET behaves like a unipolar FET. A high value of Rexto can be used for modelling the above effect. Also in [30] it has been shown that the saturation velocity is almost independent of the gate-length while the carrier mobility generally decreases when reducing the gatelength. Rext (Vgs , Vgd ) = achieved by this method do not include the contact and access resistances. These resistances are the main contributions to the parasitic resistances and as a result, this method underestimates the contact resistances. For example in [2] the extracted drain and source resistances for graphene of about 3 µm width are 13.9 Ω (42 Ω.µm) and 2.4 Ω (7 Ω.µm) respectively. These values are more than an order of magnitude lower than the lowest reported contact resistance for a single layer graphene at room temperature (500 Ω.µm) [28]. Consequently we modify the the method as follows: 1. Open and short structures with identical layouts, excluding the graphene are used to extract Ls , Lg , Ld , Cpg , Cpd and Rg . 2. For extraction of Rd and Rs (R0 and Rext ), the drain is biased at a low voltage. At low drain voltage (∣Vds /Vgd ∣,∣Vds /Vgs ∣ ≪ 1) equation 8 and 12 can be writ√ 1 + (V /V0 )2 Vds and IDS = ten as IDS = µ (W /L)Q gs e 0 √ 2 µh (W /L)Q0 1 + (Vgs /V0 ) Vds respectively and as a result, the drain to source resistance becomes α µe (18) RDS = 2R0 + 2Rexto + √ 1 + (Vgs /V0 )2 for Vgs ≫ VDirac and αµh RDS = 2R0 + √ 1 + (Vgs /V0 )2 (19) for Vgs ≪ VDirac where αµe,h = L/(W µe,h Q0 ). By fitting the Rds profile with above equations Rd and Rs (R0 and Rexto ) as well as αµe,h and V0 can be extracted. 3. The capacitors of the intrinsic device can be obtained by biasing the gate voltage of the G-FET at the minimum conductivity (VGS = VDirac ) and measuring the S-parameters. Since the device is biased at VDirac , the transconductance, gm becomes zero. The parasitic elements obtained from the previous sections can be deembedded [33] and the remaining small signal model is depicted in Fig. 4. It can be shown that III. PARAMETER E XTRACTION Device parameters are extracted from both DC and Sparameters measurements. The measured S-parameters are used in order to extract the extrinsic parasitic elements as well as intrinsic capacitors (Fig. 1). For parasitic elements, S-parameters of the open and short structures are needed. In conventional FETs, the open and short structures are obtained by biasing the device under the forward-biased gate and pinched off respectively (VDS = 0, cold-FET technique) [31], [32]. However, the G-FET channel can not be pinched off. Consequently, the method must be modified for the GFET. As described in [2], ’open’ and ’short’ structures with identical layouts, excluding graphene are used to deembed the parasitic elements and calculate the intrinsic elements by Yintrinsic = ((YDU T − Yopen )−1 − (Yshort − Yopen )−1 )−1 . However, since the ’open’ and ’short’ structures do not include graphene in channel, the parasitic drain and source resistances Fig. 4. Small signal model after deembedding the parasitic elements at VGS = VDirac . Cgd = −Y12 /jω, Cgs = (Y11 + Y12 )/jω, (20) Cds = Im(Y22 + Y12 )/ω (21) C = (Cgs + Cgd )/(LW ) (22) where Y is the admittance matrix of the structure in Fig. 4. 4. The rest of the parameters are calculated as follow: Q0 = CV0 µe,h = L/(W αµe,h Q0 ) (23) (24) 5 TABLE I MODEL PARAMETERS Element Cgs Cgd Cds Cpd Cpg Ls Ld value 25 fF 22 fF 48 fF 18 fF 20 fF 31 pH 43 pH Element Lg Rg Ro Rexto µe µh V0 value 75 pH 12 Ω 28 Ω 5Ω 2000 cm2 /V s 2400 cm2 /V s 0.41 V IV. R ESULTS The model has been implemented in a circuit simulation software (ADS) and verified by the experimental DC and RF characterization of a G-FET. The G-FET fabrication process is the same as the method described in [7]. By the method described in the previous section, the model’s parameters have been extracted and presented in table I (Wch = 20µm and Lch = 1µm) with fitting parameters of m = 1, V1 = 1/3 V and V2 = 1/4 V. Fig. 5. Fig. 6. Fitted data (solid line) versus measured data (⧫) for the drain-source resistance (VDirac = 1 V, VDS = 0.1 V ). SEM image of the G-FET, Lg = 1µm , Wg = 20µm. The measured G-FET drain-source resistance versus the fitted model is shown in Fig. 6. The transconductance and IDS − VDS characteristic curves at different VGS are depicted in figures 7 and 8 respectively. As can be seen the model is in good agreement with the measurement. For the IDS − VDS characterization, a pulsed I-V measurement was used to avoid the self heating effect at high drain voltages [34]. The simulated voltage range has been extended beyond the experiment range in order to show the formation of p-n junction along the channel. For example at VGS − VDirac = 2 V, when VDS reaches 2 V a p-n junction starts to form. Moreover, measured and modeled S-parameters are shown in Fig. 9. The measurement is performed from 1 GHz to 30 GHz with an Agilent PNA with on-wafer probing. As described in the previous section, for extracting the capacitors of the intrinsic device, S-parameters of the G-FET biased at VGS = VDirac are needed. At this bias point, gm = 0, and consequently S12 = S21 . Therefore, VDirac can be found by tuning the gate voltage until S12 and S21 completely coincide. For example for VDS = 0.1 and 0.5 V in this device, we have VDirac = 1 and 1.22 V respectively (Fig. 9 a,b). Fig. 9 c,d show the measured and modeled S-parameters for the maximum gm with VDS = 0.1 and 0.5 V. In order to investigate the model accuracy more in detail, a power spectrum analysis is performed. This is done by Fig. 7. Model (solid line) versus measured data (⧫) for the transconductance (VDirac = 1 V, VDS = 0.1 V ). sweeping the gate bias voltage and superimposing a low frequency (10 MHz) signal to the gate, and measuring the power spectrum at the drain. Harmonic balance analysis is used to simulate the power spectrum using the model. Fig. 10 demonstrates a good agreement between the simulated and measured power spectrum up to the third harmonic. Also, as can be seen from the first harmonic, there is no power gain in this G-FET. V. D ESIGN EXAMPLE Here, we show how the model can be used for designing and analyzing a G-FET based circuit. Fig. 11 shows the circuit schematic of a subharmonic resistive G-FET mixer [7]. The input-signal frequencies are 2 GHz (fRF ) and 1.01 GHz (fLO ) and the desired output-signal frequency (fIF = ∣fRF −2×fLO ∣) is 20 MHz. By running a harmonic balance load-pull simulation, the optimum RF and IF embedding impedances for a given LO power, PLO can be found. Fig. 12 depicts contour plots showing simulated conversion loss (CL) for RF and IF 6 Fig. 10. Model (solid line) and measured power spectrum. fin = 10 MHz, Pin = 0 dBm, VDS = 0.5 V. Fig. 8. Model (solid line) versus measured data (⧫) for IDS − VDS characteristic curves at VGS − VDirac =-3 to 3 volts. The measurement is based on the pulsed-IV method. Fig. 11. A subharmonic G-FET mixer circuit [7]. Fig. 9. S-parameters Model (●) versus measurement (○), a) VGS = 1 V , VDS = 0.1 V, gm = 0, b) VGS = 1.22 V , VDS = 0.5 V, gm = 0 c) VGS = 1.5 V , VDS = 0.1 V, gm = 27 µS/µm , d) VGS = 1.75 V , VDS = 0.5 V, gm = 137 µS/µm impedance levels at PLO = 15 dBm. The plot is achieved by the following steps. First, the global optimum embedding impedances are estimated by sweeping one impedance (e.g. ZRF ) in the Γ-plane while keeping the other impedance (e.g. ZIF ) at Z0 = 50 Ω. Then, ZRF,opt (ZIF,opt ) is obtained by sweeping ZRF (ZIF ) in the Γ-plane with ZIF (ZRF ) kept at the estimated global optimum value. Fig. 13 shows ZRF,opt , ZIF,opt , the simulated and measured mixer CL versus PLO . There is a good agreement between the measured and simulated CL. The optimum impedances are located near the Fig. 12. Simulated CL contour plots for RF and IF embedding impedances. The inner contour is for 22.8 dB CL and the outer contour is for 30 dB CL (0.5-dB step), PLO = 15 dBm and, ZLO and all other mixing terms are terminated with Z0 (VGS = VDirac = 1 V). real axis of the Γ-plane, therefore the real values are plotted. The optimum RF and IF embedding impedances essentially 7 coincide and it is due to the low RF frequency. It can be seen that by using the optimum embedding impedances at low PLO levels, a lower CL can be reached. Also, it is seen that at high PLO , the CL starts to increase. This is because at high PLO , the time duration of ’off state’ becomes much shorter than that of ’on state’. In conventional subharmonic FET mixers, the device can be biased such that the time duration of the ’on state’ and ’off state’ becomes close together, and consequently the mixer CL monotonically decreases by increasing PLO and saturates [35], [36]. The demonstrated mixer exhibits a high CL which is attributed to the low on-off current ratio (≃ 3) of the G-FET. Further simulation shows that by utilizing a GFET with a on-off ratio of 10-20 in the mixer, a CL of 15-17 dB is achievable [37]. A higher on-off ratio can be obtained by reducing the contact resistance levels as well as creating a bandgap in the G-FET channel [38]. Fig. 13. Simulated and measured CL, as well as the optimum embedding impedances versus PLO (ZLO = Z0 ). The measured data is compensated for filter losses and the device parameters are the same as those in table I (VGS = VDirac = 1 V). VI. C ONCLUSION We have proposed and evaluated a closed form large signal model for the G-FETs. The model is semiempirical and is derived for a single layer zero-bandgap graphene. The model accepts different carrier mobilities for electrons and holes and it can predict the asymmetric transfer characteristic of a GFET. Using a semiempirical approach reduces the complexity of the calculations and enables us to have an analytical model suitable for circuit-level simulations. The model has few fitting parameters which can straightforwardly be extracted using a novel method. This new extraction method gives a more accurate estimation of the drain and source contact resistances than previously reported methods. The model is experimentally verified at DC and RF. The DC measurements agree well with the model and also the power spectrum analysis shows good agreement up to the 3rd order. As a practical example for the model, a harmonic-balance load pull analysis is performed to extract the optimum embedding impedances of a G-FET subharmonic mixer. In this case the measured and simulated mixer CL are compared, and the simulated result follows the measured data. This example shows how the model can be used to analyse and design G-FET circuits. For future model development, the effect of self-heating should be considered. Moreover, for ultra-thin gate dielectrics, where the effect of the graphene quantum capacitance is no longer negligible, the model should be modified. ACKNOWLEDGMENT The authors would like to thank Prof. Erik Kollberg, Prof. Olof Engström and Micheal Andersson (Chalmers University of Technology) for fruitful discussions. R EFERENCES [1] Y.-M. Lin, C. Dimitrakopoulos, K. A. Jenkins, D. B. Farmer, H.-Y. Chiu, A. Grill and Ph. Avouris ”100-GHz transistors from wafer-scale epitaxial graphene,” Science, 327, 662, 2010. [2] Lei Liao, Yung-Chen Lin, Mingqiang Bao, Rui Cheng, Jingwei Bai, Yuan Liu, Yongquan Qu, Kang L. Wang, Yu Huang and Xiangfeng Duan, ”High-speed graphene transistors with a self-aligned nanowire gate,” Nature, 467, 305-308, 2010. [3] H. Wang, D. Nezich, J. Kong, and T. Palacios, ”Graphene frequency multipliers,” IEEE Electron Device Letters, vol. 30, no. 5, pp. 547-549, May 2009. [4] H. Wang, A. Hsu, K. K. Kim, J. Kong, and T. Palacios, ”Gigahertz ambipolar frequency multiplier based on CVD graphene,” in IEDM Tech. Dig., 2010, pp. 572-575. [5] H. Wang, A. Hsu, J. Wu, J. Kong, and T. Palacios, ”Graphene-based ambipolar RF mixers,” IEEE Electron Device Letters, vol. 31, no. 9, pp. 906-908, Sep. 2010. [6] L. Liao, J. Bai, R. Cheng, H. Zhou, L. Liu, Y. Liu, Y. Huang, and X, Duan, ”Scalable Fabrication of Self-Aligned Graphene Transistors and Circuits on Glass,” Nano Letters, DOI: 10.1021/nl201922c, June 7, 2011. [7] O. Habibpour, S. Cherednichenko, J. Vukusic and J. Stake, ”A subharmonic graphene FET mixer,” IEEE Electron Device Letters, vol. 33, no. 1, pp. 71-73, 2012. [8] Y. M. Lin, A. V. Garcia, S. J. Han, D. B. Farmer, I. Meric, Y. Sun, Y. Wu, C. Dimitrakopoulos, A. Grill, P. Avouris, K. A. Jenkins ”Wafer-Scale Graphene Integrated Circuit,” Science, 332, 1294, 2011. [9] M. Cheli, G. Fiori, and G. Iannaccone, ”A semianalytical model of bilayer graphene field effect transistor,” IEEE Transactions on Electron Devices, vol. 56, no. 12, pp. 2979-2986, Dec. 2009. [10] S. Thiele, and F. Schwierz, ”Modeling of the steady state characteristics of large-area graphene field-effect transistors,” Journal of Applied Physics, vol. 110, 034506, 2011. [11] J. G. Champlain, ”A first principles theoretical examination of graphenebased field effect transistors,” Journal of Applied Physics, vol. 109, 084515, 2011. [12] S. O. Koswatta, A. Valdes-Garcia, M. B. Steiner, L. Yu-Ming Lin and P. Avouris, ”Ultimate RF Performance Potential of Carbon Electronics”, IEEE Transactions on Microwave Theory and Techniques, vol. 59, no. 10, pp. 2739-2750, 2011. [13] I. Meric, M. Y. Han, A. F. Young, B. Ozyilmaz, P. Kim, and K. L. Shepard, ”Current saturation in zero-bandgap, top-gated graphene fieldeffect transistors,” Nature Nanotechnology, vol. 3, no. 11, pp. 654-659, Nov. 2008. [14] D. Jimenez and O. Moldovan, ”Explicit Drain-Current Model of Graphene Field-Effect Transistors Targeting Analog and Radio-Frequency Applications”, IEEE Transactions on Electron Devices, vol. 58, no. 11, pp. 4049-4052, 2011. [15] D. Jimenez and O. Moldovan, ”Modeling of the Output and Transfer Characteristics of Graphene Field-Effect Transistors”, IEEE Transactions on Electron Devices Nanotechnology, vol. 10, no. 5, pp. 1113-1119, 2011. [16] O. Habibpour, S. Cherednichenko, J. Vukusic and J. Stake, ”Mobility improvement and microwave characterization of a graphene field effect transistor with silicon nitride gate dielectrics,” IEEE Electron Device Letters, vol. 32, pp. 871-873, 2011. [17] H. Wang, A. Hsu, J. Kong,D. A. Antoniadis, and T. Palacios,”Compact Virtual-Source CurrentVoltage Model for Top- and Back-Gated Graphene Field-Effect Transistors”, IEEE Transactions on electron devices, vol. 58, no. 5, pp. 1523-1533, May 2011. 8 [18] B. Huard, N. Stander, J. A. Sulpizio, and D. Goldhaber-Gordon, ”Evidence of the role of contacts on the observed electron-hole asymmetry in graphene,” Physical Review B, 78, 121402, 2008. [19] S. Kim, J. Nah, I. Jo, D. Shahrjerdi, L. Colombo, Z. Yao, E. Tutuc. and S. K. Banejee, ”Realization of a High Mobility Dual-gated Graphene Field Effect Transistor with Al2O3 dielectric,” Applied Physics Letters, vol. 94, 062107, 2009. [20] K. Yhland, N. Rorsman, M. Garcia, H.F. Merkel, ”A Symmetrical Nonlinear HFET/MESFET Model Suitable for Intermodulation Analysis of Amplifiers and Resistive Mixers,” IEEE Transactions on Microwave Theory and Techniques, vol. 48, no. 1, pp. 15-22, Jan 2000. [21] S. Adam, E. H. Hwang, V. M. Galitski, and S. Das Sarma, ”A selfconsistent theory for graphene transport,” Proc. Natl. Acad. Sci. U.S.A. 104, 18392 (2007). [22] J. Chauhan and J. Guoa, ”High-field transport and velocity saturation in graphene,” Applied Physics Letters, vol. 95, 023120, 2009. [23] V. E. Dorgan, M. H. Bae, and E. Pop, ”Mobility and saturation velocity in graphene on SiO2,” Applied Physics Letters, vol. 97, 023120, 2010. [24] J. Xia, F. Chen, J. Li and N. Tao, ”Measurement of the quantum capacitance of graphene,” Nature Nanotechnology, vol. 5, pp. 505 - 509, 2009. [25] S. A. Maas, Nonlinear Microwave and RF circuits, 2nd ed. Artech House, 2003. [26] W. J. Liu, H. Y. Yu, S. H. Xu, Q. Zhang, X. Zou, J. L. Wang, K. L. Pey, J. Wei, H. L. Zhu, and M. F. Li, ”Understanding Asymmetric Transportation Behavior in Graphene Field-Effect Transistors Using Scanning Kelvin Probe Microscopy,” IEEE Electron Device Letters, vol. 32, no. 2, pp. 128-130, Feb. 2011. [27] E. H. Lee, K. Balasubramanian, R. T. Weitz, M. Burghard, and K. Kern, ”Contact and edge effects in graphene devices,” Nature Nanotechnology, 3, pp. 486 - 490, 2008. [28] K. Nagashio, T. Nishimura, K. Kita and A. Toriumi, ”Metal/graphene contact as a performance Killer of ultra-high mobility graphene analysis of intrinsic mobility and contact resistance,” IEEE International Electron Devices Meeting, pp. 565-568, Dec. 2009 [29] S. J. Han, Z. Chen, A. A Bo and Y. Sun, ”Channel-Length-Dependent Transport Behaviors of Graphene Field-Effect Transistors,” IEEE Electron Device Letters, vol. 32, no. 6, pp. 812-814, May, 2011. [30] I. Meric, C. R. Dean, A. F. Young, N. Baklitskaya, N. J. Tremblay, C. Nuckolls, Ph. Kim, Philip and K. L. Shepard, ”Channel Length Scaling in Graphene Field-Effect Transistors Studied with Pulsed CurrentVoltage Measurements,” Nano Letters, vol. 11, no. 3, pp. 1093-1097, 2011. [31] G. Dambrine, A. Cappy, F. Heliodore, and E. Playez, ”A new method for determining the FET small-signal equivalent circuit,” IEEE Transactions on Microwave Theory Tech., vol. 36, pp. 1151-1159, July 1988. [32] R. Tayrani, J. E. Gerber, T. Daniel, R. Pengelly, S. Rayrnond and U. L. Rohde, ”A new and reliable direct parasitic extraction method for MESFETs and HEMTs,” 23rd European Microwave Conference, 1993. [33] R. Bauer and Paul Penfield, ”De-embedding and unterminating,” IEEE Transactions on Microwave Theory and Techniques, vol. 22, no. 3, pp. 1151-1159, March 1974. [34] M. Freitag, M. Steiner, Y. Martin, V.i Perebeinos, Zh. Chen, J. C. Tsang, and Ph. Avouris, ”Energy Dissipation in Graphene Field-Effect Transistors,” Nano Letters, vol. 9, no. 5, pp. 1883-1888, 2009. [35] H. Zirath, ”A subharmonically Pumped Resistive Dual-HEMT Mixer,” IEEE MTT-S Digest, pp. 875-878, 1991. [36] K. Yhland ”Simplified Analysis of Resistive Mixers,” IEEE Microwave and Wireless components letters, vol. 17, no. 8, pp. 604-606, 2007. [37] O. Habibpour, ”Fabrication, Characterisation and Modelling of Subharmonic Graphene FET Mixers,” Licentiate dissertation, Chalmers Univ. Technol., Göteborg, Sweden, 2011. [38] M. Sprinkle, M. Ruan, Y. Hu, J. Hankinson, M. Rubio-Roy, B.Zhang, X. Wu, C. Berger, and W. A. de Heer, ”Scalable templated growth of graphene nanoribbons on SiC,” Nature Nanotechnology, vol. 5, pp. 727731, 2010. Omid Habibpour (S’08) received the B.S degree in electrical engineering from Sharif University of Technology, Tehran, Iran, in 2002 and the M.S. degree in electrical engineering from Amirkabir University of technology, Tehran, Iran, in 2004. He is currently working toward the Ph.D. degree in graphene electronics with the Terahertz and Millimetre Wave Laboratory, Department of Microtechnology and Nanoscience, Chalmers University of Technology, Göteborg, Sweden. Josip Vukusic received the diploma and Ph.D. degree in photonics from Chalmers University of Technology, Göteborg, Sweden, in 1997 and 2003, respectively. Since 2004, he has been with the Terahertz and Millimetre Wave Laboratory working on terahertz (THz) technology. His is currently involved in modeling, fabrication, and characterization of frequency multipliers and photomixers for THz generation. Jan Stake (S’95- M’00 - SM’06) was born in Uddevalla, Sweden, in 1971. He received the M.Sc. degree in electrical engineering and the Ph.D. degree in microwave electronics from Chalmers University of Technology, Göteborg, Sweden, in 1994 and 1999, respectively. In 1997 he was a Research Assistant at the University of Virginia, Charlottesville, USA. From 1999 to 2001, he was a Research Fellow in the millimetre wave group at the Rutherford Appleton Laboratory, UK. He then joined Saab Combitech Systems AB as a Senior RF/microwave Engineer until 2003. From 2000 to 2006, he held different academic positions at Chalmers and was also Head of the Nanofabrication Laboratory at MC2 between 2003 and 2006. During the summer 2007, he was a Visiting Professor in the Submillimeter Wave Advanced Technology (SWAT) group at Caltech/JPL, Pasadena, USA. He is currently Professor and Head of the Terahertz and Millimetre Wave Laboratory at the department of Microtechnology and Nanoscience (MC2), Chalmers, Gteborg, Sweden. His research involves sources and detectors for terahertz frequencies, high frequency semiconductor devices, graphene electronics, terahertz measurement techniques and applications. He is also co-founder of Wasa Millimeter Wave AB. Prof. Stake serves as Topical Editor for the IEEE Transactions on Terahertz Science ad Technology.