1. Introduction

Three-dimensional integrated circuits (3D ICs) are generally considered to be one of the most prominent solutions that offer a way beyond Moore's law. A through-silicon via (TSV) is a popular choice for the vertical signal connection in 3D ICs and it is compatible with current process technologies. Some IC products,such as CMOS image sensors,have been manufactured using the TSV based 3D IC technology. The 3D IC technology will be applied to more commercial products.

The accurate extraction of resistance and capacitance (RC) parasitics of TSVs is essential in studying the performance and power consumption of 3D ICs. It is important for both the design-time circuit optimization,and the performance verification guaranteeing a reasonable production yield. Because TSVs across two IC tiers penetrate the silicon substrate,there is a kind of cylindrical metal-oxide-semiconductor (MOS) capacitor around the TSV,in addition to the conventional electrostatic capacitances among metal wires. As shown in Figure 1,$C_{\rm TSV}$ between a TSV and its surrounding silicon stands for the MOS capacitance,while other capacitance components belong to the conventional electrostatic capacitance. To accurately model the parasitic effects of TSV,all these capacitance components should be considered. This requires a comprehensive treatment to the semiconductor effect and electrostatic field. Field-solver simulators,such as the Sdevice of Synopsys^[1],is able to perform a complete extraction of the TSV structure,and produce accurate capacitance results. However,this kind of field solver employs an electromagnetic finite element method (FEM),and costs huge runtime and memory usage.

In Reference ^[2],the RLC parameters of the TSV were modeled as functions of physical parameters and metal characteristics. However,only one TSV was considered for modeling the semiconductor effect in MOS capacitance,with effects from surrounding TSVs ignored. In Reference ^[3],an equivalent circuit model was built for a 2-TSV layout structure. The method cannot be extended to the scenario with multiple TSVs. A model considering multiple TSVs was proposed based on the multi-conductor transmission line theory in Reference ^[4]. However,the silicon depletion region around the TSV was ignored in that work. In References ^{[5, 6, 7, 8]},the focus is on the electrostatic capacitance among ITVs and horizontal wires,instead of the complete RC parasitics including the MOS capacitance. We see that the existing works either consider the MOS capacitance and electrostatic capacitance separately,or are not suitable for arbitrary circuit layout with multiple TSVs and metal wires. The only possible tool for capturing both semiconductor and electrostatic capacitances is the simulator-like Sdevice^[1],which solves 3D semiconductor and electromagnetic Poisson equations. However,for a general TSV layout in 3D IC,this kind of simulator suffers from the huge computing expense.

In this work,we combine both extraction techniques for MOS capacitance and for electrostatic capacitance,establishing a comprehensive modeling and extraction framework for the general TSV layout. This framework includes a RC equivalent circuit suitable for low and medium signal frequencies,the analytical extraction technique for MOS capacitance,and a random walk based extraction technique for electrostatic capacitance. The random walk based technique is advantageous over the conventional deterministic methods for capacitance extraction,like the FEM and boundary element method,because it does not rely on generating and solving a matrix equation. Compared with the deterministic methods,the floating random walk method is very suitable for large-scale structures^{[8, 10]}. The proposed hybrid method exhibits good modeling accuracy at 1 GHz signal frequency,and demonstrates up to 47$\times$ speedup over the Sdevice based simulation method for multi-TSV test cases.

2. Background

2.1 The MOS capacitance of a single TSV

In this subsection,we discuss the MOS capacitance of a TSV in the situation where only a single TSV exists in the circuit layout^[2]. The TSV serves as the interconnection between the bottom tier and the top tier,as shown in Figure 2. The copper metal,silicon oxide liner and the silicon substrate form a cylindrical MOS structure. With given voltage imposed on the TSV,there is a depletion region in the silicon around the TSV. Therefore,the MOS capacitance includes the contributions from the oxide linear region and the depletion region. Among them,the liner capacitance A structure with one TSV. (a) The cross-section view. (b) The top view. Fig.2

$\begin{equation} C_{\rm ox}=\frac{2\pi \varepsilon_{\rm ox}l_{\rm TSV}}{\ln \dfrac{R_{\rm ox}}{R_{\rm metal}}}, \end{equation}$

(1)

