Monte Carlo and Worst-Case Circuit Analysis using LTSpice

SPICE is a handy tool for evaluating circuits without having to first breadboard them, and through its "directives," it provides a powerful method for analyzing how a circuit might perform with components exhibiting real-world tolerances. 

One such method of "real-world" analysis is Monte Carlo analysis, which, with each new analysis run, randomly varies parameters (within their user-defined limits) to give the user a useful picture of actual circuit performance.

However, as a circuit designer, I'm most often interested in worst-case performance.  That is, I want to know how a circuit performs at the extremes of component values, to ensure that  I've met whatever design specification I'm designing to.  And although Monte Carlo analysis can tell me what the performance is at these limits, if the circuit contains many components, it can take quite a lot of runs before its random selection of parameter values happens to simultaneously correspond to the worst-case limits of all of the components (and it's quite possible that I'll never see the true worst-case limits -- after all, it's a matter of chance).

To truly evaluate performance at a circuit's worst-case limits, we can perform a "worst-case" analysis in lieu of a Monte Carlo analysis.  This analysis has an added advantage, too, in that not as many runs are required to ensure that we've truly evaluated all of the components over all possible variations in their tolerances -- after all, we don't need to measure between the tolerances limits (which is what Monte Carlo analysis does).  For example, if one of the components is a 10K ohm resistor with a 5% tolerance, worst-case analysis will only use resistance values of 10.5K and 9.5K for this component  and not any other value between these limits (whereas Monte Carlo analysis would randomly use any value between, and including, these limits).

Linear Technology's LTSpice can handle both Monte Carlo analysis as well as worst-case analysis.  Unfortunately, it's not obvious how to do this from their "Help Topics."  For example, "Monte Carlo", when entered into LTSpice's search field, returns no results.  Not much use, that!

Nevertheless, LTSpice does indeed have a pre-defined Monte Carlo function.  This is the "mc" function, and a search of the Help Topics for "mc" will point to the .PARAM topic, and under this heading we find the function mc (x,y), which, when invoked, returns a "random number between x*(1+y) and x*(1-y) with uniform distribution."

To use this function, rather than define a resistor's value as, say, 10K, we define it as "{mc(10K,0.05)}", where 10K is its nominal value, and 0.05 is its tolerance (5%).  (An example will follow below).

OK, so the predefined mc function handles Monte Carlo analysis, but there is no pre-defined "worst case" function.   We need to create this ourselves, but it's not too difficult.  This can be done using Spice's ".function" directive.  Here's an example of a function for worst-case analysis (we'll use this later, too):

     .function wc(nom,tola) if (run == 1, nom, if(flat(1)>0,nom*(1+tola),nom*(1-tola)))


In this definition:

• .function is the Spice directive for defining a function.
• wc is the name I've given this instance of this worst-case function.
• nom is the "nominal" value (e.g. 10K for a 10K resistor).
• tola is the tolerance (defined elsewhere with a .param directive (e.g. 0.1 for a 10% tolerance)).
• run is a variable (defined elsewhere in the .step directive) identifying the current run count).
• flat(1) is a Spice function that returns a random number between -1 and 1. 

(Note that the definition of .function and flat can both be found in LTSpice Help).

So how do we use this function we've just defined?

Take as an example a 10K ohm resistor.  Normally, we'd just label its value in the LTSpice schematic as "10K".  However, to vary its value between its worst-case values, we instead use a more complex label for its value.  Rather than entering "10K", in this case we'll enter "{wc_a(10K,tola)}" into the component's value field.  (Note the use of the curlicue brackets).

So what happens when we run the analysis?

If this is the first analysis run (run is defined elsewhere in the .step directive and simply is the current run being performed), then, because run = 1, the function returns the value assigned to nom, in this case, 10K.

But for each new run after the first run, and for each component defined with a wc_a function, the function flat(1) is re-evaluated for that component.  If the random result returned from flat(1) is greater than 0, then the tolerance is added to the component's nominal value (e.g. for a resistor whose nominal (nom) value is 10K, if tola is 0.1 (10% tolerance), the resistor's value is set to 11K).  Otherwise, the tolerance is subtracted from the component's nominal value (e.g. the 10K resistor is set to 9K).

Here's a demonstration of both Monte Carlo and Worst Case analysis.  Consider this basic circuit: 

(Click on image to enlarge)

Given these component values, a frequency sweep of the input from 1 Hz to 1KHz shows that the circuit has the following gain and phase transfer function when measured at its Vout node:

(Click on image to enlarge)

But what happens when the component values vary over their tolerance range?  Let's suppose that the resistors have 10% tolerance and the capacitors have 20% tolerance.  Let's perform a Monte Carlo analysis on this circuit, given these tolerance values. The same circuit, but now set up with its Monte-Carlo functions and .param directives, looks like this:

(Click on image to enlarge)

(Note that within the "mc" function, I'm not setting the tolerance field to an actual number (although I could have done it this way, too).  Instead, I'm using a separate directive (.param) to define the tolerances (in this case, 10% and 20%) globally.)

Running the Monte Carlo analysis 1000 times (via the directive ".step param run 1 1000 1") gives us the following spread of gain and phase plots:

(Click on image to enlarge)

But we can't be sure that we've truly evaluated worst-case performance.  So let's instead use our new "worst-case" function for a worst-case evaluation.

The schematic, with its new Spice directives and functions, now looks like this:

(Click on image to enlarge)

