I have been slowly learning about sensitivity testing using the Neyer D-Optimal method. There is an excellant chance that I have misunderstood something. On the off chance that someone with a better understanding should happen by this page, I would appreciate any feedback that would help my understanding. The reason for conducting these R&D projects isn't to produce a NARAM report, (although that is a very nice side effect) it is to learn something new and perhaps even useful.
Cody Steele presented a B division Research and Development report at NARAM 46 that caught my attention. When I saw the title of "All-Fire/No-Fire Tests on Estes Igniters Using the Neyer D-Optimal Method" I knew that I had to attend the oral presentation. The reason is that I had previously run across Mr.Neyers work and in particular his paper on the D-Optimal method. There was also a better than average chance that I was the only one in the audience that had ever read this paper.
The oral presentation was a dissapointment. I know that it was a B division report and that expectations should be lowered but it was still not very impressive. The first problem was that Cody demonstrated a profound misunderstanding of how current limited power supplies work. He tested the igniters at 6V, 9V, and 12V and even went so far as to use different launch controllers for each voltage.
The fallacy here is that the voltage will only matter if the resistance in the circuit is large enough so that this voltage is reached on the power supply output before the current limit. Since the resistance of the igniter was about 1 ohm and the wiring added not much more, at the test currents of less than 2 amps the power supply always went into current limit before reaching even 6 volts output. So it was unsurprising that the results for each voltage were similar.
The next and biggest problem was that the test protocol stated was nowhere close to the protocol of the D-Optimal method. Cody's stated method was to apply a current, then increase the current and apply the higher current until the igniter ignited. Compare this to the flow chart in Figure 2. of Mr. Neyer's paper. (Note: It gets worse. There are at least two errors in the flow chart.)
Another problem was that the computation of the sample mean and variance is not a simple thing. At the time of the oral presentation it had been a while since I read the report so I had forgotten the details but I knew that it was not straightforward. So during the Q&A portion of the presentation I asked Cody if he purchased or had donated Mr. Neyer's software for conducting the test. The answer was "no". Later I asked Cody how he performed the calculation and he said that he used Microsoft Excel. I was unaware that Excel performed such an esoteric statistical calculation.
So the idea of repeating the tests to confirm the results has been bouncing around for a while. Now that NARTS has provided a CD with all of the NARAM 46 reports on it, the written report confirms the problems of the oral presentation.
Preparing to repeat these tests breaks down into two areas: hardware and software.
For the hardware I need a constant current source that I can adjust to the various stimulous levels required. In addition it would be very handy if it could provide a time limited pulse.
Not wanting to spend a lot of cash, I recalled a circuit for a constant current source I had seen on the Pass DIY web site. That circuit was for the Zen amplifier and regulated current on the "high" side. Since I preferred operating on the low side, a simple transformation was required. I happened to already have all of the parts required to build the current source. I also had a PIC programmer and a supply of 12C508 microcontrollers to make a pulser. Here is the resulting schematic. I plan on using a 12V battery as the power source so the maximum current will be limited because of the 1 ohm resistor. A practical upper limit for testing 1 ohm igniters will be 6 amps. If I need to go higher than that I am limited to 20 volts unless I add a 20V zener to protect the gate of the power FET.
One of the big decisions is just how long to make the current pulse. Cody used a stated time of one second and I suspect this was measured by hand. Not very repeatable. In addition, in most cases if you are launching a rocket with a single motor you aren't too worried about how long the igniter takes to function just so long as it does. A one second duration is fine for that. But if you want to launch a cluster then if the igniter functions after even half a second the rocket is likely to be gone. So a shorter time is required.
After banging around some simulations in Rocksim I settled on a pulse width of 100 ms. For the simulations I created a rocket with a mass exactly equal to one motor and then determined how long it took it to travel 6" when loaded with the same motor type.
I also dug out an old DATAQ DI-194 data acquistion board. This will only do 240 samples a second but it will still provide a decent record of the current applied to the igniter. I hacked up some sample code provided by DATAQ to capture data from the tests and save it to disk.
So I now have a current source but that was the easy part.
Figuring out how the D-Optimal test works and how to do it myself has been a lot harder than building the current source.
The first step was figuring out how to compute the "Maximum Likelihood Estimate" or MLE. Because of the nature of the test, you can't simply average the currents of the igniters that fired.
The Likelihood is a function of the mean and variance (or its square root, standard deviation if you prefer) and the probability distribution function of the variable. In only a few cases will the MLE have a closed form solution and in this case we must resort to iteration to find it. The method is simple. Pick starting values for the mean and variance and compute the likelihood function at that point. Then add a small delta to the mean and/or variance and compute the Likelihood function again. The point with the higher Likelihood is the new starting point. Repeat until adding the delta (positive and negative) always results in a lower likelihood. That is pretty simple, what is messy is the Likelihood function. Equation 1 of Mr. Neyer's report is:
This is a nasty looking thing. A rough translation into something resembling English is: The Likelihood is a function of the mean and standard deviation and is given by the product over the samples of the combination of the number of trials (n) taken the number of successes (n * p) at a time, times the probability of success at this stimulus level ( P(z) ) raised to the power of the number of successes (n*p), times the probability of a failure ( Q(z) ) raised to the power of the number of failures (n * q).
Like I said, nasty. But as it turns out, it gets a lot simpler. In the course of running the tests, all of the stimulus levels will almost certainly be unique. This results in the value of 'n' always being 1. Also there will either be a single success or a single failure. The result is that the combination will either be 1 thing taken 1 at a time (which is 1) or 1 thing taken 0 at a time (which is also 1) so that drops out of the computation. the probabilities will be raised to powers that are either 1 or 0. Which makes that very easy as well.
One other complication, or perhaps simplification, that comes up is that we could be taking the product of lots of small numbers. For numerical reasons related to the finite precision of computer arithmetic, this is not a good thing. The normal workaround to this is to compute the log of the likelihood function. This has the effect of changing the product operator into a summation.
Now it is simple. For each sample, if it was a success, calculate the probability of success at that level and add its log to the total. If the sample is a failure, add the log of the probability of a failure at that stimulus level to the total.
That simplified nicely. Now there is the question of just what probability distribution function I should use.
Each individual igniter has some minimum current required to activate it. Because it is not possible to make every igniter absolutely identically, there is some variation in this level. Which is why you need to perform fancy tests. To make things worse, you can't just increase the applied current until the igniter functions because currents lower than that required to cause it to function will almost certainly alter its properties. Including the current required to make it function. So after each test you throw the test item away no matter the outcome.
Fear not. While the igniter is no longer useful for purposes of the test, it is still perfectly useful for igniting rocket motors. Even if it ignited during testing. It turns out that at the current levels being used the bridge wire almost never burns in two. The bare bridge wire without the pyrogen is nearly as good an igniter as one with the pyrogen. Just so long as it is touching the BP propellant, it will work. However, it would not be a good idea to use any of these igniters for clusters.
But I digress. The current required to activate any given igniter is subject to random variations in its construction. we don't know what these are nor do we much care. One of the interesting things about statistics is that if you combine several random processes into one, the combination tends towards a "normal" distribution. So we will assume that the igniter currents follow a normal distribution, with a few wrinkles.
When current is applied to the igniter for a set period of time, we are not so much applying a fixed current to the igniter as a fixed amount of energy. The temperature rise in the igniter and therefore the chance of igniting the pyrogen is a function of the energy applied. The power applied to the igniter can be computed by multiplying its resistance by the square of the current (I squared R). Power is simply energy per unit time, Joules per second in this case. Multiply by the time interval and you get the energy in Joules. The upshot of this is the chance of activating the igniter will vary with the square of the current. So the usual response is to use the log normal distribution.
Another wrinkle is that while each individual igniter's miminum current is randomly distributed and described by the log normal distribution, the chances of activating an igniter at any given current is the chance at that current plus the chances for all lower currents. So we have to use the log normal cumulative distribution function.
The code to compute this function in C is:
double probability( double x, double m, double s)
{
double z;
z = (log(x)-m)/s;
return( 0.5 * (1.0 + erf(z/M_SQRT2)));
}
In getting ready to do these tests I performed a simple check of my equipment. Mr. Neyer provides a "crippleware" version of his software for free and it is limited to eight samples and that is what I did. I used a one second (no stopwatch involved) pulse width. This was just to get a feel for how the process worked. I used igniters pulled from packages of Apogee A3-6 motors purchased at a NARAM several years ago. For the cost of the igniters I also got uncertified motors. :-)
I started with a range on the average of 0.5 to 2.0 and the guess on standard deviation as 0.05. The data for this practice run was:
| Stimulus Level (Amps) | Pass/Fail |
| 1.25 | F |
| 1.625 | P |
| 1.4375 | F |
| 1.53125 | P |
| 1.48438 | P |
| 1.39537 | P |
| 1.3045 | P |
| 1.12697 | F |
When I computed the sample mean and standard deviation using my code I got an average current of 1.34 amps with a standard deviation of 0.13 amps. the value of the log likelihood function at this point was -3.37. This data series is too short to be of much use other than seeing how the process works. I plan on a using somewhere between twenty and thirty samples for the actual tests. In order to see if there is any variation over time, I will run two all-fire tests. The first on random igniters pulled from my range box and the second on brand new, or at least newly purchased igniters.
One of the tricky things about the MLE is that you must have an overlap in the data to get valid results. In other words, the highest failure must be greater than the lowest success. Fortunately this data has a failure at 1.4375 amps and two successes at lower currents. This is yet another problem with Cody Steele's report. Three of the six sets of data have no overlap.
The next step is to figure out how Figure 2. works. One of the things holding me up at the moment is the Fisher information matrix. So far this looks really nasty with second order partial derivatives of the likelihood function but if I am lucky, things will get simpler after I understand it better. But that could take a while.
I hear someone asking "Why don't you just buy the software?". Because when someone advertises a product which has a pretty narrow market and fails to mention a price, it scares me. I am reminded of the phrase "If you have to ask the price, you can't afford it."