© Strojni{ki vestnik 46(2000)1,14-23 ISSN 0039-2480 UDK 519.2:621.039.58 Pregledni znanstveni ~lanek (1.02) © Journal of Mechanical Engineering 46(2000)1,14-23 ISSN 0039-2480 UDC 519.2:621.039.58 Review scientific paper (1.02) Dolo~anje odzivne povr{ine z optimalnim statisti~nim cenilnikom Response Surface Generation with Optimal Statistical Estimator Andrej Pro{ek - Borut Mavko Na področju jedrske tehnike se odzivna površina uporablja za reševanje problemov, povezanih z jedrsko varnostjo. Glavni namen te študije je bil zasnovati orodje za avtomatsko določanje odzivne površine zapletenih in nelinearnih pojavov ter ga preskusiti pri določanju odzivne površine za najvišjo temperaturo srajčke med malo izlivno nezgodo v jedrski elektrarni. Za določanje odzivne površine smo uporabili optimalni statistični cenilnik, ki smo ga priredili za uporabo v večdimenzionalnem prostoru. V postopek smo vgradili dva statistična kazalca, ki povesta, kako točno smo napovedali posamezne točke in možnost nastavljanja širine Gaussove krivulje. Za preskus delovanja optimalnega statističnega cenilnika smo uporabili rezultate male izlivne nezgode 59 različnih primerov. Uporaba optimalnega statističnega cenilca za določanje odzivne površine je pokazala vrsto prednosti pred regresijsko analizo, ki se večinoma uporablja v svetu. Rezultati so pokazali, da z optimalnim statističnim cenilcem lahko dovolj natančno napovemo odzivno površino za najvišjo temperaturo srajčke, kar z regresijsko analizo ni bilo mogoče. Še več, z avtomatiziranim določanjem odzivne površine se je odprla široka možnost uporabe optimalnega statističnega cenilnika za oceno negotovosti poljubne vrste časovno odvisnih pojavov in nezgod. © 2000 Strojniški vestnik. Vse pravice pridržane. (Ključne besede: površine odzivne, optimalni statistični cenilec, najvišja temperatura srajčke, varnost jedrska) In the field of nuclear engineering the response surface is used to solve some problems related to nuclear safety. The main purpose of the study was to develop a tool suitable for response surface generation of complex and non-linear phenomena and to demonstrate its applicability for the response surface generation of peak cladding temperature during a small-break loss-ofcoolant accident in a nuclear power plant. The optimal statistical estimator, adapted for use in multi-dimensional space, was used for response surface generation. For assessing the adequacy and predictive capability of the optimal statistical estimator two statistics were built in, and the possibility to set the width of the Gaussian curve. The performance of the optimal statistical estimator was tested with the results from 59 different calculations of the small-break loss-of-coolant accident. The application of the optimal statistical estimator shows several advantages when compared to the more commonly used regression analysis. The results showed that the response surface for the peak cladding temperature was adequately predicted by the optimal statistical estimator but not with regression analysis. Furthermore, an ability to automate the response surface generation provides the possibility of using the optimal statistical estimator for an uncertainty evaluation of any kind of time dependent phenomena and transients. © 2000 Journal of Mechanical Engineering. All rights reserved. (Keywords: response surface, optimal statistical estimator, peak cladding temperature, nuclear safety) 0 UVOD Na področju jedrske tehnike se odzivna površina uporablja za reševanje različnih problemov, povezanih z jedrsko varnostjo. Odzivna površina nadomesti termohidravlični računalniški program, ko za statistično analizo potrebujemo na tisoče izračunov. Običajno so odzivno površino računali za eno samo vrednost, npr. najvišjo temperaturo srajčke 0 INTRODUCTION In the field of nuclear engineering the re-sponse surface is used to solve various problems re-lated to nuclear safety. When thousands of complex computer code runs are needed for statistical analy-sis, the response surface is used to replace the computer code. Usually the response surface is generated for a single value parameter, for example: peak grin^SfcflMISDSD VBgfFMK stran 14 A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation goriva, najnižji nivo v sredici reaktorja ali največji tlak sistema. Uporaba odzivne površine na področju jedrske tehnike se v svetu zvečuje. Že leta 1989 so v ZDA za izračun negotovosti najvišje temperature srajčke med veliko izlivno nezgodo uporabili odzivno površino [1]. Za določanje odzivne površine so uporabili regresijsko analizo (prilagajanje s polinomi). Leta 1990 so, ponovno v ZDA, uporabili odzivno površino za določitev nastavitvene točke, pri kateri se odpirajo varnostni ventili na sekundarni strani jedrske elektrarne. Za prilagajanje z regresijsko analizo so uporabili linearne in kvadratne člene. Novo uporabo odzivne površine je zaslediti v letu 1992 za veliko in malo izlivno nezgodo izlivno nezgodo ([2] do [4]). Leta 1996 je bila odzivna površina uporabljena za določitev negotovosti najvišje temperature srajčke med veliko izlivno nezgodo [5]. Odzivno površino uporablja za svoje preračune negotovosti najvišje temperature srajčke od leta 1996 dalje tudi podjetje Westinghouse [6]. Čeprav je že bilo opravljenega dosti dela, je vsem študijam skupno to, da lahko uporabijo odzivno površino le za eno vrednost parametra, ne pa tudi za časovni potek. Glavni namen našega dela je bil samodejno določiti odzivno površino za poljubne diskretne točke, ki popisujejo časovni potek parametra. To je pomembno pri uporabi postopkov, ki za ocenjevanje negotovosti uporabljajo odzivno površino, saj regresijska analiza ni uporabna za vse vrste parametrov in s tem scenarijev, za katere bi želeli oceniti negotovost. 1 DOLOČANJE ODZIVNE POVRŠINE Da bi lahko statistično sklepali o negotovosti rezultatov, sama analiza občutljivosti ni zadostna, ker je število izračunov premajhno. Ker zaradi omejenih računalniških zmogljivosti ni mogoče izvesti več tisoč izračunov, se iz določenega števila izračunov določi odzivno površino in s postopkom Monte Carlo naključno spreminja vhodne parametre. Na podlagi več deset tisoč vrednosti, dobljenimi z odzivno površino, se da statistično sklepati o negotovosti rezultata. Do zdaj se je v svetu za določanje odzivne površine najpogosteje uporabljala regresijska metoda ([1] do [5]). Regresijska metoda je primerna predvsem za opis pojavov, pri katerih je odvisnost pretežno linearna. Stopnja polinomov določa število potrebnih preračunov s termohidravličnim programom, samo določanje odzivne površine je zamudno in neprimerno za računalniško avtomatizacijo. Zaradi naštetih omejitev, predvsem nezmožnosti obravnavanja močno nelinearnih pojavov, regresijska metoda npr. ni bila primerna za uporabo na mali izlivni nezgodi. Poskusili smo z linearno interpolacijo, ki je bila uporabna le za dve cladding temperature, lowest reactor core level or peak system pressure. The response surface is increasingly used in the field of nuclear engineering. In 1989 the response surface was used in the USA for the uncertainty evaluation of the peak cladding temperature [1]. Regression analysis (polynomial fit) was used for the response surface generation. In 1990, again in the USA, the re-sponse surface was developed for the determination of high-pressure setpoints for the safety valves on the secondary side of a nuclear power plant. Linear and quadratic terms were used for a polynomial fit. In 1992, other applications of the response surface were in a large-break and a small-break loss-of-coolant accident ([2] to [4]). In 1996, the response surface was used to evaluate uncertainties in the peak cladding temperature during a large-break loss-of-coolant accident. Westinghouse have also begun to use the response surface to evaluate uncertainties in the peak cladding temperature since 1996. Although much work has already been done, in all the studies so far the response surfaces devel-oped could only be used for single-value parameters, with no possibility for continuous-valued parameters. The main purpose of this work was to automate the response surface generation for any set of discrete points, characterising the time trend of the parameter. This is very important for application of uncer-tainty methods, which use the response surface, be-cause regression analysis is not applicable to all kinds of parameters and scenaria to be evaluated for un-certainty. 1 RESPONSE SURFACE GENERATION In order to statistically quantify the uncer-tainties, the number of calculations in the sensitivity analysis is too small. Because computer capabilities introduce a limit of the thousand calculations, the practice is to develop the response surface from se-lected calculations and with the Monte Carlo method, randomly select input parameters. The uncertainty is evaluated on the basis of ten thousand values, pre-dicted by the response surface. In the past, the re-sponse surface was generated mostly with regression analysis ([1] to [5]). Regression analysis is most suitable for a description of parameters with mostly linear behav-iour. The order of the polynomial determines the number of calculations to be performed with the thermalhydraulic code, the response surface generation is time consuming and not adequate for computer automation. Because of these limitations, especially the capability of the response surface to describe highly non-linear phenomena, regression analysis was not applicable to a small-break loss-of-coolant accident. Therefore, we tried with a linear interpolation method | gfin=i(gurMini5nLn 00-1_____ stran 15 I^BSSIfTMlGC A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation neodvisni spremenljivki [7]. Linearna interpolacija torej ni rešitev za večje število parametrov. Da bi lahko določili negotovost izračuna male izlivne nezgode, smo zato potrebovali novo orodje. Postopek optimalnega statističnega cenilnika, uporabljen že leta 1992 [2] se je izkazal za primernega, če ga posplošimo za uporabo v večdimenzionalnem prostoru. Postopek optimalnega statističnega cenilnika namreč ni odvisen od števila izhodnih parametrov, uporablja zelo nelinearne funkcije in je preprost za računalniško uporabo. 2 POSTOPEK OPTIMALNEGA STATISTIČNEGA CENILNIKA Za določitev odzivne površine smo uporabili postopek optimalnega statističnega cenilca (OSC), ki je bil formuliran in uporabljen pri modeliranju ultrazvoka [8]. Na področju jedrske tehnike smo postopek prvič uporabili leta 1992 za primerjavo rezultatov regresijske analize. Uporabili smo osnovni niz enačb iz [8]. V letu 1998 smo jih dodatno posplošili za večdimenzijski primer uporabe ([9] in [10]). Odzivno površino določimo iz znanih izračunanih ali merjenih vrednosti. Predpostavimo, da parameter, ki ga želimo nadomestiti z odzivno površino, lahko opišemo s točkami x, ki sestavljajo popoln vektor X=(x1, x2,...,x). Ta vektor lahko sestavimo iz dveh delnih vektorjev G=(x1, x2,...,xM; 0) in H=(0; xM+1, xM+2,...,x). Simbol 0 označuje podatke, ki niso podani, I pa število vhodnih in izhodnih točk. Popoln vektor lahko zapišemo kot sestav: using two independent input parameters [7]. However, the linear interpolation does not have a solution when we have more than two input parameters. The need for a new tool was identified for the uncertainty evaluation of small-break loss-of-cool-ant results. The optimal statistical estimator developed and applied in 1992 seems to be correct if adapted for use in the multidimensional space. Namely, the optimal statistical estimator is independent of the number of output parameters, it uses highly non-linear functions and is simple for a computer application. 2 OPTIMAL STATISTICAL ESTIMATOR METHOD For the response surface generation the optimal statistical estimator (OSC) that was formulated and used in the modelling of ultrasonic data [8], was used. In the nuclear field the method was for the first time adopted in 1992 for the comparison of OSC with regression analysis. The basic equations derived in [8] were used. In 1998, the method was further improved for multi-dimensional space ([9] and [10]). The response surface is predicted from the calculated or the measured values. In the following we assume that the parameter under consideration can be characterised by a finite number of data x, which are represented by the complete data vector X=(x1, x2,...,x). The complete data vector is often composed of two partial data vectors G=(x1, x2,...,xM; 0) and H=( 0; xM+1, xM+2,...,xI). Here the symbol 0 indicates the corresponding data, which are not specified, and I the number of input and output parameters. The complete data vector can then be expressed by a composition: X = G®H = (x1,x2,...,xM,xM+1,x x M +2, ..., xI ) (1). V našem primeru vektor G pomeni M vhodnih vrednosti (npr. vhodni negotovi parametri), H pa (I-M) izhodnih vrednosti (npr. najvišjo temperaturo srajčke in največji tlak sistema). Za optimalni statistični cenilnik HO so v [8] izpeljali izraz, ki predstavlja linearno kombinacijo izračunanih ali izmerjenih vrednosti H in koeficientov C: In our case, the vector G represents M input data points (for example, the values of input parameters) and H represents (I-M) output data points (for example, the peak cladding temperature and peak system pressure). The optimal statistical estimator HO , derived in ref. [8], is expressed as a linear combination of given values H and coefficients C : kjer je: H0(G) = N CnHn (2), where: Hn=(xn(M+1),xn(M+2),...,xnI) N število izračunanih ali merjenih vrednosti in: N is the number of calculated or measured values and: C da(G -Gn) N Gn) (3). grin^SfcflMISDSD VBgfFMK stran 16 A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation Koeficienti C so merilo podobnosti med danim vektorjem G in vektorjem vhodnih parametrov G za n-ti izračun. Pri tem je približek funkcije d Gaussova funkcija: The coefficient C represents a measure of the similarity between a given vector of input data G and the vector of input data G, for n-th calculation. The approxi-mation of the d function is the Gaussian function: da ( G-Gn ) = n 1 2p ¦ s exp 1 M 2 2 s in and Gn=(xn1,xn2,...,xn ) (4) (5). Sedaj se pojavi problem, kako najbolje izbrati širino si. Glavni namen glajenja je raztegniti vpliv posameznih merjenih točk v njihovo okolico, npr. do sosednjih točk. Če želimo z N vzorci pokriti prostor približno enakomerno, potem moramo za si v i-ti dimenziji izbrati : si= Si f, i kjer je Si razlika med največjo in najmanjšo vrednostjo iz množice vhodnih podatkov xni (n=1, 2,...,N) za i-to dimenzijo, Ni pa je število korakov med točkami v isti dimenziji. Faktor f je faktor za širino Gaussove krivulje, ki ga izberemo na podlagi prej opravljenega testa optimalnega statističnega cenilnika. S faktorjem za širino Gaussove krivulje se nastavi prispevek posameznih točk h končnemu rezultatu. Z uporabo izpeljanega optimalnega statističnega cenilnika H0 lahko končni vektor zapišemo kot: The problem now appears as how best to select the width si. The main purpose of smoothing is to stretch the influence of a particular input data point into its surroundings, e.g., approximately to its neighbours. If we want to cover the volume by sam-ples uniformly, we define si for the dimension i: i = 1,2,...,M (6), where Si is the distance between the minimum and maximum value of the set of input data points xni (n=1, 2,...,N) in the i-th dimension and Ni is the number of intervals between data points in the i-th dimension. The factor f, for the width of the Gaussian curve, is to be selected by the user, based on the desired and previously tested performance of the optimal statistical estimator. The contribution of each data point to the final output parameter estimation can be adjusted by this factor. Using the derived optimal statistical estimator H0 , the complete estimated vector can be defined as: y = (g®h0 Y(X) (7). Funkcijo Y(X) je treba modelirati z računalniškim programom. Na vektor Y vplivajo vhodne vrednosti, ki jih neposredno priredimo izhodu, medtem ko komplementarne vrednosti določimo z optimalnim statističnim cenilnikom H0(G) . Ena izmed najpomembnejših lastnosti tega cenilca je, da vsebuje zelo nelinearne funkcije C. Dobljeni sistem je tako v bistvu nelinearen, čeprav je vektor H0 izražen kot linearna kombinacija vrednosti H in koeficientov C . Za natančnost prileganja odzivne površine izračunanim ali izmerjenim točkam smo uporabili srednji kvadratni pogrešek za m-ti parameter: N The function Y(X) is modelled by the computer code. Vector Y is influenced by input parameters, directly transferred to the output, while the complementary values are determined by the optimal statistical estimator H0(G) . One of the most important characteristics of this estimator is that involves the highly non-linear function C . n For assessing the adequacy and predictive capability of the optimal statistical estimator, the root mean square error for m-th parameter was used: RMSm = Z(x )2 in koeficient določitve R2: N R2m Y(x -x )2 / j \ estn ,m avg,m / N / ,\xnm n=1 x avg,m )2 ; m=(M+1), (M+2),...,I and coefficient of determination R2: ; m=(M+1), (M+2),...,I (8) (9). Sin^ObJJPsflDslJSD I stran 17 glTMDDC A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation Pri tem je xm n-ta, s programom izračunana vrednost za m-ti izhodni parameter, xest m n-ta vrednost, ocenjena z optimalnim statističnim cenilnikom za m-ti izhodni parameter in xavgm povprečna vrednost, s programom izračunanih N točk za m-ti izhodni parameter. Ujemanje med izračunanimi vrednostmi in vrednostmi, ki jih napove optimalni statistični cenilnik, je najboljše, ko se srednji kvadratni pogrešek manjša proti vrednosti RMS = 0 in ko se koeficient določitve bliža vrednosti '= 1. Da bi Rm zračunali vrednost izhodnih parametrov, vrednost vhodnih parametrov (x1,x2,...,xM) vsakič naključno spreminjamo (ali jih vnesemo) in potem z optimalnim statističnim cenilnikom ocenimo izhodno vrednost z uporabo enačb ((2) do (5)). Pri vsaki oceni izračunamo nove vrednosti funkcij C, medtem ko so vrednosti H izračunane točke. Če primerjamo OSC z regresijsko analizo, vidimo, da OSC ni odvisen od števila podanih točk za opis odzivne površine, medtem ko je pri regresijski analizi to število točk odvisno od reda polinoma. Če primerjamo naravo pojavov, potem je prednost OSC pred regresijsko analizo v tem, da z OSC lahko popišemo nelinearne in zapletene funkcijske odvisnosti, z regresijsko analizo pa ne. Nenazadnje je OSC algoritem zelo primeren za računalniško uporabo, pri regresijski analizi pa se običajno uporabljajo statistični paketi. 3 IZRAČUN NAJVIŠJE TEMPERATURE SRAJČKE ZA MALO IZLIVNO NEZGODO Da bi pokazali delovanje postopka OSC, smo uporabili rezultate analize občutljivosti najvišje temperature srajčke na izbrane vhodne parametre med malo izlivno nezgodo [11]. Za analizo male izlivne nezgode je bil uporabljen termohidravlični računalniški program RELAP5/MOD3.2. Uporabljeni vhodni model za program RELAP5/ MOD3.2 je bil model dvozančne tlačnovodne jedrske elektrarne. Ko smo si izbrali scenarij in elektrarno, smo identificirali pomembne pojave, ki vplivajo na najvišjo temperaturo srajčke in jih razvrstili po pomembnosti za jedrsko varnost. Pojavom smo pripisali parametre, s katerimi se jih da popisati v računalniškem programu RELAP5. Izbrani parametri (normirani) so bili enofazni iztočni koeficient (SDC), dvofazni iztočni koeficient (TPDC), toplotna prestopnost (HTC), medfazno trenje (IDC) in cepitveni delež zaostale toplote (FPYF). Te parametre smo v analizi občutljivosti spreminjali in dobili različne izračunane vrednosti najvišje temperature srajčke (PCTRELAP5), prikazane v preglednici 1. Za te izračunane najvišje temperature srajčke smo določili odzivno površino, ki smo jo potrebovali za statistično določitev negotovosti z uporabo postopka Monte Carlo. Here xnm is the n-th code calculated value of the m-th output parameter, xestn,m is n-th estimated value with the optimal statistical estimator and xavg,m is the mean of the N code calculated values of the m-th output parameter. The predictive capability of the optimal statistical estimator, assessing with the two proposed statistics, is perfect when RMSm=0 and Rm2 = 1. To produce output results the values of the input parameters (x1,x2,...,xM) were randomly sampled (or input by the user) each time and then the corre-sponding unknown output values were estimated by the optimal statistical estimator using Eqs. ((2) to (5)). Each time a new coefficient Cn is calculated, while the values of Hn are calculated points obtained by computer code. When comparing the OSC to the regression analysis, we see that for the regression analy-sis with the polynomial the amount of data is pre-scribed to describe the response surface while in the OSE more data mean more information that can be extracted. When comparing the nature of the phenomena, the advantage of OSC with respect to the regression analysis, is in its ability to predict very complex and highly non-linear functions. Fi-nally, the algorithm for OSC is suitable for computer automation, while for regression analysis statis-tical packages are used. 3 CALCULATION OF PEAK CLADDING TEMPERATURE FOR SMALL-BREAK LOSS-OF-COOLANT ACCIDENT To demonstrate the OSC method, the results of the sensitivity analysis of the peak cladding temperature on selected input parameters during the small-break loss-of-coolant accident were used [11]. For the small-break loss-of-coolant accident analy-sis the RELAP5/MOD3.2 thermal-hydraulic computer code was used. The input model for the RELAP5/ MOD3.2 code was a model of a two-loop pressu-rised water reactor. After the scenario selection, all the important phenomena influencing on the peak cladding temperature were identified and ranked by their importance in nuclear safety. Then key RELAP5 code parameters were selected to represent the im-portant phenomena. The selected parameters, all normalised, were subcooled discharge coefficient (SDC), two-phase discharge coefficient (TPDC), heat transfer coefficient (HTC), interphase drag coefficient (IDC), and fission product yield factor (FPYF). In the sensitivity analysis these parameters were varied to obtain the corresponding calculated peak cladding temperatures (PCTRELAP5), which are shown in Table 1. For the calculated peak cladding temperatures the response surface was generated and sampled with the Monte Carlo method to generate an approximate distribution that characterises the uncertainty. grin^SfcflMISDSD VH^tTPsDDIK stran 18 A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation Preglednica 1. Izračunana in napovedana najvišja temperatura srajčke v odvisnosti od vhodnih parametrov Table 1. Calculated and predicted peak cladding temperature as a function of input parameters Izračun PCTrelap5 PCTosc Calculation SDC TPDC HTC IDC FPYF (K) (K) 1 0,917 1,435 1 1 1 1073 1073,1 2 0,833 1,436 1 1 1 1028 1028,8 3 1,001 1,437 1 1 1 1154 1146,3 4 0,917 1 1 1 1 1176 1175,9 5 0,917 1,311 1 1 1 1093 1092,4 6 0,917 1,559 1 1 1 1056 1056,5 7 0,917 1,435 0,75 1 1 1117 1117 8 0,917 1,435 1,25 1 1 1098 1098 9 0,917 1,435 1 0,9174 1 1073 1073 10 0,917 1,435 1 1,0826 1 1067 1067 11 0,917 1,435 1 1 0,9 1056 1056 12 0,917 1,435 1 1 1,1 1054 1054 13 1,001 1 1 1 1 1058 1058 14 0,833 1 1 1 1 1083 1083 15 1,001 1,311 1 1 1 1025 1028,6 16 1,001 1,559 1 1 1 1004 1008,2 17 1,001 1 0,75 1 1 1066 1066 18 1,001 1 1,25 1 1 1032 1032 19 1,001 1 1 0,9174 1 1058 1058 20 1,001 1 1 1,0826 1 1062 1062 21 1,001 1 1 1 0,9 1013 1013 22 1,001 1 1 1 1,1 1076 1076 23 1,001 1 0,75 1 0,9 1082 1082 24 0,917 1,435 0,75 1 1,1 1070 1072,1 25 0,917 1,435 1,25 1 0,9 1012 1016,2 26 0,917 1,435 1,25 1 1,1 1124 1123,2 27 0,917 1 0,75 1 1,1 1129 1129 28 0,833 1,435 1 1 0,9 1030 1030 29 0,833 1,311 1 1 1 1055 1054,3 30 0,833 1,435 1,25 1 1,1 1233 1233 31 0,833 1,435 1,25 1 0,9 1015 1015 32 0,833 1,435 0,75 1 0,9 1101 1101 33 0,917 1 0,75 1 0,9 1048 1048 34 0,917 1 1,25 1 0,9 1040 1040 35 0,917 1 1,25 1 1,1 1095 1095 36 0,917 1,311 0,75 1 0,9 1068 1068 37 0,917 1,311 0,75 1 1,1 1143 1140,9 38 0,917 1,311 1,25 1 0,9 1159 1154,8 39 0,917 1,311 1,25 1 1,1 1096 1096,8 40 1,001 1,311 0,75 1 0,9 982,1 983,5 41 1,001 1,311 0,75 1 1,1 1092 1092,1 42 1,001 1,311 1,25 1 0,9 974,9 975,1 43 1,001 1,311 1,25 1 1,1 1054 1055 44 1,001 1,435 0,75 1 0,9 1030 1028,7 45 1,001 1,435 0,75 1 1,1 1096 1095,9 46 1,001 1,435 1,25 1 0,9 981,3 981,1 47 1,001 1,435 1,25 1 1,1 1090 1089 48 0,917 1 0,75 0,9174 1,1 1129 1129 49 0,833 1,435 1 0,9174 0,9 1030 1030 50 0,833 1,435 0,75 0,9174 0,9 1101 1101 51 0,917 1 0,75 0,9174 0,9 1048 1048 52 0,917 1,311 0,75 0,9174 1,1 1143 1143 53 1,001 1,311 0,75 0,9174 1,1 1092 1092 54 0,917 1 0,75 1,0826 1,1 1142 1142 55 0,833 1,435 1 1,0826 0,9 1096 1096 56 0,833 1,435 0,75 1,0826 0,9 1089 1089 57 0,917 1 0,75 1,0826 0,9 1048 1048 58 0,917 1,311 0,75 1,0826 1,1 1110 1110 59 1,001 1,311 0,75 1,0826 1,1 1099 1099 | gfin=i(gurMini5nLn 00-1_____ stran 19 I^BSSIfTMlGC A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation 4 REZULTATI DOLOČANJA ODZIVNE POVRŠINE ZA NAJVIŠJO TEMPERATURO SRAJČKE Preglednica 1 kaže izračunane vrednosti najvišje temperature srajčke s programom RELAP5/ MOD3.2 (PCT AP5) in napovedane vrednosti najvišje temperature srajčke z optimalnim statističnim cenilnikom (PCTOSC) v odvisnosti od petih vhodnih parametrov (SDC, TPDC, HTC, IDC, FPYF) 59 primerov. V našem primeru vhodni vektor G (n = 1, 2,...,59) sestavlja pet vhodnih parametrov, n izhodni vektor (H, n = 1, 2,...,59) pa ena komponenta, tj. najvišja temperatura srajčke PCTRELAP. Faktor širine za Gaussovo krivuljo f se izbere po kriteriju, da vpliv izračunanih točk razširimo do npr. sosednje točke in da je pri tem točnost prilagajanja še zadosti velika (R2 > 0,95). V našem primeru je za izbrani f = 0,25 imel srednji kvadratni pogrešek vrednost 1,55 K in koeficient določitve R2 vrednost 0,97. Z vstavitvijo vrednosti H, izračunanih s programom RELAP5/ MOD3.2, v enačbo (2) lahko napovemo najvišje temperature srajčke za poljubno kombinacijo vhodnih parametrov G znotraj podanih meja, če koeficiente C izračunamo po enačbi (3). Za nove izbrane vrednosti vhodnih parametrov G se določijo novi koeficienti C . Naslednje vprašanje, ki se nam zastavi je n. kakšne so vrednosti odzivne površine v točkah, za katere nimamo podanih izračunanih vrednosti? Ker je obravnavana odzivna površina petdimenzionalna, ni mogoča grafična predstavitev na eni sliki. Slike od 1 do 5 kažejo krivulje, pri katerih v osnovnem primeru spreminjamo samo en parameter. S krožcem so označene točke, izračunane s programom RELAP5/ MOD3.2, ki smo jih med seboj povezali s črtkano pikčasto linijo. Slike kažejo tudi vpliv faktorja za širino Gaussove krivulje na odzivno površino, dobljeno z OSC. Na slikah je prikazana tudi krivulja, dobljena z regresijsko analizo, ki je označena z “regr.”. Za določitev negotovosti najvišje temperature srajčke smo naključno spreminjali vhodne podatke in z OSC izračunali najvišje temperature srajčke. Ta postopek smo ponovili 100000 krat, in za rezultat dobili, da obstaja 95% verjetnost, da najvišja temperatura srajčke ne bo presegla 1157 K, srednja vrednost znaša 1085 K, medtem ko je negotovost razlika med 95. odstotkom in srednjo vrednostjo, in znaša 72 K. Srednja vrednost najvišje temperature srajčke z njej pripadajočo negotovostjo je precej pod dopustno mejo 1478 K. 5 RAZPRAVA Pravih vrednosti najvišjih temperatur srajčk med izračunanimi točkami ne poznamo. En možen približek je linearna odvisnost. OSC kot približek potegne krivuljo s prevojem, gledano v eni izmeri. Manjši ko je faktor f, bolj stopničast je prevoj. Želimo 4 RESULTS OF RESPONSE SURFACE DEVELOPMENT FOR PEAK CLADDING TEMPERATURE Table 1 shows the values of the calculated peak cladding temperatures (PCTRELAP5) and the pre-dicted peak cladding temperatures (PCTOSC) by OSC as a function of five input parameters for 59 cases. In our case, the input vector Gn (n = 1, 2,...,59) are five input parameters and the output vector Hn (n = 1, 2,...,59) has only one value, i.e. PCTRELAP5. The factor for the width of the Gaussian curve (f) is se-lected on the basis of criteria to stretch the influ-ence of a particular input data point into its sur-roundings, e.g. approximately to its neighbours and that the accuracy of fit is adequate (R2 > 0.95). In our case the root mean square error was equal to 1.55 K and the coefficient of determination was 0.97 for f = 0.25. By inserting the values for Hn, which were calculated with RELAP5/MOD3.2, into eq. (2) we can predict the peak cladding temperatures for any combination of input parameters G within parameter boundaries, if coefficients Cn are calcu-lated by eq. (3). The next question is, what are the values of the response surface in the points where the cal-culated values are not given? Because the con-structed response surface is five dimensional, visual presentation for all dimensions is not possible. Fig-ures 1 to 5 show the response surface with one parameter varied in base case calculation. The points calculated by the RELAP5/MOD3.2 code are marked with circles, connected with dash-dotted lines. In the figures is shown the influence of the factor for the width of the Gaussian curve on the response surface predicted by OSC. In the figures is also shown the regression curve which is labelled with “regr.”. To evaluate uncertainties in the peak clad-ding temperatures, input parameters were randomly selected by the Monte Carlo method and the peak cladding temperatures were estimated by OSC. This procedure was repeated 100000 times, and the result indicates that there is the 95% probability that the peak cladding temperature will not exceed 1157 K. The mean value is expected to be 1085 K, and the uncertainty of PCT is the difference between 95% and mean value, which is 72 K. The peak cladding temperature with its uncertainty is well above the cri-terion, 1478 K. 5 DISCUSSION The peak cladding temperatures between the calculated points are unknown values. One possible fit is a linear dependence. The OSC fit is an inflected curve in the one-dimensional case. The smaller is the factor f, the closer is the inflected curve to the grin^SfcflMISDSD ^BSfiTTMlliC | stran 20 A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation si sicer čimbolj blage prevoje, vendar je to odvisno od narave pojavov. Teže, ko jih je popisati, manjši faktor za širino Gaussove krivulje moramo uporabiti, da še dobimo zadovoljivo ujemanje v podanih izračunanih točkah. Posebej bi želeli opozoriti, da so približki na slikah 1 do 5 rezultat petdimenzionalne odzivne površine. Ker odzivno površino potrebujemo za integracijo Monte Carlo, za napoved verjetnosti zadostuje optimalni statistični cenilnik s stopničasto funkcijo. 1160 r--------------------------------- step function. We wanted smooth curves but smooth-ness is dependent on the nature of the phenomenon. The more complex is the phenomenon, the smaller is the factor for the width of Gaussian curve, which is needed to adequately fit the calculated points. It is worth noting that the fits in Figures 1 to 5 are the result of a five-dimensional response surface. Because the response surface is needed for the Monte Carlo integration, the robust response surface with step function transitions still satisfies for the probability evaluation. 1140 1120 1100 1080 1060 1040 1020 0.82 0.86 0.9 0.94 SPDC 0.98 1.02 Sl. 1. Najvišja temperatura srajčke v odvisnosti od enofaznega iztočnega koeficienta Fig. 1. Peak cladding temperature as a function of subcooled discharge coefficient 1200 1160 1120 1080 1040 1000 1 1.1 1.2 1.3 TPDC 1.4 1.5 1.6 Sl. 2. Najvišja temperatura srajčke v odvisnosti od dvofaznega iztočnega koeficienta Fig. 2. Peak cladding temperature as a function of two phase discharge coefficient 1120 1110 1100 1090 1080 1070 1060 0.75 0.85 0.95 1.05 HTC 1.15 1.25 Sl. 3. Najvišja temperatura srajčke v odvisnosti od toplotne prestopnosti Fig. 3. Peak cladding temperature as a function of heat transfer coefficient isfFIsJBJbJJIMlSlCšD I stran 21 glTMDDC A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation 1081 1079 1077 1075 1073 1071 1069 1067 1065 — -•¦ RELAP5 OSC, f=0.25 OSC, f=0.3 OSC, f=0.4 regr. 'v «. \\ " *\ \-^v 0.9 0.95 1 IDC 1.05 1.1 Sl. 4. Najvišja temperatura srajčke v odvisnosti od medfaznega trenja Fig. 4. Peak cladding temperature as a function of interphase drag coefficient 1080 1075 1070 1065 1060 1055 1050 0.9 0.95 / ^""T* r^5^ /y'^ *\ \ \^^ // Vs — -•¦ - RELAP5 'S? J? OSC, f=0.25 OSC, f=0.3 OSC, f=0.4 \ 's regr. ^—^» 1 FPYF 1.05 1.1 Sl. 5. Najvišja temperatura srajčke v odvisnosti od cepitvenega deleža zaostale toplote Fig. 5. Peak cladding temperature as a function of fission product yield factor V našem primeru optimalni statistični cenilnik skupaj s postopkom Monte Carlo uporabimo za določanje negotovosti termohidravličnih računalniških programov. Dobra lastnost OSC je, da je funkcija med dvema poznanima točkama monotona. To tudi pomeni, da bodo vse vrednosti ležale med najvišjo in najnižjo izračunano vrednostjo s programom RELAP5/MOD3.2. Druga dobra lastnost je, da se da določanje odzivne površine avtomatizirati. Za določitev odzivne površine ni predpisano število potrebnih točk. Več točk ko imamo, večje je zaupanje v rezultate. S primerno izbiro faktorja za širino Gaussove krivulje lahko dosežemo želeno točnost ujemanja napovedanih točk s podanimi izračunanimi točkami za še tako zapletene odvisnosti. Z zoževanjem širine Gaussove krivulje so prehodi med podanimi izračunanimi točkami vedno bolj stopničasti, zato optimalni statistični cenilnik deluje bolj grobo. V takem primeru je priporočljivo povečati število izračunanih točk. 6 SKLEP Prirejeni optimalni statistični cenilnik za potrebe računanja negotovosti termohidravličnih računalniških In our case, the optimal statistical estimator with the Monte Carlo method is used for the uncer-tainty evaluation of thermal-hydraulic computer codes. One good characteristic of the OSC is that between the two code-calculated values, the function is monotonic. This also means that all the values pre-dicted by OSC will be between the minimum and maximum RELAP5/MOD3.2 code calculated value. A second valuable characteristic of the OSC is that the response surface generation can be automated. The number of calculated points is not prescribed, but the higher the number of calculated points, the higher the confidence level. With the proper selection of the width of the Gaussian curve, the desired accuracy of fit in the calculated values can be obtained for very com-plex phenomena. By decreasing the width of the Gaussian curve the transitions between the points are increasingly stepwise, as a result, the response sur-face performance is crude. In such cases, it is recom-mended that more calculated points are provided. 6 CONCLUSION The adapted optimal statistical estimator for the uncertainty evaluation of the thermal-hydraulic grin^SfcflMISDSD VH^tTPsDDIK stran 22 A. Pro{ek - B. Mavko: Dolo~anje odzivne povr{ine - Response Surface Generation programov je pokazal, da se da določiti odzivno površino za pojave in procese z zapleteno in nelinearno odvisnostjo. Še več, določanje odzivne površine z OSC se da z računalniškim programom avtomatizirati. Avtomatiziran postopek določanja odzivne površine zelo razširi področje uporabe metod za določanje negotovosti preračunov s termo-hidravličnimi programi od izlivnih nezgod na poljubne nezgode, kar do zdaj ni bilo mogoče. codes showed that the response surface can be de-veloped for complex and non-linear phenomena and processes. Furthermore, the response surface generation can be automated. The automated procedure for response sur-face generation extends the use of uncertainty evaluation methods for thermal-hydraulic codes from a loss-of-coolant accident for which the uncertainty was evaluated for any accident. 7 LITERATURA 7 REFERENCES [I] Boyack, B.E et al. (1990) Quantifying reactor safety margin parts 1 to 6. Nuclear engineering and design, Amsterdam, 119, 1-117. [2] Prošek, A., B. Mavko, A. Stritar (1992) Ocena velike izlivne nezgode v jedrski elektrarni z analizo negotovosti. Zbornik Kuhljevi dnevi ’92, Slovensko društvo za mehaniko, Ljubljana, 178-185. [3] Mavko, B., A. Stritar, A. Prošek (1993) Aplication of code scaling, applicability and uncertainty methodology to large break LOCA analysis of two-loop PWE. Nuclear engineering and design, Amsterdam, 143, 95-109. [4] Ortiz, M.G., L.S. Ghan (1992) Uncertainty analysis of minimum vessel liquid inventory during a small-break LOCA in a B&W plant - an application of the CSAU methodology using the RELAP5/MOD3 computer code. NUREG/CR-5818, EGG.2665, Idaho National Laboratory. [5] Haskin, E.F., Bevan, B.D., C. Ding (1996) Efficient uncertainty analyses using fast probability integration. Nuclear engineering and design, Amsterdam, 166, 225-248. [6] Rombouts, D., Denil, D., Simon, C, C Matthys (1998) Westinghouse advanced safety analysis technology for plant power upratings. Proc. of nuclear energy in central Europe ’98, Nuclear Society of Slovenia, Ljubljana, 439-444. [7] Mavko, B., A. Prošek (1997) Peak cladding temperature response surface generation based on simulations of a small-break loss-of-coolant accident scenario. Proc. of 4th regional meeting: Nuclear energy in central Europe, Nuclear Society of Slovenia, Ljubljana, 605-612. [8] Grabec, I., W. Sachse (1991) Automatic modeling of physical phenomena: Applicaton to ultrasonic data. J.Appl. Phys., 69(9), 6233-6244. [9] Prošek, A. (1998) Ocena negotovosti realističnih simulacij potekov nezgod v jedrskih elektrarnah. Doktorska disertacija, univerza v Ljubljani, Fakultete za matematiko in fiziko, Ljubljana. [10] Prošek, A., B. Mavko (1999) Evaluating Code Uncertainty - II: An optimal statistical estimator method to evaluate the uncertainties of calculated time trends. Nuclear Technology, La Grange Park, Vol. 126, 186-195. [II] Prošek, A., B. Mavko (1999) Evaluating Code Uncertainty - I: Using the CSAU method for uncertainty analysis of a two-loop PWR SBLOCA, Nuclear Technology, La Grange Park, Vol. 126, 170-185. Naslov avtorjev: dr. Andrej Prošek profdr. Borut Mavko Inštitut Jožef Stefan Jamova 39 1000 Ljubljana Authors’ Address: Dr. Andrej Prošek Prof.Dr. Borut Mavko Jožef Stefan Institute Jamova 39 1000 Ljubljana, Slovenia Prejeto: Received: 21.4.1999 Sprejeto: Accepted: 29.2.2000