The Monte-Carlo method without sorting for uncertainty propagation analysis in PRA

The Monte-Carlo method without sorting for uncertainty propagation analysis in PRA

Reliability Eng~ineering 10 ~1985) 233 243 The Monte-Carlo Method Without Sorting for Uncertainty Propagation Analysis in PRA Soon Heung Chang, Joo ...

401KB Sizes 1 Downloads 58 Views

Reliability Eng~ineering 10 ~1985) 233 243

The Monte-Carlo Method Without Sorting for Uncertainty Propagation Analysis in PRA

Soon Heung Chang, Joo Young Park and Myung Ki Kim Nuclear Engineering Department, Korea Advanced Institute of Science and Technology, PO Box 150, Chongryang, Seoul, South Korea (Received: 21 April 1984)

ABSTRACT The Monte-Carlo method has been widely used Jor the analysis of an uncertainty propagation in Jault tree analysis. The conventional Monte-Carlo method is time-consuming because of the sorting process. In this stud),, the Monte-Carlo method without sorting is developed by using the segmentation between the sujficient upper bound and the sufficient lower bound of the evaluated top event jrequencies into meshes and interpolation within the mesh. It is" shown that the results obtainedJrom the new method are almost equivalent to those from the conventional method, and that the new method reduces greatly the computational running time.

1.

INTRODUCTION

Fault-tree analysis is one of the most important procedures in the safety study of the nuclear power plant. The top event in a fault tree has uncertainty because of the imperfection of a fault tree construction and the uncertainties of the basic event probabilities.1 Many methods have been developed to obtain the uncertainty of the top event probability 233 Reliability Engineering 0143-8174/85/$03-30 © Elsevier Applied Science Publishers Ltd, England, 1985. Printed in Great Britain

Soon HeungChang, Joo Young Park, Myung Ki Kim

234

Specify parameters of component distribution

rI

Generate random number for each component characteristic distribution

Combine random component values according to system function to yield value of system characteristic distribution

NO

System characteristic function of component characteristics

Sampling complete'? YES Sort distribution values

Calculate cumulative ft. and density function

Fig. 1. Flow chart of conventional Monte-Carlo method.

Monte-Carlo method without sorting in PRA

235

propagated from the uncertainties of the basic event probabilities. Wellknown methods in these analyses can be divided into three groups, as follows. (1) Methods to obtain the cumulative distribution function of the top event using the discretizations of basic event distributions. (2) Methods to obtain the top event distribution from the moments of the top event probabilities which are obtained using the moments of the basic event probabilities. (3) Methods using the Monte-Carlo technique and sorting. Each method has characteristics of its own. In the first group are the histogram method and D P D method. 2 --4 The histogram method provides the desired accuracy in an uncertainty analysis at the expense of high CPU-time requirements. The D P D method may also require comparatively high CPU-time to produce the high accuracy when there are many basic events. The moment method and the method of using fitted response surfaces belong to the second group. 5'6 In the moment method, an empirical distribution is obtained from the first and second moments of the top event. This method has been known to give comparatively reasonable results. The method of using fitted response surfaces has the problem of accuracy and requires high C PU-time in the case of many basic events. However, this method has the advantage of getting the information on which basic event is more sensitive to the uncertainty of the top event probabilities. Methods in the third group are the most powerful that can be used in any case if the basic event distribution is once assumed. The general procedures in these methods are shown in Fig. 1. However, the sorting which this group used is a very time-consuming process. Based on the review of several uncertainty propagation methods, a new method combining the three groups is developed. By using the new method, the computer code M O M O D is developed. The output and the run-time of M O M O D is compared with those of SAMPLE, which is a computer code using the conventional Monte-Carlo method and was used in WASH-1400.1 2.

METHOD DESCRIPTION

The new method developed can be described as the following six steps. First, generate the randomly-sampled probability for each basic event

236

Soon Heung Chang, Joo Young Park, Myung Ki Kim

based on the characterized distribution for the Monte-Carlo simulations. The new method described in this study can treat the basic events which are characterized by any combination of the following probability distributions: (1) (2) (3) (4) (5) (6) (7) (8)

normal, log-normal, log-uniform, uniform, exponential, gamma, Weibull and discrete distribution defined by user.