where $R_{\rm metal}$,$R_{\rm ox}$,$l_{\rm TSV}$,and $\varepsilon_{\rm ox}$ are the metal radius,the oxide liner radius,the length of the TSV,and the permittivity of oxide,respectively (see Figure 2(b)).

To get the radius for the depletion region ($R_{\rm dep})$ and the capacitance $C_{\rm dep}$ we have to solve a 1-D Poisson equation in a cylindrical coordinate system^[2]. Firstly,we obtain the formula for the Si-SiO$_{2}$ surface potential ($\phi_{\rm s})$:

$\begin{equation} \phi_{\rm s}=\frac{qN_{\rm a}(R_{\rm dep}^2-R_{\rm ox}^2)}{4\varepsilon_{\rm si}}- \frac{qN_{\rm a}R_{\rm dep}^2}{2\varepsilon_{\rm si}} \ln \frac{R_{\rm ox}}{R_{\rm dep}}, \end{equation}$

(2)

where $N_{\rm a}$ is the doping concentration of the P-Si substrate and $\varepsilon_{\rm si}$ is the permittivity of silicon. $R_{\rm dep}$ varies with the change of voltage imposed on the TSV ($V_{\rm TSV})$. With the TSV bias voltages from low to high,there are three distinct regions of operation in the silicon substrate: accumulation,depletion,and inversion. The flat-band voltage ($V_{\rm fb})$ of the TSV is defined as the TSV voltage at which $\phi_{\rm s}$ equals to zero. At that situation,$R_{\rm dep}$ equals to $R_{\rm ox}$,which can be derived from Equation (2). When the TSV voltage becomes larger than $V_{\rm fb}$,the operation region is converted from accumulation to depletion. The flat-band voltage can be calculated with^[9]:

$\begin{equation} V_{\rm fb}=\phi_{\rm ms}-\frac{Q_{\rm s}}{C_{\rm ox}}, \end{equation}$

(3)

where $\phi_{\rm ms}$ is the work function difference of copper and silicon^[9],and $Q_{\rm s}$ is the charges inside the oxide liner.

$\begin{equation} Q_{\rm s}=2\pi R_{\rm ox}Q_{\rm ot}l_{\rm TSV}, \end{equation}$

(4)

where $Q_{\rm ot}$ is the oxide charge density^[9]. As $V_{\rm TSV}$ increases beyond a threshold voltage ($V_{\rm th})$,the MOS capacitance transfers from the depletion operation region to the inversion operation region. $V_{\rm th}$ is defined as the TSV voltage at which $\phi_{\rm s}$ equals to the double of Fermi potential ($\phi_{\rm Bp}$)^[9].

$\begin{equation} \phi_{\rm Bp}=\frac{kT}{q}\ln \frac{N_{\rm a}}{n_{\rm i}}, \end{equation}$

(5)

where $T$ and $k$ represent the temperature and Boltzmann constant,respectively. $n_{\rm i}$ is the intrinsic carrier concentration,and $q$ is the electrical charge of an electron. With Equations (2) and (5),we have:

$\begin{equation} \frac{qN_{\rm a}(R^2_{\rm dep}-R_{\rm ox}^2)}{4\varepsilon_{\rm si}}-\frac{qN_{\rm a}R^2_{\rm dep}}{2\varepsilon_{\rm si}} \ln \frac{R_{\rm ox}}{R_{\rm dep}}-2\frac{kT}{q}\ln \frac{N_{\rm a}}{n_{\rm i}}=0. \end{equation}$

(6)

The value of $R_{\rm dep}$ corresponding to $V_{\rm th}$ can be solved with Equation (6),and it is denoted by $R_{\rm max}$. Then,$V_{\rm th}$ can be obtained.

$\begin{equation} V_{\rm th}=V_{\rm fb}+2\phi_{\rm Bp}+\frac{qN_{\rm a}(R^2_{\rm max}-R^2_{\rm ox})}{2 \pi \varepsilon_{\rm ox}} \ln \frac{R_{\rm ox}}{R_{\rm metal}}. \end{equation}$

(7)

For most ICs,the signal voltage on the TSV is between 0 and VDD. It is easy to see from Equation (3) that the flat band voltage $V_{\rm fb}$ is negative. Thus,the whole accumulation operation region and part of the depletion operation region do not need to be taken into account. Within the working voltage region,the TSV MOS capacitance is the series combination of the oxide and depletion capacitances,calculated by:

$\begin{equation} C_{\rm TSV}=\frac{C_{\rm ox}C_{\rm dep}}{C_{\rm ox}+C_{\rm dep}}, \end{equation}$

(8)

where

$\begin{equation} C_{\rm dep}=\frac{2\pi \varepsilon_{\rm si}l_{\rm TSV}}{\ln \dfrac{R_{\rm dep}}{R_{\rm ox}}} \end{equation}$

(9)

As long as $V_{\rm TSV}$ belongs to the range of the depletion operation region,$R_{\rm dep}$ can be solved from

$\begin{equation} V_{\rm TSV}=V_{\rm fb}+\phi_{\rm s}+\frac{qN_{\rm a}(R^2_{\rm dep}-R^2_{\rm ox})}{2 \pi \varepsilon_{\rm ox}} \ln \frac{R_{\rm ox}}{R_{\rm metal}}. \end{equation}$

(10)

The depletion capacitance $C_{\rm dep}$ becomes the minimum at the inversion operation region,because $R_{\rm dep}$ reaches its maximum value $R_{\rm max}$. This minimum depletion capacitance ($C_{\rm dep,\,min})$ is:

$\begin{equation} C_{\rm dep,\,min}=\frac{2\pi \varepsilon_{\rm si}l_{\rm TSV}}{\ln \dfrac{R_{\rm max}}{R_{\rm ox}}}, \end{equation}$

(11)

The above derivation shows that the MOS capacitance is depending on the voltage on the TSV. This is a distinct difference from the conventional electrostatic capacitance. To validate the accuracy of the above formulas,simulation experiments are carried out with the Sdevice whose results are compared with those from Equations (1)-(11). For a structure with a copper TSV embedded in a P-Si silicon layer,we set the geometry parameters and the doping concentration as those listed in Table 1.

Considering $V_{\rm TSV}$ ranging from 0 to 3 V,the MOS capacitance of the TSV is calculated from Equations (1)-(11) and simulated with the Sdevice. For this case,$V_{\rm fb}$ equals to $-0.24$~V according to Equation (3),$V_{\rm th}$ equals to 0.97 V according to Equations~(3),(6),and (7). $C_{\rm TSV}$ remains $C_{\rm dep,\,min}$ and equals to 3.75 $\times$ 10$^{-14}$ F for the inversion operation region according to Equations (6),(8) and (11). $C_{\rm TSV}$ is calculated from Equations (3),(8)-(10) for the depletion operation region. A comparison of the analytical method and Sdevice results is plotted in Figure 3. It shows that the maximum discrepancy between both methods is less than 2.1{\%}. This verifies the accuracy of the analytical model for calculating the MOS capacitance of a single cylindrical TSV.

2.2 Electrostatic coupling capacitances among TSVs

The electrostatic coupling capacitances among TSVs and conventional wires may be comparable to the MOS capacitance and cannot be ignored^{[5, 7, 8]}. There are many capacitance extraction techniques for conventional interconnect wires^[11]. However,for the cylindrical TSVs especially in the context of general layout,specific field solver techniques are necessary to guarantee a satisfying accuracy^{[8, 12]}.

Recently,an efficient floating random walk (FRW) based approach was proposed to extract the electrostatic capacitances of TSVs^[8]. The approach produces highly accurate results,and is ten times faster than existing capacitance solvers. The basic ideas in this approach are briefly introduced as follows.

The fundamental formula of the FRW algorithm is

$\begin{equation} \phi(r)=\oint_s P(r,r^{(1)})\phi(r^{(1)}){\rm d}r^{(1)}, \end{equation}$

(12)

where $\phi(r)$ is the electric potential at point $r$,and $S$ is a closed surface surrounding $r$. $P(r,r^{(1)})$ is called the surface Green's function. If $\phi(r^{(1)})$ is unknown,we apply Equation (12) recursively to obtain the following nested formula:

$\begin{equation} \begin{split} {} \phi(r)= \oint_{s(1)}P^{(1)}(r,r^{(1)})\oint_{s(2)}P^{(2)}(r^{(1)},r^{(2)})\cdots \\[2mm] \oint_{s(k+1)}P^{(k+1)}(r^{(k)},r^{(k+1)})\phi(r^{(k+1)}){\rm d}r^{(k+1)}\cdots {\rm d}r^{(2)}{\rm d}r^{(1)}, \end{split} \end{equation}$

(13)

where $S_{(i)}$ ($i =$ 1,$\cdots$,$k+1)$ is the surface of the $i$-th cube centered at $r^{(i-1)}$. $P^{(i)}(r^{(i-1)}$,$r^{(i)})$,($i=$ 1,$\cdots$,$k+1)$,are the surface Green's functions relating to the potentials at $r^{(i-1)}$ to $r^{(i)}$. This can be interpreted as a floating random walk procedure,as shown in Figure 4. To calculate the coupling capacitance between conductor blocks,the charge of the conductor is related with the electric potential through the Gauss theorem.

Therefore,with the random walks starting from the Gaussian surface of a specified conductor (master),we can calculate the coupling capacitance between the master and other conductors (see Figure 4)^[10]. The FRW based capacitance extraction is a discretization-free method,and thus enjoys the advantages of better scalability for large structures,tunable accuracy,and higher parallelism,etc.

In Reference ^[8],the FRW algorithm was extended to handle the cylindrical TSV structures. The idea is to allow the transition cube which is suitable for the conventional Manhattan shape of interconnects to rotate an angle for better touching the TSV surface (see Figure 5). This increases the probability of terminating the random walk quickly,and therefore improves the computational efficiency. Also,a special space management technique was devised to accelerate the calculation of size of maximum FRW hop for a general,large layout including wires and TSVs.

3. A comprehensive parasitic model and extraction methodology for TSVs

The purpose of this work is to combine the existing works on MOS capacitance calculation (Section 2.1) and the electrostatic capacitance extraction (Section 2.2) for an efficient and comprehensive modeling and extraction framework for the general TSV layout. Our focus is to derive the equivalent circuit for analyzing the signal integrity on a single TSV net (which we call “victim TSV”) and to propose an efficient extraction technique. Because the signal frequency of ICs is not very high in most scenarios,we only consider the low and medium frequency range. Therefore,in this work we do not consider the inductance effect among the TSVs.

Due to the lossy nature of the silicon substrate,the substrate resistances need to be considered as well. For a homogeneous silicon substrate,the resistance element is in parallel to the capacitance element,and can be calculated with:

$\begin{equation} R_{\rm si}=\frac{\varepsilon_{\rm si}}{\sigma_{\rm si}C_{\rm si}}, \end{equation}$

(14)

where $\sigma_{\rm si}$ is the electrical conductivity of silicon. Therefore,we will focus our discussion on the capacitance elements modeling bulk silicon effect,MOS effect and other electrical coupling effects,and the resulted RC equivalent circuit model.

3.1 Two-TSV structure

We first consider a simple layout including only two TSVs,and present the equivalent RC circuit model. As shown in Figure 6,the impact from other TSV (“aggressor”) on the victim TSV is considered. Because of the size of the TSV,its DC resistance is negligible. Therefore in Figure 6,each TSV corresponds to a circuit node,and R and C elements reflect the couplings among two TSVs and the grounded substrate. In the equivalent circuit,$C_{\rm tsv}$'s and $C_{\rm si}$'s represent the MOS capacitances and the electrostatic capacitances in the silicon substrate,respectively. The parallel $R_{\rm si}$'s models the conduction effect of silicon and can be calculated with Equation (14). $C_{\rm tsv}$'s can be calculated with the techniques in Section 2.1,while $C_{\rm si}$'s can be extracted with the FRW based field solver^[8].

