Statistical timing analysis methods for circuits are described which compensate for circuit elements having correlated timing delays with a high degree of computational efficiency. An quadratic timing model is used to represent each delay element along a circuit path, wherein each element's delay has a first-order relationship to local variations and a second-order relationship to global variations. Propagation of the modeled delays through the circuit is efficiently done via straightforward ADD operations where an input propagates through another element in a circuit path, and via a MAX operation (or an approximation thereof) where two or more inputs merge at an intersection. The inputs to the MAX operator can be tested for gaussianity, and can be processed by the MAX operation (or its approximation) if they are substantially gaussian. Otherwise, they may be stored in a tuple for processing at later points along the circuit path.