System and method for statistical modeling and statistical timing analysis of integrated circuits

[0001] This invention relates to the field of computer-aided design (CAD) tools for integrated circuits (ICs). More specifically, the invention relates to using statistical methods for delay modeling and static timing analysis of integrated circuits.

[0002] FIG. 1 shows a typical prior art flow (100) for nominal modeling and prediction of the timing properties of a digital integrated circuit or integrated circuit macro (portion of an integrated circuit). Box 110 shows the manufacturing process which is characterized by a set of process parameters. These parameters are well-known and varied. A parameter is any physical characteristic of a manufacturing process that determines the behavior of the product being manufactured. Non-limiting examples are temperature of the oven, junction depth and concentration of doping impurities.

[0003] These parameters are used to create a device model as shown in Box 120, which captures the electrical behavior of individual devices (such as transistors) and wires. Gates, such as a 2-input NAND gate, are comprised of multiple transistors and wires. The device models of Box 120 are used in Box 130 to create timing models of gates. Device models and timing models are also well-known. Timing models of individual gates often have timing attributes like rising delay, falling delay, rising slew and falling slew associated with them (see IEEE standard for integrated circuit (IC) delay and power calculation system, IEEE Standard 1481-1999, 1999, pages 1-390, available at hffp://fp.ieeexplore.ieee.org/ie15/6837/18380/00846710.pdfisNumber18380&prodstandards& arnumber00846710).

[0004] Finally, the gate timing models are used to predict timing properties of an entire chip or macro of the chip by means of static timing analysis (see R. B. Hitchcock, Sr., G. L. Smith and D. D. Cheng, Timing analysis of computer hardware, IBM Journal of Research and Development pages 100-105, January, 1982), as shown in Box 140. While FIG. 1 depicts the most common steps in modeling processes, devices, gates and integrated circuits, it is to be understood that this is a non-limiting modeling method.

[0005] A typical gate delay model is explained in reference to FIG. 2. FIG. 2(a) depicts a particular transition of a 3-input NAND gate in which input A rises and output Y falls, while B and C are held at logic 1. The NAND gate example has an input slew rate, s_in, for input A and an output a capacitance, C_out. FIG. 2(b) shows the waveforms of A and Y. The difference between the 50% crossing times of A and Y is the delay of the gate, d_out. The time taken by Y to make its transition is the output slew of the gate, s^out. The gate delay model comprises a model for the delay and slew for every valid input-to-output transition (also called pin-to-pin transition). Each such transition will be represented by an edge of the timing graph during static timing analysis. FIG. 2(c) shows the set of edges contributed by the 3-input NAND gate to the timing graph. Note that when A is switching, B and C are held at logic 1 because otherwise the output of the gate would not switch due to the logical behavior of the NAND gate. Repeated circuit simulations of the gate being modeled are performed a priori to measure the output delay and output slew for each possible input-to-output transition. Then the measurements corresponding to each transition are typically fitted to analytic delay and slew models by adjusting the coefficients of the models to minimize the error between the measured delays/slews and those predicted by the model. In this non-limiting example, a typical delay model is shown below:

[0006] where a_1,2,3,4 and b_1,2,3,4 are coefficients of the delay model, C_out is the output capacitance of the gate, s_in is the slew of the input transition, d_out is the delay of the gate and s_out is the output slew of the gate. The delay/slew models are typically functions of input slew and output capacitive load (the model arguments) as shown in the equations above. These modeling methods are well known. For example, see IEEE standard for integrated circuit (IC) delay and power calculation system, IEEE Standard 1481-1999, 1999, pages 1-390, available at http://fp.ieeexplore.ieee.org/iel5/6837/18380/00846710.pdfisNumber18380&prodstandards& arnumber00846710 which is herein incorporated by reference in its entirety.

[0007] The gate delay models obtained by the above known procedures are used during static timing analysis to predict the timing properties of an integrated circuit or macro of an integrated circuit. The delay of the critical path of the circuit in late-mode timing analysis dictates the fastest clock frequency at which the circuit can be operated. Since static timing analysis predicts the delay of the critical path, it gives us a method of determining the fastest clock frequency at which the circuit will safely work.

[0008] The procedure outlined above works very well to predict the nominal performance of an integrated circuit, i.e., its performance in the absence of variations. Unfortunately, variations are inevitable in any practical situation. There are variations both due to imperfect control of the manufacturing process and due to environmental variations (such as temperature and power-supply voltage). Some of these variations can cause catastrophic defects, resulting in an unrecoverable failure of the circuit. Others cause the timing performance of the circuit to vary, and these variations are called parametric variations. Parametric variations cause manufactured chips to have a range of critical path delays, and therefore to be operable at a range of maximum clock frequencies. Some integrated circuits such as microprocessors are sorted into bins depending on their working frequency; chips in the high-frequency bin command the highest price and hence the highest profit. These are called sorted designs. Other chips are required to work at a certain frequency and chips below that performance level are useless. These are called non-sorted designs.

[0009] In any case, predicting the parametric variations in the critical path delay is very important. The critical path delay determines the fastest frequency of operation, and its probability distribution in turn determines the yield and profitability of an integrated circuit. Analyzing the parametric variations is a key step in estimating the benefit of design changes to improve the yield (percentage of working parts) of an integrated circuit. Design for Manufacturability implies changing the design so as to obtain the highest number of working chips at the highest possible performance so as to maximize profitability.

[0010] Variations due to manufacturing are numerous and complex. The variations stem from the lack of perfect control of machines in the fabrication line. Dust contamination in the clean room, mechanical tolerances, lens aberrations, inhomogeneous gas flows, variations in gas pressure, insufficient control of temperature, deposition failures with spinning and lack of perfect surface planarity with polishing are sources of parameter variations. A non-exhaustive set of 35 phenomena giving rise to manufacturing variations is listed on page 12 of Chapter 2 of Simon, P., Yield Modeling for Deep Sub-Micron IC Design, Ph.D. dissertation, ISBN 90-75341-28-8, Eindhoven University of Technology, Eindhoven, The Netherlands, December 2001, which is herein incorporated by reference in its entirety.

[0011] Variations can be lot-to-lot, wafer-to-wafer, die-to-die and intra-die. The variations can also be functions of the position of a chip on the wafer, or the position of a gate on a die, for example. These variations exhibit complex correlations. The variability in these parameters is clearly stochastic in nature.

[0012] One method of predicting the effect of parameter variations is by Monte Carlo analysis in which the space of all possible variations is sampled and a complete analysis of the type shown in FIG. 1 (flow 100) is conducted at each sample point. With a sufficient number of such samples, the distribution of chip performance can be predicted with a measurable confidence. To reduce the number of analyses required, the variations are typically analyzed by so-called best case, nominal and worst-case analyses. The processing and environmental conditions are all assumed to be correlated to produce the best (fastest) possible circuits for the best-case analysis. Likewise, worst-case assumptions are made to predict the slowest performance of the circuit. In each analysis, in reference to FIG. 1, a new set of device models 120 is created, a new set of gate delay models 130 is created and a new static timing analysis 140 is conducted. Thus multiple analyses are conducted to predict the distribution of the performance of the final manufactured product.

[0013] Prior art methods such as the Monte Carlo analysis and the best/nominal/worst-case analyses mentioned above suffer from several weaknesses. First, they either ignore or restrict the nature of correlations between the sources of variation. Consequently, the results they predict suffer from inaccuracy, and are often needlessly pessimistic. Second, they can suffer from extremely long computer run times. Third, many methods assume that the delay of each gate or wire is available as a probability distribution function with arbitrary correlations between the various delays, but do not provide an efficient means for deriving such models. Fourth, prior art methods often compute discrete points on the probability distribution of the performance of the circuit as opposed to computing the entire probability distribution. The prior art fails to gain knowledge of the entire distribution that would be advantageous in design for manufacturability, and yield/profit modeling and optimization of the integrated circuit or integrated circuit macro. Finally, prior art methods do not provide a full methodology or flow for statistical modeling and timing of integrated circuits, going from process variations all the way to predicting variations in the performance of integrated circuits.

[0014] An object of this invention is an improved system and method for modeling and predicting the timing of (integrated) circuits in the face of manufacturing and environmental variations.

[0015] It is another object of this invention to provide a system and method for creating statistical timing models of individual gates in a gate library.

[0016] It is still another object of this invention to compute efficiently the sensitivity of the coefficients of analytic delay models with respect to sources of parameter variations.

[0017] It is yet another object to use timing models of individual gates along with relevant sensitivity information to carry out a statistical timing analysis of an integrated circuit or macro of an integrated circuit.

[0018] It is a further object to compute the probability distribution of the slack of the critical path of an integrated circuit or integrated circuit macro.

[0019] It is a more particular object to predict the yield of an integrated circuit or integrated circuit macro at a given performance specification.

[0020] It is another particular object to predict the profit from an integrated circuit or integrated circuit macro.

[0021] It is yet another particular object to predict the maximum performance of an integrated circuit or integrated circuit macro given a required yield.

[0022] It is a further particular object to predict the impact of a change in the design or manufacturing process or equipment on the statistical performance, yield and profitability of a circuit.

[0023] This invention describes a system, method, and program product that provide a comprehensive methodology for statistical modeling and timing of integrated circuits and integrated circuit macros. More specifically, it provides a means for efficiently computing the sensitivities of gate delay models (and their coefficients) to sources of parameter variations. These sensitivities are used to determine the probability distribution of the delay and slew of each gate and wire, as well as the correlations between these delays and slews. In some preferred embodiments, these timing models are used in an inventive statistical static timing analysis method to predict the statistical performance of an integrated circuit or integrated circuit macro.

[0024] The foregoing and other objects, aspects, and advantages will be better understood from the following non-limiting detailed description of preferred embodiments of the invention with reference to the drawings that include the following:

[0025] FIG. 1 is a process flow diagram showing typical prior art methodology of nominal modeling and static timing of an integrated circuit or integrated circuit macro.

[0026] FIG. 2(a) is a non-limiting, prior art example diagram of a single transition of a sample 3-input NAND gate.

[0027] FIG. 2(b) shows the waveforms of an example single transition and the definition of output delay and slew for the simple 3-input NAND gate of FIG. 2(a).

[0028] FIG. 2(c) is a diagram showing the set of edges contributed by the sample gate of FIG. 2(a) to the example timing graph of a circuit in which this gate occurs.

[0029] FIG. 3 is a process flow diagram showing one preferred embodiment of the inventive methodology for statistical modeling and static timing of a circuit, e.g., an integrated circuit or integrated circuit macro.

[0030] FIG. 4 is a flow chart of a preferred process of determining the sensitivity of gate timing model coefficients of a single transition of a gate with respect to all sources of variation.

[0031] FIG. 5 is a flow chart of a preferred embodiment of a novel statistical static timing procedure.

[0032] While the present invention can be applied to any general circuit analysis, the preferred embodiment of the invention involves a method of creating statistical timing models for electrical circuit gates such as NAND and NOR gates, and more preferably then using these statistical timing models in a process of statistical static timing analysis to predict the statistical timing properties of an integrated circuit or integrated circuit macro. Thus this disclosure will describe the use of the invention in integrated circuit design without loss of generality.

[0033] FIG. 3 shows a process flow diagram for statistical modeling and static timing 300. It is assumed that statistical process models (box 310), and transistor models (box 320) are available from well-known prior-art techniques (see, for example, C. Yu, T. Maung, C. J. Spanos, D. S. Boning, J. E. Chung, H. -Y. Liu, K. -J. Chang and D. J. Bartelink, Use of short-loop electrical measurements for yield improvement, IEEE Transactions on Semiconductor Manufacturing, vol. 8, number 2, May 1995, pages 150-159). This invention focuses on the creation of gate delay models along with relevant sensitivity information (box 330 and FIG. 4, below) and statistical static timing analysis (box 340 and process 500, below).

[0034] We begin by describing the steps involved in creating gate delay models with relevant sensitivity information (box 330). These will later be employed to conduct statistical timing analysis (box 340 and process 500). Typically (as indicated by Box 130 of FIG. 1), the actual behavior of the circuit, in a preferred embodiment, the set of timing properties of all members of a gate library, is characterized by repeatedly running a circuit simulation program, using the steps described below. First, an analytical form is decided upon to model (or predict) each timing behavior, such as delay and output slew. The analytical form is usually a polynomial. The model arguments are typically input slew and output capacitance, but could also include ambient temperature and power-supply voltage, for example. The predictive model is a function of these model arguments and tunable coefficients. The model arguments are then sampled at a multiplicity of sample value points, and a circuit simulation program is employed at each sample value point to determine the actual timing behavior of the circuit. Curve fitting or regression fitting is then used on all the collected data. In the fitting procedure, the tunable coefficients are adjusted to minimize the difference between the actual behavior and predicted behavior over all of the sample value points. At the end of this procedure, each actual timing behavior can be modeled or predicted as an analytic function of model arguments with known, constant coefficients. The above procedure is well known and widely practiced.

[0035] In practice, timing properties are modeled by different methods, and the models themselves are p created and used in different ways, and the models are functions of different quantities, some of which are subject to variability and others not It is to be understood that the present inventive method can easily be applied to any of these existing methodologies by a practitioner of ordinary skill in the art.

[0036] The ensuing description presents a method of determining the timing models along with relevant sensitivity information depicted in box 330 of FIG. 3.

[0037] Assume that a timing property such as the delay of a pin-to-pin transition of a particular gate's predicted behavior is modeled by the following predictive model:

[0038] where d_pred is the predicted timing property (or generally a measure of predictive behavior), a_i are a set of m tunable coefficients of the model, s_in is the input slew, C_out is the output load (see FIG. 2(a)), and d_pred is an analytic, possibly nonlinear function, i.e., the timing model. For multiple, say n, preselected sample value points of input slew and output load (set of one or more model arguments), the gate is simulated with a circuit simulator to obtain a set of n measurements that we call d_actual.

[0039] The model d_pred is now tuned or adjusted (the tunable coefficients, a_i, are changed) to approximate the actual behavior, d_actual, that is obtained from the simulator (or other source of actual behavior information). To do this, any curve fitting method could be used. In the preferred embodiment, the process of fitting involves adjusting the coefficients a_i in order to minimize the least squared error. Specifically, the error $E=\underset{j=1}{\overset{n}{}}{\left[d_{\mathrm{actual}}^j-d_{\mathrm{pred}}^j\right]}^2$

[0040] where the superscript, j, indicates the measurement number, i.e., the index identifying the sample point. By minimizing this error, the coefficients of the gate timing model, d_pred, are obtained. Now, we are interested in finding the complete sensitivity of d_pred to all relevant sources of parameter variation, be they manufacturing variations or environmental variations.

[0041] Assume that all sources of parameter variation are collected in a vector P with variations denoted by P. Then the first-order change in the predicted delay d_pred can be expressed as $d_{\mathrm{pred}}={\left[\frac{d_{\mathrm{pred}}}{P}\right]}^{T}P$

[0042] where the complete differential can be expanded as $d_{\mathrm{pred}}={\left[\left\{\underset{i=1}{\overset{m}{}}\frac{d_{\mathrm{pred}}}{a_i}\frac{a_i}{P}\right\}+\frac{d_{\mathrm{pred}}}{s_{in}}\frac{s_{in}}{P}+\frac{d_{\mathrm{pred}}}{C_{\mathrm{out}}}\frac{C_{\mathrm{out}}}{P}\right]}^{T}P.$

[0043] The term in the square brackets above is the total sensitivity of d_pred to P. Obviously, to compute this total sensitivity, we require a_i/P, i1,2, . . . , m, which are the sensitivities of the tunable coefficients of the delay model to the sources of variation. One inventive method of determining and applying these sensitivities is described below. Note that computing these sensitivities by perturbation or finite differencing is both prohibitively expensive and numerically dangerous. Each source of variation would have to be perturbed sequentially, and the measuring and fitting procedure repeated for each source of variation. Numerical problems could be encountered because of the noisiness of the circuit simulator and because of finite-precision computer arithmetic. The amount of the perturbation is also crucial. If the perturbation is too small, numerical problems can result because of division by small quantities. If the perturbation is too large, we are computing a poor approximation to the sensitivity. Thus we would like to determine the required sensitivities without resorting to perturbation or finite differencing.

[0044] Differentiating the above expression for E with respect to a particular coefficient a_i, we obtain $\frac{E}{a_i}=\underset{j=1}{\overset{n}{}}2\left[d_{\mathrm{actual}}^j-d_{\mathrm{pred}}^j\right]\frac{d_{\mathrm{pred}}^j}{a_i}=0.$

[0045] The fitting procedure ensures that E is at a minimum with respect to each a_i and hence the gradient of E with respect to each a_i is zero.

[0046] We now differentiate the above equation with respect to a parameter p, a source of variation (such as the effective length of a transistor), to obtain $\underset{j=1}{\overset{n}{}}2\left[\frac{d_{\mathrm{actual}}^j}{p}-\frac{d_{\mathrm{pred}}^j}{p}\right]\frac{d_{\mathrm{pred}}^j}{a_i}=0.$

[0047] This represents the conditions under which the tunable coefficients yield a minimum model error irrespective of any source of variation. No matter what the parameter variations, this relation will be satisfied to guarantee minimum (least squares) error.

[0048] Note that $\frac{d_{\mathrm{pred}}^j}{a_i}$

[0049] is not a function of p since we assume that even though p shows variations, the sample slew and output load values at which d_actual is measured are not changed, i.e., $\frac{d_{\mathrm{pred}}^j}{a_i}$

[0050] depends on the values of the sample points which are independent of the sources of variation. Only the tunable coefficients of the model depend on the sources of variation. Also, since the delay model is a polynomial with coefficients a_i, i1,2, . . . , m, $\frac{d_{\mathrm{pred}}^j}{a_i}$

[0051] is not a function of the a_i coefficients either. This is because the tunable coefficients appear linearly in the polynomial model. Since d_pred is a function of a_i, which are in turn functions of p, we can write by simple chain-ruling: $\underset{j=1}{\overset{n}{}}2\left[\frac{d_{\mathrm{actual}}^j}{p}-\left\{\underset{k=1}{\overset{m}{}}\frac{d_{\mathrm{pred}}^j}{a_k}\frac{a_k}{p}\right\}\right]\frac{d_{\mathrm{pred}}^j}{a_i}=0.$

[0052] Note that the left hand side consists of a summation of terms in which we differentiate at each of the n pairs of sample value points for s^j_in and C^j_out which are considered constant and independent of the parameter variations.

[0053] Rearranging terms and dividing by 2, we obtain $\underset{j=1}{\overset{n}{}}\left[\frac{d_{\mathrm{pred}}^j}{a_i}\left\{\underset{k=1}{\overset{m}{}}\left(\frac{d_{\mathrm{pred}}^j}{a_k}\frac{a_k}{p}\right)\right\}\right]=\underset{j=1}{\overset{n}{}}\left[\frac{d_{\mathrm{actual}}^j}{p}\frac{d_{\mathrm{pred}}^j}{a_i}\right],i=1,2,,m.$

[0054] Since we can perform the above differentiation and manipulation with respect to each of m coefficients, we obtain m equations like the one above. The left-hand side of these equations is simply a linear combination of the $\frac{a_k}{p}$

[0055] terms, which are our unknowns. The coefficients of these terms, $\frac{d_{\mathrm{pred}}^j}{a_k},$

[0056] the sensitivity of the model (or predicted behavior) with respect to each tunable coefficient, a_k, are easily evaluated as constants at each of the n sample value points.

[0057] Assume that the circuit simulator is able to provide the sensitivity $\frac{d_{\mathrm{actual}}^j}{p}$

[0058] of each measurement (actual behavior) with respect to p, the parameter which shows variability. Then the right-hand sides of the equations are known. Thus we end up with m linear equations in m unknowns that can be solved to obtain the sensitivity of the tunable coefficients of the delay model with respect to one variable parameter of interest, $\frac{a_i}{p}.$

[0059] The physical significance of the above equation is as follows: if a parameter, p, has a small variation, p, the delay measurements will change. If the fitting procedure were to be repeated with the changed measurements, it would lead to changed tunable coefficients, a_i. The relationship between the change in the coefficient values and the variation of the parameter is given by the above equation, without having to repeat either the measurements or the fitting procedure. This is in contrast with prior-art perturbation methods.

[0060] Now we consider the situation with multiple variable parameters. The above analysis can easily be repeated for a second parameter. Interestingly, the left-hand-side coefficients are unchanged, and only the right-hand-side vector is changed since the sensitivities of d_actual to the second parameter will be different. As before, even in the case of multiple parameter variations, the sensitivities of the delay coefficients to the multiple parameter variations are governed by the equation above without having to repeat either the measurements or the fitting procedure. Thus if the linear equations for the first parameter were solved by the well-known LU (Lower-triangular/Upper-triangular) factorization method, the LU factors obtained can be reused for the second and every subsequent source of variability. Only forward and backward substitution steps are required to obtain the sensitivities of the coefficients of the delay equations to the second and every subsequent parameter It is well known that the factorization step is much more computationally expensive than forward and backward substitution. Hence, the computational effort required for each subsequent parameter is a small additional overhead after the initial LU factorization. The inventive method thus allows a significant reduction in computational effort over any prior-art method for the case of multiple parameter variations. LU factorization is a well-known procedure, described for example in Chapter 3 of L. T. Pillage, R. A. Rohrer and C. Visweswariah, Electronic Circuit and System Simulation Methods, McGraw Hill, 1995.

[0061] Finally, we note that the preferred method of computing sensitivities of the measurements with respect to the multiple sources of variation in the circuit simulator is by the adjoint method (see S. W. Director and R. A. Rohrer, The generalized adjoint network and network sensitivities, IEEE Transactions on Circuit Theory, vol. CT-16, number 3, pages 318-323, August, 1969), in which the sensitivity of a measurement with respect to multiple parameters can be computed by means of a single adjoint analysis.

[0062] A further preferred method of efficiently determining the measurements (actual timing behavior) and the sensitivities of the measurements to the sources of variation is described below.

[0063] Consider a library that has a number of gate types (e.g., inverter, NAND, NOR), each of which is available with various combinations of transistor sizes (i.e., a number of power levels, ratios and taper ratios). We will call each unique gate in the library (representing a unique combination of circuit topology and transistor sizes) a cell. The characterization process requires us to model the timing of each edge of the timing graph contributed by each cell. For each such edge, we must simulate the cell at a multiplicity of preselected sample value points of input slew and output capacitance, and at each sample value measure the delay, output slew, and sensitivities of delay and output slew to all sources of variation. One way to do this is to invoke a separate circuit simulation for each sample value, carrying out an adjoint analysis in addition to each nominal circuit simulation. Another way to do this is to concatenate in time the simulation at multiple sample point values to form groups of simulations. For example, the four single-input transitions of a 2-input NAND gate can be measured in a single simulation by cycling the inputs through patterns such as A11101, B10111. Further, the simulations for multiple input slews can then be further concatenated in time. Concatenation in time saves the computer time required by the circuit simulator to set up the circuit equations and perform a DC analysis from scratch for each new simulation. Thus the setup and DC analysis computer time is leveraged over all the measurements in a group.

[0064] Concatenation in time also vastly improves the efficiency of computing the sensitivities of the measurements to sources of variation. The right-hand-side of the final set of equations above can be thought of as a linear combination of the sensitivities of individual measurements to sources of variation, or in fact the sensitivity of a linear combination of individual measurements. The coefficients of the linear combination are known once the nominal simulations have been performed. Even better, this linear combination is the same for all sources of variation. We will call this linear combination of individual measurements a scalar function. What we desire is the sensitivity not of individual measurements but of the scalar function, with respect to each source of variation.

[0065] The sensitivity of any scalar and differentiable function of multiple measurements with respect to any number of sources of variation can be computed in a single adjoint analysis (see A. R. Conn, R. A. Haring and C. Visweswariah, Method of efficient gradient computation, U.S. Pat. No. 5,886,908, issued March 1999). Thus grouping all the measurements to be simulated and time-concatenating within each group saves computer time both in recording measurements and in determining the sensitivities of the measurements.

[0066] Of course, there are some practical limitations to this type of grouping. Cells of different gate types cannot be grouped. On occasion, the overhead of an extra pattern or vector of simulation may be required to bring the inputs and internal nodes of the gate to the logical state required to make the next required transition. Measurements at multiple sample values of input slew can easily be concatenated in time. Depending on the limitations of the circuit simulator employed, there may be other constraints in deciding which measurements to concatenate.

[0067] Within these limitations, the maximal grouping will buy the most gains in computer efficiency during the library characterization process.

[0068] The steps involved in computing sensitivities of delay model coefficients of a single pin-to-pin transition of a single library gate with respect to variation parameters are shown in the flow diagram 400 of FIG. 4. In step 410, the circuit simulator is employed to take measurements at numerous sample value points of s_in and C_out (the model arguments), and to compute the sensitivities of each of these measurements of actual timing behavior to the sources of parameter variation, preferably by the adjoint method and preferably using grouping and time-concatenation. In box 420, curve fitting is performed by well-known prior art methods to minimize the least squares error between the predicted and measured delays. In box 430, a set of linear equations is formed by the inventive method with sensitivities of the tunable delay coefficients to the first variation parameter as the unknowns. These equations are then solved by the well-known LU factorization method. In box 440, the LU factors are repeatedly used with different right-hand side values to compute sensitivities with respect to subsequent variable parameters efficiently.

[0069] The steps of the flow chart of FIG. 4 are repeated for every pin-to-pin transition of every member of the circuit library to obtain the sensitivities of all the delay coefficients to each of the sources of parameter variation. Depending on the details of the delay model, the above procedure must be repeated for each aspect of the timing model such as output delay and output slew. Referring back to FIG. 3, all of the procedure of FIG. 4 applied to every pin-to-pin transition of every member of the circuit library enables us to determine the gate delay models along with relevant sensitivities as indicated in Box 330 of FIG. 3. In other words, repeatedly applying the steps of FIG. 4 results in accomplishment of Box 330 of FIG. 3.

[0070] Now we describe an alternative embodiment of this invention, which is a method of statistical static timing analysis, shown as Box 340 of FIG. 3. In this part of the embodiment, we assume that timing models along with relevant sensitivity information are available for each edge of the timing graph. These models are preferably but not necessarily computed by the inventive method described above in FIG. 4. However, any prior-art timing model could be used including timing models expressed in the DCL language (see IEEE standard for integrated circuit (IC) delay and power calculation system, IEEE Standard 1481-1999, 1999, pages 1-390, available at http://fb.ieeexplore.ieee.org/ie15/6837/18380/00846710.pdfisNumber18380&prodstandards& arnumber00846710).

[0071] In the ensuing detailed description, we will assume that the circuit is a combinational circuit whose timing can be represented by a directed acyclic graph (DAG) and that we are interested in late-mode delays only. The nodes of the graph are the rising and falling signals of the circuit, and all pin-to-pin input-to-output transition delays of all gates and wires are the edges of the graph. It is to be understood that while these simplifications are made to clarify the exposition, the extension to early-mode analysis, to circuits with different kinds of latches, to dynamic circuits, to circuits with loops of transparent latches, and to situations dealing with slack instead of delay can be accomplished by one of ordinary skill in the art.

[0072] FIG. 5 shows a flow chart (500) of the steps involved in the inventive method of a statistical timing analysis process. In other words, the details of Box 340 of FIG. 3 can be found in flow chart 500 of FIG. 5. The individual steps are explained in detail in the subsequent paragraphs. Note that any circuit delay model that provides sensitivities of the tunable coefficients with a respect to the sources of variation (330) can be used with process 500. However, the process 330 described above is preferred.

[0073] The essence of the procedure (500) is to recursively consider two paths through the circuit at a time and determine the probability density function (PDF) of the maximum delay of the two paths, taking into account the covariance of the delay of the two paths. To generally address the effects of any given source of variation, covariance must be accounted for because two or more of the paths are likely to be commonly influenced by the source of variation. There are two main reasons for this correlation. The first is that the two paths may share one or more wires or gates, and so the path delays show correlation. The second is that all paths are exposed to the same or similar manufacturing and environmental variations. For example, all paths will be exposed to the same oven temperature during manufacturing. Thus the paths would be expected to exhibit behavior with a common influence (correlation) as a result. As this procedure continues, an edge criticality vector Q of the probability of each edge being critical is maintained.

[0074] The process begins with a nominal static timing analysis by prior art methods (510). These static timing analysis methods provide the list of nominally most-critical paths in order of criticality and the nominal slack of each of these paths. The next step (520) is to determine the matrix (a square matrix of size number of edges by number of edges) which represents the variance and covariance of all the edges of the timing graph which was constructed during the nominal static timing analysis. Let us denote by d_pred the change in predicted delay as a consequence of the statistical manufacturing and environmental parameter variations (collectively the sources of variation). Step 520 obtains statistical parameters like the variance of d_pred of a single edge in the delay graph. For instance, given the covariance matrix of all parameter variations, denoted by V, the variance of d_pred, denoted by ²(d_pred), can be obtained by applying the gradient vectors appropriately as shown below: ${}^2\left(d_{\mathrm{pred}}\right)={\left[\frac{d_{\mathrm{pred}}}{P}\right]}^{T}V\left[\frac{d_{\mathrm{pred}}}{P}\right].$

[0075] Assuming two different edges (associated with the same or different gates), say edge i and edge k, and assuming that the covariance matrix V covers all parameter variations relevant to both edges, then the covariance of d_predⁱ and d_pred^k denoted by the statistical expectation E $E\left[d_{\mathrm{pred}}^i{d}_{\mathrm{pred}}^k\right]$

[0076] is readily determined by $E\left[d_{\mathrm{pred}}^i{d}_{\mathrm{pred}}^k\right]={\left[\frac{d_{\mathrm{pred}}^i}{P}\right]}^{T}V\left[\frac{d_{\mathrm{pred}}^k}{P}\right].$

[0077] Thus the variances and correlations that make up the matrix are determined in box 520.

[0078] The matrix C is built (step 530) that is a sparse path-edge matrix that contains only 1s and 0s. Each path is represented by a row that contains a 1 in columns corresponding to edges in the path and a zero in columns corresponding to all other edges. The paths are sorted by criticality (530). As a practical matter, such a large matrix is never built or computed; rather, it is used here for notational convenience. The mean delays of paths are simply

[0079] where d_path is a vector of mean path delays and d_edge the vector of mean edge delays.

[0080] We take the nominally two most critical paths of the circuit, since they have the highest probability of contributing to the probability density function of the most critical path However, the method will work (perhaps less efficiently) irrespective of the path order chosen. We call the rows of C corresponding to these paths Q₁ and Q₂, respectively (step 530). Using these two rows, we can compute the 22 submatrix of CC^T (step 540). The edge delays are assumed to be distributed as multivariate normal distributions. Thus the means, variances and the covariance of the two paths being considered are now known. It is to be noted that because of the sparsity and special structure of C, the product CC^T can be computed efficiently.

[0081] The next step (550) is to compute the probability density function of the worst delay of the two paths being considered. The probability that the first path is binding (i.e., dominates the timing of the other path) and has a delay of is $\frac{1}{\sqrt{2}}\frac{1}{{}_1}{}^{-\frac{{\left(-{}_1\right)}^2}{{}_1^2}}\mathrm{erf}\left(\frac{\frac{-{}_2}{{}_2}-\frac{-{}_1}{{}_1}}{\sqrt{1-{}^2}}\right),$

[0082] where $\mathrm{COV}\left(\mathrm{path1},\mathrm{path2}\right)=\left[\begin{array}{cc}{}_1^2& {}_1{}_2\\ {}_1{}_2& {}_2^2\end{array}\right]$

[0083] and ₁ and ₂ are the mean delays of the two paths and the covariance matrix determined in step 540 allows us to determine ₁, ₂ and . Likewise, we can write the probability that path 2 dominates and has a particular delay. The summation of these two probability distributions over all values of gives us the required probability distribution function. The resulting probability distribution function is normalized by computing its mean and variance, and then approximated by a normal distribution with the same mean and variance values (step 550). Further, the integral of the probability that path 1 dominates path 2 over all values of gives us the total binding probability that path 1 dominates (step 550). The physical significance of the computed PDF is that if the circuit were to consist of only the paths considered thus far, the circuit would have a critical path delay as predicted by the probability density function. Further, the binding probabilities tell us the probabilities that the critical path delay is limited by the first path or the second, and, of course, the two binding probabilities add up to 1. The above probability calculations are well-known, see for example P. Z. Peebles, Jr., Probability, random variables, and random signal principles, McGraw Hill, Second Edition, 1987. A check is made (step 560). If the paths of the circuit are exhausted or the change in the probability density function is negligibly small compared to the desired accuracy, the process 500 terminates (step 570). If not, the process continues to step 580.

[0084] The next step (580) is to update our vector of edge criticality probabilities by setting

[0085] to represent the cumulative edge criticality probabilities thus far. There is no approximation involved in the above formula since the probability that an edge is critical is the sum of the probabilities that each path that includes this edge is critical. Q₁ represents the vector of edge criticality probabilities considering only the set of paths processed so far, and b₁ is the probability that one of the paths from this set is critical (or binding). Similarly, Q₂ represents the vector of edge criticality probabilities of the next path being considered (a vector of 1s and 0s), and b₂ is the probability that that path is critical.

[0086] We then repeat the above recursive procedure by computing the required means, variances and covariances (step 590) and then returning to step 550. At some point, we will either run out of paths or the change in the probability distribution function will be negligibly small (step 560). At this point, the procedure terminates (step 570), and thus the statistical timing analysis (box 340 of FIG. 3) will have been accomplished.

[0087] The above method can advantageously be combined with uncertainty-aware tuning to render it more efficient. If the circuit undergoes uncertainty-aware tuning before the statistical static timing analysis, then there will be relatively few critical paths that must be considered since the wall of equally critical paths is avoided. Thus the statistical timing analysis procedure can be terminated earlier, and the inaccuracy due to the approximation of making the PDF Gaussian at each path-merging recursion will be reduced. Uncertainty-aware tuning is described in X. Bai, C. Visweswariah, P. N. Strenski and D. J. Hathaway, Uncertainty-aware circuit optimization, Proceedings of the Design Automation Conference, June 2002, New Orleans, La. and also in C. Visweswariah, X. Bai, D. J. Hathaway and P. N. Strenski, Parameter variation tolerant method for circuit design optimization, U.S. patent application filed May 2002, Docket FIS9-2002-0034-US1.

[0088] Once we have the probability distribution function of the most critical delay of the circuit, the information can be used in several interesting ways. We can predict the yield of the circuit at a given performance requirement. Since we have the full PDF of the critical path delay, we simply accumulate the probability that the circuit has the required performance or better to obtain the yield. More specifically, for instance, given a maximum delay d_max that can be tolerated, we compute the percentage of circuits exhibiting a delay not larger than d_max which would be the required yield. We can also calculate the performance requirement in order to guarantee a certain yield. To do this, we simply compute the cumulative probability distribution by integrating the PDF and then sample the cumulative probability distribution at the required yield to obtain the performance of the circuit that would guarantee the required yield.

[0089] We can also obtain the final probability of each edge being critical. To do this, we simply inspect the Q₁ vector of binding probabilities at the end of the statistical timing. We can also compute the probability that a particular path dominates the timing of the circuit. To do so, we simply accumulate an updated binding probability for each considered path as the statistical timing analysis proceeds. The accumulation implies the product of the binding probabilities that any group of paths that includes the path of interest is binding. All of the above probabilities are extremely useful in guiding circuit optimization techniques to improve performance and/or yield and/or profitability.

[0090] It is to be understood that the same techniques can be applied with slight modification to early mode timing in which we seek to find the earliest time at which a signal changes from the stable logical value of the previous cycle of operation. In early mode, the fastest path dominates the timing of the circuit. The previously described procedure can be repeated with two minor changes. Given two paths, path 1 and path 2, the probability that the first path is binding (i.e., dominates the timing of the other path) and has a delay of is $\begin{array}{cccc}& & & \frac{1}{\sqrt{2}}\frac{1}{{}_1}{}^{-\frac{{\left(-{}_1\right)}^2}{{}_1^2}}\left[1-\mathrm{erf}\left(\frac{\frac{-{}_2}{{}_2}-\frac{-{}_1}{{}_1}}{\sqrt{1-{}^2}}\right)\right],\\ \mathrm{where}& & & \\ & & & \mathrm{COV}\left(\mathrm{path1},\mathrm{path2}\right)=\left[\begin{array}{cc}{}_1^2& {}_1{}_2\\ {}_1{}_2& {}_2^2\end{array}\right]\end{array}$

[0091] and ₁ and ₂ are the mean early-mode delays of the two paths and the covariance matrix determined in step 540 allows us to determine ₁, ₂ and . Likewise, we can write the probability that path 2 dominates and has a particular delay. The summation of these two probability distributions over all values of gives us the required probability distribution function. The binding probabilities are then computed in a similar fashion to the previous procedure, except that the binding probability of a path now represents the probability that the path being considered has the lowest delay among all the paths considered thus far. So, other than computing the PDF and binding probabilities in a slightly different manner, the procedure for statistical static timing in early mode is similar to that for late mode.