Referring to Fig. 6.6, a rectangular region Ш = {(x, y)a < x < b, 0 < y < c} is superimposed to enclose the area Ф = {(x, y)a < x < b, 0 < y = g(x) < c} represented by Eq. (6.49). By the hit-and-miss method, the rectangular region Ш containing the area under g(x), that is, Ф, is hung on the wall, and one is to throw n darts on it. Assume that the darts fly in a random fashion and that all n darts hit within the rectangular region. The area under g(x), then, can be estimated as the proportion of n darts hitting the target multiplied by the known area of rectangular region Ш, that is,
G = A(6.50)
where G is the estimate of the true area G under g(x), A = c(b — a) is the area of the rectangular region, and nh is the number of darts hitting the target out of a total of n trials.
The hit-and-miss method can be implemented numerically on a computer. The two coordinates (Xi, Yi) on the rectangular region Ш, which represents the location where the ith dart lands, are treated as two independent random variables that can be generated from two uniform distributions. That is, Xi is generated from U(a, b) and Yi from U(0, c). When Yi < g(Xi), the dart hits its target; otherwise, the dart misses the target. A simple hit-and-miss algorithm is given as follows:
1. Generate 2n uniform random variates from U(0, 1). Form them arbitrarily into n pairs, that is, (u1, u[), (u2, u’2),…, (un, u’n).
2. Compute xi = a + (b — a)ui and g(xi), for i = 1,2,…, n.
3. Count the number of cases nh that g(xt) > cu.
4. Estimate the integral G by Eq. (6.50).
Note that G is an estimator of the integral G; it is therefore also a random variable. It can be shown that G is unbiased, namely,
E(G) = A x = Ap = A (A = G (6.51)
where nh/n, the proportion of n darts hitting the target, is an unbiased estimator of the true probability of hits, and p simply is the ratio of the area under g(x) to the area of the rectangular region. Furthermore, the standard error associated with the estimator G is
|
|
|
|
|
|
As can be seen, the precision associated with G, represented by its inverse of standard deviation, using the hit-and-miss Monte Carlo integration method increases with n1/2.
A practical question is how many trials have to be carried out so that the estimated G satisfies a specified accuracy requirement. In other words, one would like to determine a minimum number of trials n such that the following relationship holds:
P(|G — G<e) > a (6.53)
in which є is the specified maximum error between G and G, and a is the minimum probability that G would be within є around the exact solution. Applying the Chebyshev inequality, the minimum number of trials required to achieve Eq. (6.53) can be determined as (Rubinstein, 1981)
(1 — p)p[c(b — a)]2 (1 — p)pA2
n >———— / 2——— = n * 2 (6.54)
(1 — а)є2 (1 — а)є2
Note that the required number of trials n increases as the specified error level є decreases and as the confidence level a increases. In addition, for the specified є and a, Eq. (6.54) indicates that the required number of trials n can be reduced by letting p approach 1. This implies that selecting an enclosed region Ш as close to Ф as possible would reduce the required number of trials. However, consideration must be given to the ease of generating random variates for U’ in the algorithm.
When the number of trials n is sufficiently large, the random variable T,
G — G
T = ^^ (6.55)
sG
approximately, has the standard normal distribution, that is, T ~ N(0, 1), with sG being the sample estimator of aG, that is,
|
|
|
|
|
|
Hence the (1 — 2a)-percent (a < 0.5) confidence interval for G then can be obtained as
G ± sgza
with za = Ф 1(1 — a).
Example 6.6 Suppose that the time to failure of a pump in a water distribution system follows an exponential distribution with the parameter в = 0.0008/h (i. e., 7 failures per year). The PDF of the time to failure of the pump can be expressed as
ft(t) = 0.0008e-0 0008t for t > 0
Determine the failure probability of the pump within its first 200 hours of operation by the hit-and-miss algorithm with n = 2000. Also compute the standard deviation associated with the estimated failure probability and derive the 95 percent confidence interval containing the exact failure probability.
Solution The probability that the pump would fail within 200 hours can be computed as
r-200 ,■ 200
ft (t) dt = 0.0008e—00008t dt
00 which is the area under the PDF between 0 and 200 hours (Fig. 6.7). Using the hit — and-miss Monte Carlo method, a rectangular area with a height of 0.0010 over the interval [0, 200] is imposed to contain the area representing the pump failure probability.
The area of the rectangle can be easily determined as A = 0.001(200) = 0.2. The hit-and-miss algorithm then can be outlined in the following steps:
1. Initialize i = 0 and щ = 0.
2. Let i = i + 1, and generate a pair of standard uniform random variates (ui, ui) from U(0, 1).
3. Let ti = 200ui, and compute ft(ti) = 0.0008e—00008ti, y = 0.001ui.
4. If ft(ti) > y, nh = nh + 1. If i = 2000, go to step 5; otherwise, go to step 1.
5. Estimate the pump failure probability as pf = A(nh/n) = 0.2(nh/n).
Using the preceding algorithm, 2000 simulations were made, and the estimated pump failure probability is pif = 0.2(nh/n) = 0.2(1500/2000) = 0.15. Comparing with the exact failure probability pf = 1 — exp(-0.16) = 0.147856, the estimated
|
failure probability by the hit-and-miss method, with n = 2000 and the rectangular area chosen, has a 1.45 percent error relative to the exact solution.
 |
 |
 |
 |
The associated standard error can be computed according to Eq. (6.56) as
Assuming normality for the estimated pump failure probability, the 95 percent confidence interval containing the exact failure probability pf is
pf ± Z0.975SPf = (0.1462,0.1538)
where zq.975 = 1.96.