Random number generation from each basic event distribution was done by the subroutine packages in IMSL. v So this part of M O M O D can be modified according to the computer used and the number of subroutines supplied. Second, from the given relationship between the top event and basic events as the sum of minimal cut sets and sampled probabilities of the basic events, obtain the probability of the top event Qr. The first and second steps are repeated for the sufficient Monte-Carlo runs. The above steps are typical procedures which the conventional MonteCarlo method takes before the sorting process. In other words, after the second step, the new method developed takes a different procedure from that of the conventional Monte-Carlo method, sorting. Third, the mean rn and the variance V of the probability of the top event Qr are used in order to obtain the sufficient confidence bounds in Qr. Assuming that the top event has an ideal log-normal distribution for this purpose, the distribution parameters p and a and the sufficient confidence bounds, AAA and ZZZ, can be obtained from the following equations: a2 = log ( ~ 2 + 1)

(1)

0-2

/.t = log (m) + -2-

(2)

AAA = exp (/~ - 3.0a)

(3)

ZZZ = exp (/~ + 3.0a)

(4)

Monte-Carlo method without sorting in PRA

237

The values AAA and ZZZ, which are 99.73 7o confidence bounds in an ideal log-normal distribution, are used as the sufficient confidence bounds in the new method. Fourth, the distance between AAA and ZZZ is divided into the meshes. There are two choices in the segmentation. One method is a linear scale segmentation with equal distances, and the other is a log-scale segmentation. Monte-Carlo simulations for obtaining the probabilities of the top event

Obtain the mean m and the variance V of the probabilities of the top event

Obtain the sufficient confidence bounds AAA, ZZZ using the values of m and V

Segment the distance between AAA and ZZZ into meshes

Determine the number of probabilities of the top event in each mesh

I Obtain the cumulative distribution function of the top event and plot Fig. 2.

Flow chart of the method developed.

Soon Heung Chang, Joo Young Park, Myung Ki Kim

238

Fifth, after segmentation the mesh that each obtained top event probability belongs to is determined, and the number of top event probabilities obtained in each mesh is counted. Finally, from the information obtained from the above step, an approximate cumulative distribution function of the top event can be obtained. The method of obtaining the cumulative distribution function of the top event is interpolation. After all these procedures, the relative frequency and the cumulative distribution function of the top event is plotted. The procedures in M O M O D which are described above are shown schematically in Fig. 2.

3.

S A M P L E C A L C U L A T I O N S A N D RESULTS

The code developed, M O M O D , was tested in four cases. In the first case, the basic events were characterized by the log-normal distribution, which is most widely used in nuclear power plant failure data. 1 This case dealt with the fault tree in Fig. 3 and the failure data in Table 1. The cumulative distribution functions of the top event obtained from S A M P L E and M O M O D , and the computer running time, are shown in Table 2. In the second case, with the same fault-tree used in the above case, the case when the failure data were characterized by the normal distributions in Table 3 was considered. The results are shown in Table 4. T System Failure

+ Fig. 3.

+ Fault tree.

TABLE 1 Failure Data

Fault

Distribution

X(I) X(2) X(3) X(4) X(5) X(6)

log-normal log-normal log-normal log-normal log-normal log-normal

Median 1.0 3.0 1-0 3.0 1.0 3.0

x x x x x x

Error factor

10- 3 10 2 10 2 10- 2 10 : 10 3

3 3 3 3 3 3

C o m m o n Mode Xc~~= )((3). X(5) QT= X(l) +

(X(2) + X(3)). (X(4) + X(5)) + X(6) +,/X(3)5. Xi5).

TABLE 2 Results of Test Case 1

Code

SAMPLE

SAMPLE

MOMOD (linear)+

MOMOD (linear)

MOMOD (log)t+

MOMOD (log)

imax*

1200

3 600

1 200

3 600

1 200

3 600

100

100

100

isegt

100 _

_

o

Percentile l 2.5 5 10 20 30 40 50 60 70 80 90 95 97.5 99 C PU-time (s) * t + tt

3.190E-3 3.868E-3 4-255E-3 4.954E-3 5.794E-3 6.602E-3 7.338E-3 8.034E-3 8.912E-3 1-003E-2 1.169E-2 1.385E-2 1.655E-2 1.847E-2 2-150E-2

3.154E-3 3.683E-3 4.154E-3 4-778E-3 5.704E-3 6.475E-3 7.192E-3 8.023E-3 8.936E-3 1.002E-2 1.149E-2 1.395E-2 1.638E-2 1.847E-2 2.171E-2

3-212E-3 3.632E-3 4.043E-3 4.707E-3 5.610E-3 6.391E-3 7.150E-3 7.864E-3 8-771E-3 9.938E-3 1.127E-2 1.369E-2 1.720E-2 2.032E-2 2.455E-2