And the analysis output, after 40 runs, looks like this:

(Click on image to enlarge)

Note the discrete intervals between plots.  This is because the worst-case analysis is only using component values that are at the +/- tolerance limits for each component, and not any intermediary values (except for the first of the 40 plots, which uses the nominal component values for its analysis).

(Note:  the above was edited on 29 May 2020 to replace two "worst-case" functions (wc_a and wc_b) with a single worst-case function, wc.)


Optimized Worst-case Analysis
(This section added 29 May 2020)

The worst-case analysis component values selected on a run-by-run basis for the simulations, above, are a function of random numbers.  Therefore, to ensure that all combinations of "worst-case" component values have been evaluated, many runs need to be evaluated (above and beyond the "optimal" number of runs equal to 2^N, where N is the number of components being varied).

In the "Comments" section of this blog post, Harry Dymond (Electrical Energy Management Group, University of Bristol, UK) on November 9, 2012, posted an excellent technique for optimizing the number of runs necessary to perform a complete worst-case analysis.

If we consider, for worst case analysis, that each component has two possible values (nominal value plus tolerance, and nominal value minus tolerance), then these two value "states" can be represented by a single binary bit.  In Harry's technique, "1" represents "nominal value plus tolerance", and "0" represents "nominal value minus tolerance.

Therefore, if we have N components in our circuit that we would like to vary the tolerances of, we can represent these components with N bits.

If we then stepped through all possible combinations of these N bits (from 0 to (2^N)-1), and at each step ran a simulation using the appropriate value of each component as defined by the state of its bit for that step, we would step through all possible combinations of the worst-case values without repeating or skipping a combination of values.

In other words, we would have optimized the number of runs to be the minimum set required for a complete worst-case analysis of our circuit

The table below demonstrates the component "bit" values for the 16 runs required for a complete worst-case analysis of 4 components:

Each component (that will be varied) is assigned a unique "Component Index" number.  This index starts at 0 and increments by 1 for each component.

For example, if we are going to vary the tolerance of four components in a circuit, these four components are assigned index values starting with 0 for the first component, 1 for the next, 2 for the third, until we reach 3 for the fourth (and last) component.

You can think of these index values as each being a "bit position" in the N-bit word (where N, in this particular example, would be 4).  For example, index 0 refers to the 2^0 bit location, index 1 refers to the 2^1 bit location, etc.  You can see this in the table, above.

This table determines how the tolerances are set for each run.  For Run = 0, all entries in the column for that run equal 0.  Therefore, all four component values will have their tolerances subtracted from their nominal values.

For the next run, Run = 1, the component whose index is 0 will have its tolerance added to its nominal value.  All other components will have their tolerances subtracted from their nominal values.

For the next run, Run = 2, now the component whose index is 1 will have its tolerance added to its nominal value.  All other components will have their tolerances subtracted from their nominal values.

And for the next run, Run = 3, now the two components whose indexes are 0 and 1 will have their tolerances added to their nominal values.  All other components will have their tolerances subtracted from their nominal values.

And so it goes, in a binary-counting fashion, until we reach the last run count of 15.


How do we do this using LTSpice?  Here's the new schematic, using the same four components that we analyzed in our earlier, non-optimized worst-case analysis, above.

The differences between this new LTSpice simulation and the earlier version are:

1.  The "wc()" function now has a third input variable: "index" (i.e. the "component index", and the function is now:

.function wc(nom,tola,index) if (run == -1, nom, if(binary_digit(run,index),nom*(1+tola),nom*(1-tola)))

2.  The wc() function for each component (in each component's "value" field) is assigned a unique index.  The index for the first component is 0, and indexes increment sequentially by 1.

3.  The "binary_digit" function is a new function.  It returns either a 1 or a 0, depending upon run-number and component index:

.function binary_digit(run,index) floor(run/(2**index))-2*floor(run/(2**(index+1)))

4.  The "run" parameter now starts at -1 (as mentioned by "anonymous" on December 13, 2012, in the comments section, below).  When "run" equals -1, the circuit simulation is run with nominal component values.

There are 4 indexes (and thus 4 bits) for this circuit.  Thus, a complete worst-case simulation will require 16  runs to test all combinations of worst-case tolerances, plus there is one additional run with components set to their nominal values.  So there are 17 runs, total.


OK, let's run the example!

First, here's the LTSpice schematic, again:

And here are the plots of the voltage at the "vout" node:

Please note that this technique (as described by Harry Dymond in his November 9, 2012 comment in my comments section below) is also described in the following article written by Linear Technology Corporation, in 2017:

https://www.embedded-computing.com/articles/getting-the-worst-case-circuit-analysis-with-a-minimal-number-of-ltspice-simulation-runs

(Note, too, that in my example, above, I've replaced Harry's use of "powerOfTwo" with the shorter word "index", used by the LTC authors.)


Other Notes and Caveats:

1.  In the "comments" section, "Thomas", on December 2, 2016, points out that worst-case analysis could miss LC circuit resonances.  If your circuit contains inductors and capacitors, Monte Carlo analysis (in lieu of Worst-case analysis) would probably be more appropriate.

2.  Further information on Worst-Case and Temperature Analysis can be found here:

http://k6jca.blogspot.com/2012/07/worst-case-and-temperature-analysis.html

(And check out the other comments in the comments section, below, too).


Resources:

LTSpice is available for free from Linear Technology.  You can find it here