Annals of Nuclear Energy 28 (2001) 1219±1235 www.elsevier.com/locate/anucene
Application of a reduced order model to BWR corewide stability analysis Miguel CecenÄas-FalcoÂn a,*, Robert M. Edwards b a Instituto de Investigaciones EleÂctricas, Av. Reforma 113 Col. Palmira, Temixco, Mor 62490, Mexico The Pennsylvania State University, Nuclear Engineering Department, 231 Sackett Building, University Park, PA 16802, USA
b
Received 27 September 2000; accepted 23 October 2000
Abstract The determination of system stability parameters from power readings is a problem usually solved by time series techniques such as autoregressive modeling. These techniques are capable of determining the system stability, but ignore the physics of the process and focus on the determination of a nth order linear model. A nonlinear reduced order system is used in conjunction with estimation techniques to present a dierent approach for stability determination. The simulation of the reduced order model shows the importance of the feedback reactivity imposed by the thermal-hydraulics; the dominant contribution to this feedback is provided by the void reactivity, being a function of power, burnup, power distribution, and in general of the operating conditions of the system. The feedback reactivity is estimated from power measurements and used in conjunction with a reduced order model to determine the system stability properties in terms of the decay ratio. # 2001 Elsevier Science Ltd. All rights reserved.
1. BWR reduced-order model A BWR reactor is a complex system that requires an elaborate computer model to reproduce the dynamics of the system with enough detail to perform accurate steady-state calculations or transient analysis. Because of the complex feedback mechanisms involved in coupled thermal-hydraulic and neutron kinetics, those models take considerable computer resources and are used mainly for o-line ana* Corresponding author. Tel.: +52-7-318-3811; fax: +52-7-318-9854. E-mail address:
[email protected] (M. CecenÄas-FalcoÂn). 0306-4549/01/$ - see front matter # 2001 Elsevier Science Ltd. All rights reserved. PII: S0306-4549(00)00118-3
1220
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
lysis. Some examples of current codes used for BWR analysis are TRAC-BF1, RELAP5 and RETRAN. Other types of codes, such as LAPUR, NUFREQ-NP and STAIF, work in the frequency domain to obtain characterizations of the dynamic response in the form of a set of numeric transfer functions. A Reduced order model for BWR dynamics is a simple set of equations that reproduce the basic BWR dynamics. It lacks the complexity of high order models used for BWR analysis, but it is faster and easy to implement and execute. A simple phenomenological model was developed in 1986 by March-Leuba (1986) to reproduce the BWR dynamics. The model includes a one-node fuel dynamics that provides a low-frequency zero, and a two-node approximation for the void reactivity feedback, which includes a pair of complex conjugate zeros. The forward loop contains a linearization of the point kinetic equations, as seen in Fig. 1. The equations for the model are: dnr dt
nr0 nr cr
dcr lnr dt
lcr
dT a1
nr dt
nr0
1
2
a2 T
3
d2 d dT a5 T a4 K a3 dt dt2 dt
4
ext DT
5
The ®rst two equations of the Reduced order model are a linearization of the point kinetics equations about a perturbed equilibrium point. The fuel and void equations contribute to the feedback path, and the feedback reactivity is composed of the Doppler and void reactivities, which are then added to any external reactivity induced by control rods.
Fig. 1. Block diagram for the reduced order model.
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1221
Applying Laplace transform to (1) and (2), eliminating cr, and considering / >>l, the transfer function for reactivity to power is nr
s nr0
s l
s s
s =
6
Applying the Laplace transform to (3) and (4), the feedback path is composed by T
s a1 nr
s s a2
7
s K
s a5 T
s s2 a3 s a4
8
The open-loop transfer function is ®nally given by the expression nr0 a1 D K Ka5
s l s2
a3 s
a4 D D G
sH
s s s
s a2
s2 a3 s a4
9
where G
s and H
s are the forward and feedback loops of the system, respectively. The two complex zeros in the quadratic term indicate a pair of complex conjugate poles in the closed-loop response, whose imaginary part is a function of the parameter K. This gives a ®rst idea of the importance of the void reactivity gain, K, in the dynamic response of the system. Fig. 1 shows that there are two feedback gains, corresponding to the Doppler and void reactivity transfer function. The Doppler contribution to the feedback reactivity is very small compared with that of the void reactivity. The total feedback reactivity gain is dominated by the gain K, and this parameter will be referred to as the feedback gain in the following sections. The LAPUR code as documented by Otaduy and March-Leuba (1990) was used to obtain the transfer function for the feedback path. LAPUR is a frequency domain code for BWR analysis, and it provides numeric transfer functions for the nuclear and thermal±hydraulic phenomena occurring in a BWR core. The transfer function for the feedback path can be ®t to match the number of poles and zeros of H
s and obtain the values for the coecients a1±a5 and K. The model parameters for a typical BWR are given in Table 1, where the system is located at the critically stable point for 44% power. The value of the gain K at this point will be denoted as K0. The Bode plot of the closed-loop system is shown in Fig. 2. A low frequency resonance at about 0.4 Hz can be observed, typical of BWRs. The resonance disappears as the feedback gain K is reduced, and becomes sharper as the gain increases.
1222
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
Table 1 Values for a typical BWR Parameter
Value
a1 a2 a3 a4 a5 K0
10.53 0.4023 4.3944 7.3518 37.0294 0.00137
Fig. 2. Bode plot of the closed loop transfer function for dierent values of K.
Fig. 3 shows the root locus for the closed loop poles as the feedback gain increases from 0.0K0 to 2.0K0, where K0 is the value for a critically stable system. The complex conjugate poles of the closed loop system reach the real axis at K=K0; for values of K greater than K0 the system goes to the unstable region. The poles lying in the real axis travel to the left as K increases. A large real negative ®fth pole is not shown in the plot. The power level also aects the system stability, the closed-loop poles travel from the stable to the unstable region as power is increased keeping the same transfer function for the thermal-hydraulics. For feedback gains greater than 10K0 the system
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1223
Fig. 3. Root Locus for the linear model.
is unstable at all power levels, and for gains less than 0.5K0 the system is stable at all power levels. The value of K=K0 determines the limit of linear stability for a given power. The linear model predicts diverging oscillations for any larger K, but the real nonlinear system will present nonlinear characteristics that may evolve into a limit cycle. The linear reduced order model presented by March-Leuba (1986) is now modi®ed to include nonlinear normalized point kinetics with one delayed group in the forward loop, and the order of the dierential equation for voids is increased by one. The set of equations describing the BWR dynamics become: dnr
nr cr dt
10
dcr
nr dt
11
dT a1
nr dt
cr l nr0
a2 T
d3 d2 d dT a a a a K T 3 4 5 6 dt dt3 dt2 dt
12
13
1224
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
ext DT
14
The values for coecients a1±a6 are computed from LAPUR transfer functions. A curve-®t from Schoukens and Pintelon (1991) to the LAPUR feedback reactivity is performed to determine the values of the coecients a1±a6. fb
s Da1 s3 Da1 a3 s2
Ka1 Da1 a4 s Ka1 a6 Da1 a5 4 nr
s s
a2 a3 s3
a2 a3 a4 s2
a2 a4 a5 s a2 a5
15
Neglecting the Doppler contribution, the block-diagram algebra for Fig. 1 leads to fb
s Ka1 s Ka1 a6 4 3 nr
s s
a2 a3 s
a2 a3 a4 s2
a2 a4 a5 s a2 a5
16
where the constant K represents the gain of the feedback loop and also the void reactivity coecient, since it multiplies the void fraction change from steady state. This one-zero and three-pole transfer function simpli®es the implementation of the model, but can not reproduce with ®delity the dynamic characteristics of the system with an accurate ®t to the Bode diagram. The transfer function (15) suggests a 3zero and 4-pole ®t to the corresponding LAPUR transfer function. Fig. 4 shows this ®t, where the crosses represent the LAPUR values and the continuous line the ®tted transfer function. The power to total reactivity feedback transfer function for 28% core ¯ow and 44% power, obtained from the curve ®t in Fig. 4 is fb
s 0:0964s3 1:3590s2 2:1483s 14:3782 4 nr
s s 3:1247s3 8:6336s2 14:6429s 3:29221
17
The kinetics constants used to complete the model are given in Table 2. The parameter K behaves as a bifurcation parameter, and plays an important role in the system stability. For values of K less than a critical value K0, the system presents damped oscillations that lead to asymptotic stability. For values equal to K0, the system reaches a limit cycle, whose amplitude increases as K increases above K0. Table 3 shows how the decay ratio (DR) changes as the bifurcation parameter K varies from 0.6K0 to 1.2K0. Stability of boiling water reactors is closely related to the core ¯ow, a process variable not explicitly represented in the Reduced order model. The reduced order model can be modi®ed to re¯ect the impact of this important variable. Making a series of LAPUR analysis for the dierent powers and core ¯ows, a set of coecients for the reduced order model can be obtained and studied. Fig. 5 shows the trend for those coecients as a function of power for dierent core ¯ows. The relatively smooth plots suggest a well-behaved function can describe the coecient values. Based on this observation, a linear interpolation for all the six coecients for
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1225
Fig. 4. Curve ®t for the power to total feedback reactivity transfer function.
Table 2 Parameter values for a typical BWR at 28% core ¯ow and 44% power Parameter
Value 2.6410 3 0.0054 0.08 4.010 5 0.02113
D l K0
Table 3 Values for the decay ratio as a function of the feedback gain K
DR
RF
0.6K0 0.8K0 1.0K0 1.2K0
0.71 0.87 1.00 1.12
0.35 0.38 0.41 0.43
1226
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
Fig. 5. Coecients of Eqs. (12) and (13) as a function of core ¯ow. The continuous, dashed, and dotted lines represent 30, 35, and 40% of core ¯ow, respectively.
a given point in the power-¯ow map was developed. As a result, a reduced order model with coecients a1±a6 adapted to the current ¯ow will improve the accuracy for the predicted void reactivity. Given power and ¯ow measurements, the model can be used in a nonlinear estimation algorithm to determine a value for the bifurcation parameter K, which is related to the stability of the system. This value of K will re¯ect the plant operating conditions, and can be used to evaluate the stability of the system. 2. Application of the Reduced order model to stability estimation The reduced order model has been used previously to evaluate stability. MunÄoz et al. (1992) used the linear version to compute Lyapunov Exponents for BWRs. Turso et al. (1994) also used the linear model in a bank of Kalman ®lters, where each ®lter has a dierent plant operating point and the ®lter that best describes the current plant conditions is chosen by a maximum a posteriori probability calculation. Autoregressive modeling represents a well established method to determine system stability. A disadvantage of the AR method is the need of a relatively large amount of data to calculate the autoregressive parameters, and the problem it poses to update the decay ratio as new measurements are obtained. The adapted reduced
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1227
order model in conjunction with a nonlinear estimation algorithm provide the basis to continuously evaluate system stability. The approach of the method presented here is based on the previously studied relationship between the system bifurcation parameter, K, and the system stability. If the bifurcation parameter is estimated using a nonlinear model, the decay ratio can be assessed from this estimated value. The nonlinear estimation method chosen was the extended Kalman ®lter, and its performance is evaluated using computer simulations and an stochastic boiling channel model (CecenÄas-FalcoÂn and Edwards, 2000). In this application, the method is focused on corewide oscillations and uses power and core ¯ow measurements. 2.1. Introduction to the Kalman ®lter The reduced order model can be used in conjunction with power measurements to estimate the bifurcation parameter value. A nonlinear estimation problem is de®ned by Eqs. (10)±(14) if noisy power measurements are taken from the plant, all coecients a1±a6 are known and K is a parameter to estimate. A Kalman ®lter can be used to solve this nonlinear estimation problem, the extended version of the ®lter is chosen because it linearizes along the state trajectories during the estimation process and obtains the best estimate for the bifurcation parameter K. This procedure yields better results than the standard linear Kalman ®lter applied to a system linearized only once about the initial equilibrium point. A general description for a nonlinear system and a measurement equation is: : x f
x; u; t w
t
18
z h
x; u; t v
t
19
where x f() w(t) z h() v(t)
= = = = = =
State vector Vector of nonlinear functions plant process noise observable or measured plant output Nonlinear measurement equation measurement noise.
The system function f and the output function h are nonlinear functions of the state vector x; but, they can be linearized about the current operating point by expanding f and h in Taylor series and neglecting high order terms. The linearized system is : x Fx w
t
20
z
21
h
xn ; t Hx v
t
1228
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
where the matrices F and H are the Jacobians of f(x,u,t) and h(x,u,t), respectively. The equations for the extended Kalman ®lter are summarized in Table 4 (Brown and Hwang, 1992). 2.2. Stability estimation using the Kalman ®lter To apply the ®lter to the BWR equations, the reduced model is further reduced by the prompt-jump approximation to eliminate a high frequency pole and is augmented with a state equation for the bifurcation parameter, which is assumed to be constant, i.e. with zero dynamics. The third-order dierential equation for the voids is rewritten as a set of three ®rst-order dierential equations de®ning the auxiliary variables x2=da/dt and x3=dx2/dt. The nonlinear equations for the Kalman ®lter become: f1
dcr dt
l cr
ext DT
f2
dT a1
nr dt
f3
d x2 dt
24
f4
dx2 x3 dt
25
f5
dx3 dt
f6
dK 0 dt
a5
nr0
a4 x2
lcr
22
a2 T
a3 x3 Ka1
nr
23
nr0 K
a6
a2 T
26
27
Table 4 Extended Kalman ®lter equations System model
: x f
x; u; t w
t
Measurement model
z h
x; t v
t
Initial conditions
x0 ; P0
Kalman gain
Kk Pk HTk
Hk Pk HTk Rk
Updated estimate
x^ k x^ k Kk
zk
Error covariance
P k
I
Project estimate
x^ k1 x^ k f
x^ k ; u; tk1
tk1
Project error covariance
Pk1
1
z^k
Kk Hk Pk k Pk Tk
Qk
tk
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1229
Even though there is not an explicit dependence with core ¯ow, the set of coecients a1±a6 are ¯ow-dependent and they are computed by the ®lter. The Kalman estimator receives power, ¯ow, and the external reactivity as input from the plant or computer simulation and produces estimates for the state vector, including K. The Kalman estimator was implemented as a Matlab script, and has the inconvenience of relatively low computational speed. However, it is possible to implement the estimator in a high level programming language, such as C or FORTRAN, and obtain improved performance. Accurate estimates were obtained for all the system states, and the bifurcation parameter converged successfully to the values used in the simulations if the sampling rate is fast enough. A period of 0.05 s was set as an adequate sampling period to achieve convergence, and an arbitrary initial value of zero was set for the bifurcation parameter. An implementation of the Reduced-order BWR model with the estimated bifurcation parameter can be excited with an external reactivity perturbation and its timedomain response analyzed. The decay ratio is calculated from the step response. As an alternative to evaluate the decay ratio, the Reduced-order model can also be linearized around the initial equilibrium point to determine the complex pole values and the decay ratio is evaluated with DR e2=!d
28
the resonant frequency in hertz is obtained from the damped frequency: RF
!d 2
29
where is the real part of the pole, and !d is the absolute value of the imaginary part. The state-space representation of the Reduced-order model used in the Kalman ®lter estimator is: 3 2 lDnr0 lnr0 0 0 2 : 3 6 0 7 cr 7 6 7 6 : a Dn a n 6 7 6 1 r0 1 r0 6 T 7 6 a1 a2 0 0 7 7 6 : 7 6 7 6 7 6 7 6 7 6 0 0 0 1 0 7 6 : 7 6 7 4 x2 5 6 0 0 0 1 7 : 7 6 0 x3 5 4 Ka1 Dnr0 Ka1 nr0 Ka1 K
a6 a2 a5 a3 a4 2 3 cr 6 T 7 6 7 6 7 6 7
30 6 7 6 7 4 x2 5 x3
1230
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
2
y 1
Dnr0
nr0
0
0
3 cr 6 T 7 6 7 7 0 6 6 7 4 x2 5 x3
31
where cr, T, x1, x2 and a are small deviations from the equilibrium point produced by a perturbation on the external reactivity ext. The decay ratio is computed with the complex conjugate eigenvalues of this system and with the estimated K. A schematic diagram of the stability estimator is shown in Fig. 6, where the inputs are power, ¯ow, and APRM readings. The power input represents the reactor power, is composed by the plant process computer and is based on all the conditioned APRM signals. As an initial evaluation of the Kalman estimator, a set of points on the power-¯ow map is selected as test points. Three points at 30% core ¯ow are chosen to represent the natural circulation line, and three points at 40% core ¯ow are chosen to represent the line at low recirculation speed. These six points lie close or inside the region of possible unstable operation. LAPUR runs showed that all of the six points are stable, even when point 3 is very close to the stability boundary. Table 5 shows the six test points.
Fig. 6. Schematic of the stability estimator.
Table 5 Study cases and stability parameters based on LAPUR Case
Power (%)
Flow (%)
DR
RF (Hz)
1 2 3 4 5 6
35 40 45 35 40 45
30 30 30 40 40 40
0.74 0.86 0.98 0.47 0.57 0.66
0.35 0.38 0.41 0.36 0.39 0.42
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1231
A reduced order model with nonlinear point kinetics was prepared for each of the six cases, and implemented with the feedback gain adjusted to match the LAPUR decay ratio. The dynamic response to a reactivity step is used as input data to the Kalman estimator. The Kalman ®lter estimates the feedback gain from the time response of the reduced order model, and two decay ratios are then computed. The ®rst DR is obtained performing a simulation with the estimated feedback gain, and the second decay ratio is obtained from the complex eigenvalues of the linearized system evaluated with the estimated feedback gain. The Kalman estimator takes about 1 min of data under these conditions to converge to the feedback gain value, as seen in Fig. 7. The results for the six cases are summarized in Table 6. The column labeled K contains the values of the feedback gain used in the generation of the power measurements, columns DR and RF contain the values for the decay ratio and resonant frequency obtained from the time response to a step perturbation. The Kalman estimator does an excellent job estimating the value of the feedback gain, even when the power measurements are generated with the one delayed group point kinetics and the Filter model uses the prompt-jump approximation. The estimated decay ratios show a deviation from the real value, even when the estimated feedback gain converges to the correct value. This deviation is produced by the dierent approach taken to compute the decay ratios. The system, or real, stability parameters are obtained from a step response of the nonlinear model, averaging the decay ratios
Fig. 7. Estimated bifurcation parameter K for case 1.
1232
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
obtained from successive peaks and ignoring the ®rst peak; the estimated values represent the asymptotic decay ratio contributed by the dominant complex conjugate poles of the linear model, obtained evaluating Eqs. (28) and (29). Fig. 8 shows a comparison of the estimated and real stability parameters. The estimator developed in this section uses a dynamic response introduced by a reactivity step. For real applications, the information of the system dynamics is contained in the power noise. The random ¯uctuations inherent in the boiling process and in the core ¯ow will produce power ¯uctuations due to variations in the moderation. These noisy ¯uctuations are enough for the ®lter to converge to an estimated feedback gain. 3. Stability estimator testing The preliminary tests to the stability estimator were based on computer-generated power measurements obtained with reduced-order models and a set of LAPUR transfer functions. A completely independent set of simulated measurements is prepared with a nuclear-coupled boiling channel model (CecenÄas-FalcoÂn and Edwards, 2000). This model solves the conservation equations for a boiling channel coupled to point or modal kinetics. Three regions are considered, one phase, subcooled and Table 6 Estimated decay ratios and resonant frequencies for the six study cases Estimated Case
K
DR
RF
K
DR
RF
1 2 3 4 5 6
0.0410 0.0360 0.0324 0.0310 0.0250 0.0223
0.7376 0.8619 0.9829 0.4687 0.5664 0.6606
0.3518 0.3836 0.4137 0.3222 0.3540 0.3871
0.0410 0.0360 0.0325 0.0310 0.0250 0.0223
0.7105 0.8290 0.9440 0.4672 0.5507 0.6376
0.3512 0.3833 0.4137 0.3220 0.3532 0.3868
Table 7 Estimated decay ratios for the six study cases using 400 seconds of noisy data generated with the stochastic nuclear-coupled boiling channel Case
1 2 3 4 5 6
Autoregressive results
Kalman estimator
DR
RF
DR
RF
0.749 0.931 0.959 0.557 0.641 0.718
0.382 0.421 0.445 0.372 0.432 0.444
0.777 0.927 0.995 0.481 0.700 0.787
0.375 0.409 0.428 0.328 0.393 0.420
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1233
Fig. 8. Comparison of the system (a) decay ratios and (b) resonant frequencies with the estimated values.
1234
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
Fig. 9. Comparison of the autoregressive and the estimated stability parameters. (a) decay ratios, (b) resonant frequencies.
M. CecenÄas-FalcoÂn, R.M. Edwards / Annals of Nuclear Energy 28 (2001) 1219±1235
1235
bulk boiling. The system is excited with white noise and the resulting power time series is used as input for the Kalman estimator. The boiling channel model was set up for the same six power/¯ow conditions de®ned in Table 5. Table 7 contains the decay ratios obtained with the boiling channel, where the void reactivity factor was de®ned as a second order equation in terms of the power to ¯ow ratio. Because of the stochastic nature of the boiling channel, decay ratios computed under the same operating conditions may not yield identical results. For each point in Table 7, a total of 10 noisy sequences of 400 s sampled at 20 Hz were generated with the nuclear-coupled boiling channel to study the statistical ¯uctuations in the decay ratios. The series were analyzed using both the autoregressive model and the Kalman-based stability estimator. Fig. 9 shows a comparison of the average decay ratios computed with both methods. The Kalman estimator provides slightly higher decay ratios than the AR method. 4. Conclusions The reduced order model is used to estimate the feedback reactivity from the thermal±hydraulics with an extended Kalman ®lter. The Kalman estimator achieved fast convergence for computer simulated data on the same model of the estimator. The capability of the estimator was tested with a set of data computed with a completely independent model. Results indicate that the estimator has a tendency to slightly overpredict the decay ratio with respect to autoregressive models, and also to underpredict the resonant frequency. As advantages for the estimator, it is recursive and can continuously provide decay ratio and resonant frequency estimations without signi®cant memory requirements. The combination of the reduced order model, estimation algorithms and data acquisition into a decay ratio estimator presents an alternative for stability monitoring in BWR plants. References Brown, R.G., Hwang, P.Y.C., 1992. Introduction to Random Signals and Applied Kalman Filtering, 2nd Edition. John Wiley and Sons. CecenÄas-FalcoÂn, M., Edwards, R.M., 2000. Stability monitoring tests using a nuclear-coupled boiling channel. Nuclear Technology 131(1). March-Leuba, J., 1986. A reduced-order model of boiling water reactor linear dynamics. Nuclear Technology 75. MunÄoz-Cobo, J.L. et al., 1992. Dynamic reconstruction and Lyapunov exponents from time series data in boiling water reactors. Application to BWR stability analysis. Annals of Nuclear Energy 19(4). Otaduy, P.J., March-Leuba, J., 1990. LAPUR User's Guide. NUREG/CR-5421 ORNL/TM-11285. Schoukens, J., Pintelon, R., 1991. Identi®cation of Linear Systems: A Practical Guide to Accurate Modeling. Pergamon Press. Turso, J.A. et al., 1994. Boiling water reactor stability analysis via Kalman ®lter-based state estimation and maximum a posteriori detection. Trans. American Nuclear Society, Vol. 71, November.