3.189E-3 3.619E-3 4.085E-3 4.761E-3 5.670E-3 6.450E-3 7.218E-3 8.000E-3 8.966E-3 1.012E-2 1-157E-2 1-408E-2 1.698E-2 1.993E-2 2.402E-2

3.189E-3 3.642E-3 4.004E-3 4.727E-3 5.614E-3 6.385E-3 7.518E-3 7.858E-3 8-765E-3 9-934E-3 1.124E-2 1.373E-2 1-721E-2 2.043E-2 2.440E-2

3.177E-3 3.623E-3 4.086E-3 4.759E-3 5.665E-3 6.455E-3 7.218E-3 8.014E-3 8.960E-3 1.013E-2 1.157E-2 1.398E-2 1-698E-2 1.994E-2 2.390E-2

18.31

98.84

5.49

13.67

5.63

14.14

imax = number of Monte-Carlo runs. iseg = number of meshes. M O M O D (linear) = M O M O D with linear segmentation. M O M O D (log) = M O M O D with log-scale segmentation.

Soon HeungChang, Joo Young Park, M y u n g K i K i m

240

TABLE 3 Failure Data

Basic event

Mean

90 '"o Error spread

X(I)

0-001

0.000 3

X(2)

0.03

0.01

Y(3)

0.01

0.003

X(4)

0.03

0.01

X(5) X(6)

0.01 0.003

0.003 0.001

Qr=Y(l)

+ (X(2) + X ( 3 ) ) . (X(4) + X(5)) + X(6) + .v/X(3)2. X(5).

TABLE 4 Results of Test Case 2

Percentile

SAMPLE

MOMOD (log)

1

4.901E-

5

5.382E - 3

3

4.90lE.-

5.400E - 3

MOMOD (linear)

3

4-899E-

3

5.404E - 3

10

5.646E - 3

5.633E - 3

5-637E -- 3

50

6.605E - 3

6.583E - 3

6.580E - 3

90

7-606E - 3

7-60 4E - 3

7.609E - 3

95

7.942E - 3

7.873E - 3

7.865E - 3

99

8.412E - 3

8.393E - 3

8.402E - 3

17- 83

4.15

4.00

C P U - t i m e (s)

TABLE 5 Results of Test Case 3

Uncertainty method MOMOD

(log)

MOMOD SAMPLE

(linear)

Moment method DPD

5 '!;,

Median

95 i'~,

2'666E - 3 2 . 6 4 2 E -- 3 2.581E- 3

5.222E - 3 5.219E - 3 5.350E- 3

1.167E - 3 1.159E - 3 1.171E- 3

2.67E - 3 2.47E - 3

5.27E - 3 5.30E - 3

t.32E - 3 1' 46E - 3

Q r = X1 + X2. X 3 q- X 4 . X 5 -~- X 6. E r r o r f a c t o r = 3.

241

Monte-Carlo method without sorting in P R A TABLE 6 Results of Test Case 3 Uncertainty method

M O M O D (log) M O M O D (linear) SAMPLE Moment method DPD

5 04,

Median

95 ~~

1.832E - 3 1.042E - 3 1.776E - 3 8.36E - 4 1.87E - 3

8.143E - 3 8.187E - 3 8.530E - 3 6.07E - 3 8.75E - 3

4-303E - 2 4.350E - 2 4.533E - 2 5-36E - 2 4-542E - 2

QT = XI + X2 • )(3 nt-X4- X5 -b X6. Error factor = 10.

I n t h e t h i r d case, u s i n g the p r e v i o u s s e n s i t i v i t y s t u d y d o n e ip. ref. 4, MOMOD

w a s c o m p a r e d w i t h the o t h e r m e t h o d s . T h e r e s u l t s a r e s h o w n

in T a b l e s 5 a n d 6. F i n a l l y , t h e case w h e n t h e b a s i c e v e n t s a r e n o t c h a r a c t e r i z e d b y t h e s a m e k i n d o f d i s t r i b u t i o n w e r e c o n s i d e r e d w i t h the s a m e f a u l t tree as t h e first case a n d t h e f a i l u r e d a t a i n T a b l e 7. T h e r e s u l t s a r e s h o w n in T a b l e 8 a n d c o m p a r e d w i t h t h o s e o f the m o d i f i e d S A M P L E . S A M P L E c a n n o t t r e a t t h e case w h e n t h e b a s i c e v e n t s a r e n o t c h a r a c t e r i z e d b y t h e s a m e k i n d o f d i s t r i b u t i o n , so S A M P L E

w a s m o d i f i e d for this p r o c e d u r e .

TABLE 7 Failure Data Fault

Distribution

X(1) X(2) X(3) X(4) X(5) X(6)

normal normal gamma normal gamma normal

Distribution parameter

mean = mean = A* = mean = A* = mean =

0.001, 90 '!o error 0.03, 90130error 3.0, B* = 0'000 4 0.03, 90"; error 3"0, B* = 0"0004 0.003, 90 ?o error

Common mode XcM = X(3). X(5). * Gamma distribution is given as follows f X a ~ I e -X.B

J(x)=IoBgF~A)

