Monte Carlo Tolerance Analysis

Main idea

When the shortcuts break down, because the formula is curved, the parts aren't normal or there are too many to track, let the computer build thousands of virtual assemblies from random parts and count how many work. That fraction is itself an estimate, with standard error $\sqrt{p(1 - p)/N}$ , so four times the runs buys half the error.

The estimate and its standard error

{\htmlData{sym=0}{\hat p}} = \frac{\#\{\text{runs in spec}\}}{\htmlData{sym=1}{{N}}}, \qquad {\htmlData{sym=2}{\operatorname{SE}(\hat p)}} = \sqrt{\frac{p(1 - p)}{\htmlData{sym=1}{{N}}}}

What the symbols mean

$\hat p$
Estimated yield
The fraction of simulated assemblies that meet the spec.
$N$
Number of runs
How many virtual assemblies were built.
$\operatorname{SE}(\hat p)$
Standard error
How much the estimate changes from one simulation to the next; in practice, plug in p̂ for p.

Planning the number of runs

N = \frac{z^2\,p(1 - p)}{{\htmlData{sym=0}{E}}^2}, \qquad \frac{\operatorname{SE}(\hat {\htmlData{sym=1}{q}})}{\htmlData{sym=1}{{q}}} = \sqrt{\frac{1 - {\htmlData{sym=1}{q}}}{{\htmlData{sym=2}{Nq}}}} \approx \frac{1}{\sqrt{{\htmlData{sym=2}{Nq}}}}

What the symbols mean

$E$
Margin of error
The ± half-width you want, at the confidence level set by z.
$q$
Failure fraction
1 − p: the rare outcome you're trying to measure.
$Nq$
Expected failures
How many failing runs you expect to see; about 100 gives ±10%.

Making random parts

X = \mu + \sigma {\htmlData{sym=0}{Z}}, \qquad X = a + (b - a)\,{\htmlData{sym=1}{U}}, \qquad X = {\htmlData{sym=2}{F^{-1}}}({\htmlData{sym=1}{U}})

What the symbols mean

$Z$
Standard normal draw
A random number from N(0, 1).
$U$
Uniform draw
A random number between 0 and 1, what computers produce directly.
$F^{-1}$
Inverse CDF
Turns a uniform draw into a draw from any distribution with CDF F.

Explanation

Monte Carlo simulation is simply building lots of pretend products on a computer, each from randomly chosen parts, and counting how many come out good. The name comes from the casino in Monaco, a nod to the dice-rolling at its heart. The formulas are short: $\hat p$ (good runs divided by total runs), its standard error $\sqrt{p(1 - p)/N}$ (the same as for any poll), the run-count formula (that standard error solved for $N$ ) and three recipes for turning random numbers into random parts. A computer does the heavy lifting; you just have to read the answer honestly.

When the shortcuts run out. An RC low-pass filter should cut off at $f_c = \dfrac{1}{2\pi RC} = 1.00$ kHz, within ±8%. Its resistor is ±5% and its capacitor ±10%, and parts in both bins are spread roughly uniformly across the band: makers sort the best parts into tighter bins, so the wide bins come out flat. The formula is a quotient, the inputs aren't normal and the tolerances are big. The linear-and-normal shortcuts of the last two lessons get shaky. So let the computer build the filters.

The recipe.

Model each input from data: its distribution, mean and SD. Probability plots tell you whether "normal" is fair.
Draw $N$ random values of every input, independently (or with the correlation you measured).
Compute the output for each of the $N$ virtual assemblies with the exact formula, no linearizing.
Count the fraction in spec: $\hat p$ .
Report $\hat p \pm z\cdot\operatorname{SE}$ , because the answer is itself random.

Why the standard error is $\sqrt{p(1 - p)/N}$ . Each run lands in spec with probability $p$ , independently of the others, so the count $X$ of good runs is binomial with mean $Np$ and variance $Np(1 - p)$ . Then $\hat p = X/N$ has

E(\hat p) = p, \qquad \operatorname{Var}(\hat p) = \frac{Np(1 - p)}{N^2} = \frac{p(1 - p)}{N}.

That's the standard error of a confidence interval for a proportion: a simulation is a survey of virtual parts. Setting $z\cdot\operatorname{SE} = E$ and solving gives the number of runs, $N = z^2p(1 - p)/E^2$ .

Rare failures need many runs. For a defect rate $q$ measured in ppm, the error that matters is relative: $\operatorname{SE}(\hat q)/q = \sqrt{(1 - q)/(Nq)} \approx 1/\sqrt{Nq}$ , one over the square root of the number of failures you expect to see. Seeing 100 failures gives about ±10%; at 50 ppm that takes 2 million runs. For very rare events engineers lean on the normal-tail formulas or smarter sampling schemes.

Where random parts come from. Computers produce uniform numbers $U$ on $(0, 1)$ . Then:

uniform on $[a, b]$ : $X = a + (b - a)U$ , with SD $(b - a)/\sqrt{12}$ , the resolution formula from Measurement Error;
normal: $X = \mu + \sigma Z$ , where $Z$ comes from a normal generator (one classic recipe is $Z = \sqrt{-2\ln U_1}\cos(2\pi U_2)$ );
anything whose CDF you can invert: $X = F^{-1}(U)$ , because $P(F^{-1}(U) \le x) = P(U \le F(x)) = F(x)$ . For an exponential lifetime that's $X = -\ln(1 - U)/\lambda$ ; for a Weibull it's $X = \eta\,(-\ln(1 - U))^{1/\beta}$ .

Simulation vs the formulas. Monte Carlo handles anything: curved formulas, flat or skewed parts, correlation, even logic like "scrap any resistor outside its own band." The price is computer time and an answer with noise in it. The formulas are instant and exact under their assumptions. Good practice uses both: if the simulation and the RSS estimate disagree by more than a few standard errors, one of the assumptions is wrong, and you've learned something about the design.

A common mistake. It's tempting to report a simulated yield as if it were exact, because a computer produced it. It isn't: rerun with fresh random numbers and you'll get a slightly different $\hat p$ , so always give its standard error. Related slips: running 1,000 trials to estimate a 100 ppm defect rate (you expect to see 0.1 failures), feeding the same stream of random numbers to two inputs (which makes them perfectly correlated), and assuming normal parts when the bins are flat.

See it

N\ \text{(simulated filters)}

500

k\ \text{(new random parts)}

1.00

T_C\ \text{(capacitor tolerance in percent)}

10.0

Each dot is one simulated RC filter: its cutoff frequency in kHz across, with a random height so the dots don't pile up, colored by whether it lands inside the ±8% spec. Draw fresh random parts with k and p̂ jumps by about one standard error; raise N and the jumps shrink like 1/√N. With a ±10% capacitor the true yield is about 75%, while the normal shortcut claims 78%. The output isn't normal, and only the simulation notices.

How precise is a simulation? The band shows where 95% of estimates land, in percentage points from the true yield, against the number of runs N. It narrows like 1/√N: going from 125 to 500 runs halves it, and so does going from 500 to 2,000. At p = 0.5 and 125 runs it's about ±8.8 points. As the true yield rises toward 0.99 the band thins, so a high yield is cheap to confirm, while its tiny failure rate is expensive to pin down.

Worked examples

Example 1

A Monte Carlo study builds $2{,}000$ virtual RC filters and finds $1{,}862$ with their cutoff in spec. Estimate the yield and give a 95% interval.

1
$\hat p = \frac{1862}{2000} = 0.931$ .

Example 2

How many runs are needed to estimate a yield near $0.99$ to within $\pm0.002$ with 95% confidence?

1
$N = \frac{z^2p(1 - p)}{E^2} = \frac{1.96^2(0.99)(0.01)}{0.002^2}$ .

Practice

Type your answer and check it. Solve them all to complete the lesson.

Problem 1+5 XP

A Monte Carlo study of $5{,}000$ simulated gearboxes finds $4{,}870$ with acceptable end play. Find the standard error of the estimated yield, in percentage points, to 3 decimal places.

Problem 2+5 XP

You want to estimate an assembly yield that's probably near $95\%$ to within $\pm0.5$ percentage points with 95% confidence ( $z = 1.96$ ). How many simulation runs do you need? (Round up.)

Problem 3+5 XP

A simulation of $10{,}000$ circuits estimates a failure rate that's really $0.1\%$ (1,000 ppm). What is the relative standard error of the estimate, $\operatorname{SE}(\hat q)/q$ , in percent? Give 3 significant figures.

Problem 4+5 XP

To simulate bearing lives that follow a Weibull distribution with $\beta = 2$ and $\eta = 1000$ h, a program solves $F(T) = U$ for each uniform random number $U$ , where $F(t) = 1 - e^{-(t/\eta)^\beta}$ . What life does $U = 0.5$ produce, in hours? Give 3 significant figures.

Problem 5+5 XP

Three spacers have lengths uniform on $[9.95, 10.05]$ mm and stack into a slot of $30.12$ mm. One Monte Carlo run draws $U = 0.82$ , $0.15$ and $0.47$ for the three spacers. What gap does this virtual assembly have, in micrometers?

Problem 6+5 XP

In which situation is a Monte Carlo analysis most clearly worth it, compared with the root-sum-square formula?

Problem 7+5 XP

A four-part stack has spec limits at $\pm2\sigma$ of the assembly (σ from RSS), so the normal model predicts a yield of $0.9545$ . A Monte Carlo run with the parts' real, flat distributions finds $38{,}292$ in spec out of $40{,}000$ . How many standard errors is the simulated yield from $0.9545$ ? (Use the standard error computed from $p = 0.9545$ .) Give 3 significant figures.

Challenge problem+25 XP

The hardest one here. Get it, and you've got this lesson.

A bracket's strength $S$ is normal with mean $50$ kN and SD $4$ kN; the load $L$ on it is normal with mean $35$ kN and SD $3$ kN, independent of $S$ . It fails when $L > S$ . How many Monte Carlo runs are needed to estimate the failure probability with a relative standard error of 10%? Give 3 significant figures.

Where it's used

Monte Carlo tolerance analysis is built into mechanical CAD packages and circuit simulators: designers run thousands of virtual builds of an engine, a lens stack or a power supply before any metal is cut. Chip designers simulate the random variation of every transistor to predict yield, and structural engineers estimate the failure probabilities of bridges and aircraft parts by sampling loads and strengths. The same logic drives risk analysis in finance and in project scheduling.

Estimating π with Random DartsStats · Computer Science How Many Resistors Meet Spec?Stats · Electrical Engineering

Propagation of Error Reliability and Hazard Functions

Monte Carlo Tolerance Analysis

Main idea

The estimate and its standard error

{\htmlData{sym=0}{\hat p}} = \frac{\#\{\text{runs in spec}\}}{\htmlData{sym=1}{{N}}}, \qquad {\htmlData{sym=2}{\operatorname{SE}(\hat p)}} = \sqrt{\frac{p(1 - p)}{\htmlData{sym=1}{{N}}}}

What the symbols mean

$\hat p$
Estimated yield
The fraction of simulated assemblies that meet the spec.
$N$
Number of runs
How many virtual assemblies were built.
$\operatorname{SE}(\hat p)$
Standard error
How much the estimate changes from one simulation to the next; in practice, plug in p̂ for p.

Planning the number of runs

N = \frac{z^2\,p(1 - p)}{{\htmlData{sym=0}{E}}^2}, \qquad \frac{\operatorname{SE}(\hat {\htmlData{sym=1}{q}})}{\htmlData{sym=1}{{q}}} = \sqrt{\frac{1 - {\htmlData{sym=1}{q}}}{{\htmlData{sym=2}{Nq}}}} \approx \frac{1}{\sqrt{{\htmlData{sym=2}{Nq}}}}

What the symbols mean

$E$
Margin of error
The ± half-width you want, at the confidence level set by z.
$q$
Failure fraction
1 − p: the rare outcome you're trying to measure.
$Nq$
Expected failures
How many failing runs you expect to see; about 100 gives ±10%.

Making random parts

X = \mu + \sigma {\htmlData{sym=0}{Z}}, \qquad X = a + (b - a)\,{\htmlData{sym=1}{U}}, \qquad X = {\htmlData{sym=2}{F^{-1}}}({\htmlData{sym=1}{U}})

What the symbols mean

$Z$
Standard normal draw
A random number from N(0, 1).
$U$
Uniform draw
A random number between 0 and 1, what computers produce directly.
$F^{-1}$
Inverse CDF
Turns a uniform draw into a draw from any distribution with CDF F.

Explanation

The recipe.

Model each input from data: its distribution, mean and SD. Probability plots tell you whether "normal" is fair.
Draw $N$ random values of every input, independently (or with the correlation you measured).
Compute the output for each of the $N$ virtual assemblies with the exact formula, no linearizing.
Count the fraction in spec: $\hat p$ .
Report $\hat p \pm z\cdot\operatorname{SE}$ , because the answer is itself random.

E(\hat p) = p, \qquad \operatorname{Var}(\hat p) = \frac{Np(1 - p)}{N^2} = \frac{p(1 - p)}{N}.

Where random parts come from. Computers produce uniform numbers $U$ on $(0, 1)$ . Then:

uniform on $[a, b]$ : $X = a + (b - a)U$ , with SD $(b - a)/\sqrt{12}$ , the resolution formula from Measurement Error;
normal: $X = \mu + \sigma Z$ , where $Z$ comes from a normal generator (one classic recipe is $Z = \sqrt{-2\ln U_1}\cos(2\pi U_2)$ );
anything whose CDF you can invert: $X = F^{-1}(U)$ , because $P(F^{-1}(U) \le x) = P(U \le F(x)) = F(x)$ . For an exponential lifetime that's $X = -\ln(1 - U)/\lambda$ ; for a Weibull it's $X = \eta\,(-\ln(1 - U))^{1/\beta}$ .

See it

N\ \text{(simulated filters)}

500

k\ \text{(new random parts)}

1.00

T_C\ \text{(capacitor tolerance in percent)}

10.0

Worked examples

Example 1

A Monte Carlo study builds $2{,}000$ virtual RC filters and finds $1{,}862$ with their cutoff in spec. Estimate the yield and give a 95% interval.

1
$\hat p = \frac{1862}{2000} = 0.931$ .

Example 2

How many runs are needed to estimate a yield near $0.99$ to within $\pm0.002$ with 95% confidence?

1
$N = \frac{z^2p(1 - p)}{E^2} = \frac{1.96^2(0.99)(0.01)}{0.002^2}$ .

Practice

Type your answer and check it. Solve them all to complete the lesson.

Problem 1+5 XP

A Monte Carlo study of $5{,}000$ simulated gearboxes finds $4{,}870$ with acceptable end play. Find the standard error of the estimated yield, in percentage points, to 3 decimal places.

Problem 2+5 XP

You want to estimate an assembly yield that's probably near $95\%$ to within $\pm0.5$ percentage points with 95% confidence ( $z = 1.96$ ). How many simulation runs do you need? (Round up.)

Problem 3+5 XP

Problem 4+5 XP

Problem 5+5 XP

Problem 6+5 XP

In which situation is a Monte Carlo analysis most clearly worth it, compared with the root-sum-square formula?

Problem 7+5 XP

Challenge problem+25 XP

The hardest one here. Get it, and you've got this lesson.

Where it's used

Estimating π with Random DartsStats · Computer Science How Many Resistors Meet Spec?Stats · Electrical Engineering