This circuit model depicts both the semiconductor effect and electrostatic coupling. With it,the lump-port elements of the circuit can be easily derived. For example,the total admittance between TSV {\#}1 and other conductors in Figure 6 can be derived:

$\begin{equation} Y_1=\frac{{\rm j}\omega C_{\rm tsv,\,1}Y_3}{{\rm j}\omega C_{\rm tsv,\,1}+Y_3}, \end{equation}$

(15)

where $\omega$ is the angular frequency of the signal,and

$\begin{equation} Y_3={\rm j}\omega C_{\rm si,\,g1}+\frac{1}{R_{\rm si,\,g1}}+\frac{\left({\rm j}\omega C_{\rm si,\,12}+\dfrac{1}{R_{\rm si,\,12}}\right)Y_4}{{\rm j}\omega C_{\rm si,\,12}+\dfrac{1}{R_{\rm si,\,12}}+Y_4}, \end{equation}$

(16)

$\begin{equation} Y_4={\rm j}\omega C_{\rm tsv,\,2}+{\rm j}\omega C_{\rm si,\,g2}+\frac{1}{R_{\rm si,\,g2}}. \end{equation}$

(17)

$Y_{3}$ and $Y_{4}$ stand for the admittances of node 3 and node 4 to the environment,respectively. Here we assume TSV {\#}2 is grounded. Then,the equivalent total capacitance of TSV {\#}1 is:

$\begin{equation} C_1=\frac{Y_1}{{\rm j} \omega}. \end{equation}$

(18)

From Equations (15)-(18) we see that at low and medium frequency,the equivalent capacitance of the TSV is frequency-dependent,and because the formulas involve the MOS capacitance,the capacitance also depends on the voltage imposed on the victim TSV.

3.2 Multi-TSV structure

The RC equivalent circuit model for the two-TSV structure can be extended for a structure including more TSVs. In Figure 7,we show a structure with 5 TSVs and its corresponding equivalent circuit. In this structure,we have one victim TSV and 4 aggressor TSVs. The symbol convention of the circuit elements is the same as that in the two-TSV example. Because our purpose is to model the signal integrity on the victim TSV,the electrostatic couplings among the aggressor TSVs are ignored in this equivalent circuit. As shown in experiment results,this simplification hardly harms the modeling accuracy. In the equivalent circuit,$C_{\rm tsv}$'s can be calculated with the analytical approach,while $C_{\rm si}$'s can be extracted with the FRW based solver. Because the method for extracting the electrostatic capacitances is very general,the layout of TSVs can be arbitrary. So,if the horizontal wires are present in the structure,we just need to add some electrostatic capacitances (extracted by the FRW based approach) to the circuit for comprehensive electrostatic/semiconductor modeling.

3.3 Extraction algorithm flow

According to the above discussion,we can easily obtain a RC equivalent circuit from the layout of the TSVs in 3D IC. As for the extraction of RC elements in the circuit,the techniques in Section 2 are employed. The algorithm for extracting the RC elements is given as follows.

With the RC elements extracted,we obtain the equivalent circuit like in Figure 7(b). Then,the lump port parameters can be calculated. For example,the total frequency-dependent capacitance can be calculated with Algorithm 2.

In Algorithm 2,we assume that the aggressor TSVs are grounded. If this is not the case,the total lump capacitance can also be derived by solving the circuit equation.

In Algorithm 1,most computational time is used to perform the FRW based electrostatic extraction. Due to the merits of the FRW based algorithm^[8],the RC extraction algorithm is versatile,and more efficient than the approach performing electrostatic/semiconductor simulation for the whole structure. In the experiments,we will demonstrate the benefits of the proposed method.

4. Numerical results

We have used the FRW program for electrostatic capacitance extraction and written a Matlab program for MOS capacitance calculation. For the structures including two,five and nine TSVs,we have tested our algorithms. The Sdevice simulator,which employs FEM for overall electrostatic/semiconductor simulation,is also run to validate the accuracy and efficiency of the proposed method. All experiments are carried out on a Linux server with an Intel Xeon E5-263 2~GHz CPU,and 32 GB memory. The accuracy criterion of the FRW extraction algorithm^[8] is set to 1{\%} 1-$\sigma$ error. In all test cases,the silicon properties are the same as those in the experiment in Section 2.1. The temperature is the room temperature.

4.1 2-TSV case

The 2-TSV structure in Section 3.1 is tested to verify the accuracy of the proposed techniques. The DC bias voltage is zero for the aggressor TSV. For each TSV,the diameter is 2.5~$\mu$m,the thickness of the oxide linear is 118.2 nm,and the length is 20 $\mu$m. The pitch is 20 $\mu$m. The electrostatic capacitance between TSVs ($C_{\rm si,\,12})$ is 2.4 fF and the electrostatic capacitances between TSVs and the substrate ($C_{\rm si,\,gj})$ are both 1.48 fF. The MOS capacitance of a TSV is shown in Figure~3.

A comparison of the analytical method and Sdevice results is plotted in Figure 8,while Sdevice simulation with high-density grid is considered accurate. From the figure,we find that the maximum discrepancy between the both technologies is 4.9{\%} for signal frequencies from 10 kHz to 500 MHz when the victim TSV bias voltage equals to $V_{\rm th}$. The discrepancy of both results is within 3{\%},when the victim TSV voltage varies from 0 to 3 V.

4.2 5-TSV case

For the 5-TSV structure in Figure 7,the parameters of each TSV and the pitch between the victim and the aggressor are the same as those in the 2-TSV case. The electrostatic capacitance between the victim and aggressor ($C_{\rm si,\,1j})$ is 1.36 fF,between the victim and substrate ($C_{\rm si,\,g1})$ it is 0.316 fF and between the aggressor and substrate ($C_{\rm si,\,gj}$,except for $C_{\rm si,\,g1})$ it is 0.851 fF. A comparison of the analytical method and Sdevice results is plotted in Figure 9. From the figure we see that the maximum discrepancy between both technologies is 4.6{\%} at the signal frequency from 10 kHz to 40 MHz when victim TSV bias voltage equals to $V_{\rm th}$. The maximum discrepancy at 500 MHz signal frequency is 5.1{\%},while the victim TSV voltage changes. The results show the model matches well with the Sdevice.

4.3. 9-TSV case

For further verification of the proposed method,we consider a 9-TSV case where eight aggressor TSVs form a square and one victim TSV exists in the center. The electrostatic capacitances between the substrate and the victim,the aggressor in the vertex and the midpoint of sides of the square are 0.08,0.657 and 0.244 fF,respectively. The electrostatic capacitances between the victim and the aggressor in the vertex and the midpoint of sides of the square are 0.383 and 1.08 fF. The comparative results of the difference are shown in Figure 10. From the figure we see that the discrepancy between both technologies is within 3{\%} at the signal frequency from 10 kHz to 800~MHz when victim TSV bias voltage equals to $V_{\rm th}$. The maximum discrepancy at 800 MHz signal frequency is 4.8{\%},while the victim TSV voltage changes. According to these structures,the proposed method is verified with the accepted loss of accuracy. For this larger case,the accuracy of the proposed method is higher.

In previous experiments,very dense discretization grids are employed in the Sdevice to guarantee that it produces the most accurate results. To compare the runtime of the proposed method fairly,we reduce the grid numbers in the Sedevice ({\#}grid) to make it produce results with a similar accuracy to the proposed method. Under 100 MHz signal frequency,the three structures are extracted with the Sdevice and the proposed method,respectively. The runtime comparisons are given in Table 2. For the 2-TSV case,the proposed method costs more time since the FRW algorithm is not beneficial for a small case. For cases with more TSVs,the runtime of the FRW is reduced. Comparing with the 2-TSV case,the runtime of the proposed method is reduced evidently for 5-TSV and 9-TSV cases. This is because more conductors around the master conductor actually make the random walks terminate earlier. On the contrary,the runtime of the Sdevice increases with the increase of TSV number because the discretization grids involved quickly increase. We see that the speedup ratio of the proposed method over Sdevice simulation increases with the number of TSVs,and reaches 47$\times$ for the 9-TSV case. With these cases,the acceleration of the proposed method is verified for a multi-TSV structure.

5. Conclusions

In this paper,we have explored the comprehensive parasitic RC model and the extraction method for TSVs in 3D IC. Combining the analytical method for calculating the MOS capacitance and the FRW based extraction for electrostatic coupling capacitances,an efficient modeling and extraction framework is proposed. For several TSV structures,the total coupling capacitance on a victim TSV is extracted with the proposed framework and the Sdevice simulator. Within the signal frequencies from 10 kHz to 1 GHz and for different TSV voltages,the proposed method produces capacitance with an error of no more than 5{\%}. As for the computational speed,the proposed method is up to 47$\times$ faster than the Sdevice simulation based method. The results show that the proposed method has good accuracy,and is very efficient especially for large structures.