elsewhere.X>0(A>0'B>0)'

spread = 0.000 3 spread = 0 0 t spread = 0-01 spread = 0-001

242

Soon Heung Chang, Joo Young Park, Myung Ki Kim TABLE 8 R e s u l t s o f Test C a s e 4

Percentile

Modified SAMPLE

MO MOD (log)

MO MOD ( linea r)

0.349 3E - 02 0-395 7E - - 0 2 0.421 0E - 02 0.501 7F 02 05924E02 06176E02

0-349 3E - 02 0'395 5E - 02 0"420 6E - 02 0"501 4E - 02 0.593 I E - 0 2 0.6173E -02

0.6654E-

0.6622E-02

+

1 5 10 50 90 95

0"348 0-395 0.420 0.501 0-593 0617

99

0.6623E

C P U - t i m e (s)

3E - 02 8E - 02 4E - 02 8E -- 02 5E-02 IE--02

1533

+02

02

4+25

4.09

In the results of the first case, it was shown that the output of M O M O D is almost equal to that of S A M P L E when the basic events were characterized by the log-normal distributions, and that the computer running time required in M O M O D was much less than that of SAM PLE. This was true especially when the larger Monte-Carlo runs were used. F r o m the results of the second case, it was shown that the output of M O M O D was reasonable when another distribution (such as normal distribution versus log-normal distribution) was employed for the failure data. In the results of the third case, it was shown that M O M O D gave results compatible with those of other methods for both the narrow spread and the large spread. It is also seen from Table 6 that M O M O D with log-scale segmentation gave the more reasonable results for the large spread compared with S A M P L E than did M O M O D with linear scale segmentation. Finally, in the results of the fourth case, it was shown that the output of M O M O D was also sufficiently accurate and effective in the computer running time saving when the basic event was characterized by the combination of the various distributions. 4.

CONCLUDING REMARKS

For the uncertainty propagations of probabilistic risk assessment+ the Monte-Carlo method without sorting is developed by using the

Monte-Carlo method without sorting in PRA

243

segmentation between the sufficient upper bound and the sufficient lower bound of the evaluated top event frequencies into meshes and interpolation within the mesh. The method developed was tested when the basic events were characterized by log-normal distribution, normal distribution and the combination of normal distribution and gamma distribution. The method developed gives results equivalent to those of the conventional Monte-Carlo method and other methods of treating the uncertainty propagation in all of test cases. The new method reduces greatly the computational running time, especially when larger MonteCarlo runs were used.

ACKNOWLEDGEMENT The financial support of ASAN acknowledged.

Foundation

for this research is

REFERENCES 1. US NRC WASH-1400 (NUREG 75/014), October, 1975. 2. Colombo, A. G. and Jaarsma, R. J. A Powerful Numerical Method to Combine Random Variables, 1EEE Trans. Reliab., R-29(2) (1980), pp. 126 9. 3. Ahmed, S., Clark, R. E. and Metcalf, D. R. A Method for Propagating Uncertainty in PRA, Nucl. Technol., 59 (1982), pp. 238-45. 4. Ahmed, S., Metcalf, D. R. and Pergram, J. W. Uncertainty Propagation in PRA: A Comparative Study, Nucl. Engng Design, 68 (1981), pp. 1--3, 5. Apostolakis, G. and Lee, Y. T. Methods for the Estimation of Confidence Bounds for the Top Event Availability of Fault Tree, Nucl. Engng Design, 41 (1977), pp. 411-19. 6. Mazumdar, M. An Approximate Method for Computation of Probability Intervals for the Top Event Probability of Fault Trees, Nucl. Engng Design, 71 (1982), pp.45 50. 7. International Mathematical & Statistical Library Inc., The I MSL Library, 1981, 8. Cox, D. C., Baybutt, P. and Kurth, R. E. Comments on the Uncertainty in Accident Consequences Calculated by Large Code due to Uncertainty in Input, Nuel. Technol., 52 (1981), pp. 439-41. 9. Jackson, P. S., Hockenbury, R. W. and Yeater, M. g. Uncertainty Analysis of System Reliability and Availability Assessment, Nucl. Engng Design, 68 (1981), pp. 5 29.