# A Novel Nonlinear BWR Stability Indicator Based on the Sample Entropy

Authors: Omar Alejandro Olvera-Guerrero; Alfonso Prieto-Guerrero; Gilberto Espinosa-Paredes
Journal: Science and Technology of Nuclear Installations
 (2018)

---

## Abstract

BWRs are thus far the simplest energy systems to transform fission energy into electrical power. However, there are still many aspects in their operation that, under certain conditions, may induce BWR unstable behavior. The default indicator to study BWR unstable behavior is the Decay Ratio (DR). However, due to the fact that BWRs show very complex responses under instability and responses that may even be chaotic, the DR might not be a suitable choice to rely on to accommodate for such intricate behavior. In this work a novel methodology based on the Sample entropy (SampEn) and the noise-assisted multivariate empirical mode decomposition (NA-MEMD) is introduced. Such methodology was developed thinking for a real time-implementation of a stability monitor. The proposed methodology was tested with a set of signals that stem from several nuclear power plants in operation today that have experienced in the past unstable events, each one of a different nature.


---

## Body

## 1. Introduction

The issue of BWR instability has been a topic of concern of many scientific and technological works for more than forty years. BWR unstable events are quite rare and may happen during BWR start up or during transients that may interfere with the conventional operating region of the reactor. BWR instability phenomena are strongly linked with the interaction between thermal hydraulic and neutron kinetics processes such interaction under certain conditions and may induce undesirable BWR behavior that could jeopardize BWR operation. Thus, a reliable prediction for the onset of BWR instability is of utter importance for BWR core safety [1].D’Auria [2] presented all BWR unstable events recorded up to 1999, which were due to various incidents (e.g., [3, 4]), while others were intentionally induced for experimental purposes (e.g., [5, 6]). Unstable oscillations in the neutron flux were recorded during these unstable events through the use of the electronic instrumentation available in the reactor. After the first events of instability were observed, the respective authorities fostered the growth of research projects to study the physical phenomena involved in BWR instability.The most common instability type for commercial BWRs is the density-wave instability [7–9]. This type of instability can be described as follows: given a flow perturbation, awaveof voids travels upwards through the channel producing a pressure drop (Δp=Δp1ϕ+Δp2ϕ, is approximately constant, where Δp1ϕ is the pressure drop in the single phase in the liquid region and Δp2ϕ is the pressure drop in the two-phase flow region) which is delayed with respect to the original perturbation. An increase in flow might induce an increase in pressure drop and a negative feedback that reduces the flow perturbation. The density wave phenomena delays this feedback and, at some frequency, the delay is equivalent to a 180∘ phase lag; so at this frequency the pressure drop feedback is positive and if the gain is large enough, the channel flow becomes unstable and oscillates at that frequency highly concentrated around 0.5 Hz, where the associated resonance frequency is related with the fluid transit time along the heated channel, which could differ from 2 seconds, i.e., approximately 1 Hz.In practice, the most commonly used indicator to study BWR stability due to wave oscillation is the decay ratio (DR), which is computed from the impulse response function obtained from a model fitting a linear second order system, accommodating the core behavior of a BWR. The DR is the current output indicator of most, if not all, of the stability monitors implemented in NPPs today [10]. The DR as a feasible BWR stability measure has been widely accepted but it has been observed that a BWR working at an operating point with a small DR can be close to instability [11]. Sometimes, the DR jumps from the stable to the far unstable region [12]. According to these works, the DR might not be a reliable option, under certain conditions. Besides, the need for linear and stationary signals might be a handicap for DR estimation. Therefore, it is necessary to explore new methodologies and indicators adapted to accommodate at their theoretical framework nonlinear and nonstationary behavior in order to study BWR unstable phenomena with as much realism as possible.We are trying to convey in this work this idea, for this we explore the Sample Entropy (SampEn) linked to the noise-assisted multivariate empirical mode decomposition (NA-MEMD) to infer whether BWR signals are stable or not. The SampEn [13] is a measure that provides an index of signalcomplexity orirregularity of a time series. In this way, the SampEn might work as a possible nonlinear BWR stability indicator. SampEn was in principle developed almost exclusively to analyze physiological time series [13] but its utility has expanded to other domains. For instance, SampEn has been tested with daily weather temperature to measure climate complexity [14]. It has been used to measure the complexity of the dynamic reconfiguration of the brain [15] to infer its association with normal aging.To properly estimate the SampEn from real BWR signals, the NA-MEMD [16] was explored. The NA-MEMD is an algorithm that decomposes non-stationary signals that stem from nonlinear sources. This technique also mitigates themode mixingphenomena [17] of the default EMD [18] technique. The EMD is the root technique that inspired the development of the NA-MEMD. The NA-MEMD produces a local and fully data-driven decomposition of a studied signal into its fast and slow oscillations. At the end of the procedure, the studied signal can be expressed as a sum of amplitude and frequency modulated (AM-FM) functions called intrinsic mode functions (IMFs), or simply calledmodes, plus a residue (the trend) of the decomposition. The combination of the NA-MEMD and the Hilbert transform is known as the Hilbert-Huang transform (HHT). The methodology we introduce here is based on the HHT and it computes an indicator linked to BWR stability, in this case the aforementioned SampEn. The NA-MEMD decomposes the analyzed BWR signal into IMFs. One or more of these extracted modes can be associated to the instability problem in BWRs. Through HHT it is possible to get the instantaneous frequency (IF) linked to each IMF. By tracking this IF and the SampEn of the IMF linked to instability, the estimation of the SampEn-based stability indicator is accomplished. The methodology proposed in here is a continuation of two previous works [19, 20], in which a Shannon Entropy estimator [21] was used in conjunction with other members of the EMD family to study the stability of artificial and real BWR signals. In the present work the Shannon Entropy estimator was replaced by the SampEn practical formula to infer whether BWR signals are stable or not.This work is organized as follows: in Section2, a full review of the SampEn and of the NA-MEMD techniques is given. In Section 3, the methodology to compute the instantaneous frequency and the proposed SampEn novel stability indicator is described. In Section 4, the validation of the methodology presented in this paper is performed doing experiments with real signals taken from the Forsmark [5] and Ringhals [6] stability benchmarks and from a typical BWR power plant during stable and unstable operation. Our major observations regarding our novel methodology based on SampEn and the NA-MEMD are thoroughly discussed in Section 5.


## 2. Preliminaries

### 2.1. The Sample Entropy

The procedure for computing SampEn was introduced by Richman and Moorman [13]. In this work, we grant a brief summary of their findings. SampEn is a measure of complexity. Let x=[x1,x2,...,xN] be a time series of length N. To compute the complexity of this time series via SampEn, follow the next steps:(1)
Build a vectorvi with m consecutive data points taken from x(1)vi=xi,xi+1,…,xi+m-1wherem is the length of sequences to be compared, also called theembedding dimension.(2)
For eachi define (2)Cim=1N-m-1∑j=1,j≠iN-mΘr-vi-vj1wherei varies in the interval (1≤i≤N-m). In here, r is the tolerance for accepting matches, r=s×std(x), where s is a scaling parameter and std(x) is the standard deviation of x. Θ· is the Heaviside function: (3)Θx=0,x<01,x≥0and·1 is the Chebyshev distance, defined as(4)vi-vj1=max⁡xi-xj,xi+1-xj+1,…,xi+m-1-xj+m-1(3)
Cim represents the proportion of vj(i≠j) whose distances to vi are less than r. Now, for each i we also define(5)Cim+1=1N-m-1∑j=1,j≠iN-mΘr-vi-vj1whereCim+1 represents the proportion corresponding to the dimension of m+1. Cim and Cim+1 have the same mold, but embedding vectors in both cases are defined in different spaces.(4)
Average across all embedding vectors, to obtain(6)Um=1N-m∑i=1N-mCimand(7)Um+1=1N-m∑i=1N-mCim+1(5)
The SampEn is computed as(8)SampEn=-ln⁡Um+1UmSampEn is the negative natural logarithm of the conditional probability that two sequences similar form points remain similar at the next point, where self matches are not included in calculating the probability. Thus, alower value of SampEn indicates more self-similarity (i.e., high order) of the studied time series whereas ahigher value of SampEn points to higher complexity of the time series.The calculation of SampEn requiresa priori determination of two unknown parameters, m and r (the length of data N is up to the user). The suggested values of r are located in the range of 0.1-0.2 times the standard deviation of the studied signal x. In this work, we looked for the r value that grants the global SampEn maximum because this maximum value leads to the correct interpretation of signals complexity [22]. So, this r of interest is fixed in our simulations at r=0.1. The value of m can be computed via the estimation of false nearest neighbor [23]. In our particular BWR stability discipline, m that was found through false nearest neighbors for our signals was most of the time equal to 2. Therefore, m is fixed at m=2 for all of our computations. SampEn is, theoretically speaking, a fraction in the interval 0<SampEn<∞. However, the next two formulas can be used to find the lower bound and the upper bound of SampEn for fixed values of m and N.The lower bound is computed as(9)2N-m-1N-m-1The upper bound is computed as(10)ln⁡N-m+ln⁡N-m-1-ln⁡2For all of our simulationsm,r and N are fixed at m=2, r=0.1×std(x) and N=300 data points. Thus the lower bound is practically ≈0 whereas the upper bound is close to 11. However, the upper bound value was never attained for any simulation.


### 2.2. The Noise-Assisted Multivariate Empirical Mode Decomposition

The multivariate empirical mode decomposition [16], commonly referred to as the MEMD, is a nonlinear filter to make the classic EMD [18] suitable for processing of multichannel signals. The behavior of the MEMD was analyzed in the presence of white Gaussian noise [20, 24] and it was found that the MEMD in essence acts as a dyadic filter bank on each channel of the multivariate input signal; such MEMD property is illustrated in Figure 1 and its steps are given below.Figure 1
Averaged spectra of IMFs (1-9) obtained from 50 realizations of eight-channel white Gaussian noise via MEMD. Reproduced from [20].The MEMD has the property of aligning the IMFs ormodesfrom different channels across the same frequency range and the NA-MEMD technique was developed to aid resolve themode mixingproblem observed in the current MEMD, two key issues for real applications. The NA-MEMD technique, which makes use of the quasi-dyadic filter bank properties of MEMD on white noise (Figure 1), is capable of reducing themode mixingproblem for classes of signals where the quasidyadic filter bank structure proves useful (such property applies for BWR signals). The NA-MEMD adds q extra channels of white noise to the input studied signal before being decomposed by the MEMD. The addition of these q extra channels helps to establish a uniformly distributed reference scale which, in turn, results in correspondingmodesexhibiting a quasi-dyadic filter bank structure. The extracted IMFs corresponding to the q extra channels of white noise are then discarded yielding a set of IMFs linked with only the original input signal of interest. The NA-MEMD procedure is summarized in the following next steps [20]:Step 1.
Create an uncorrelated Gaussian white noise time-series (q-channels) of the same length as that of the input.Step 2.
Add the noise channels (q-channels) created in Step 1 to the input multivariate (L-channels) signal, obtaining an (L+q)-channel signal.Step 3.
Process the resulting (L+q)-channel multivariate signal using the MEMD [16] algorithm, to obtain multivariate IMFs.Step 4.
From the resulting (L+q)-variate IMFs, discard theq channels corresponding to the noise, giving a set ofL-channel IMFs corresponding to the original signal.
However, it should be mentioned that the NA-MEMD for mitigating themode mixingbecomes more useful for signals in which the dyadic filter bank decomposition is relevant. This is the case for the studied BWR signals.


## 2.1. The Sample Entropy

The procedure for computing SampEn was introduced by Richman and Moorman [13]. In this work, we grant a brief summary of their findings. SampEn is a measure of complexity. Let x=[x1,x2,...,xN] be a time series of length N. To compute the complexity of this time series via SampEn, follow the next steps:(1)
Build a vectorvi with m consecutive data points taken from x(1)vi=xi,xi+1,…,xi+m-1wherem is the length of sequences to be compared, also called theembedding dimension.(2)
For eachi define (2)Cim=1N-m-1∑j=1,j≠iN-mΘr-vi-vj1wherei varies in the interval (1≤i≤N-m). In here, r is the tolerance for accepting matches, r=s×std(x), where s is a scaling parameter and std(x) is the standard deviation of x. Θ· is the Heaviside function: (3)Θx=0,x<01,x≥0and·1 is the Chebyshev distance, defined as(4)vi-vj1=max⁡xi-xj,xi+1-xj+1,…,xi+m-1-xj+m-1(3)
Cim represents the proportion of vj(i≠j) whose distances to vi are less than r. Now, for each i we also define(5)Cim+1=1N-m-1∑j=1,j≠iN-mΘr-vi-vj1whereCim+1 represents the proportion corresponding to the dimension of m+1. Cim and Cim+1 have the same mold, but embedding vectors in both cases are defined in different spaces.(4)
Average across all embedding vectors, to obtain(6)Um=1N-m∑i=1N-mCimand(7)Um+1=1N-m∑i=1N-mCim+1(5)
The SampEn is computed as(8)SampEn=-ln⁡Um+1UmSampEn is the negative natural logarithm of the conditional probability that two sequences similar form points remain similar at the next point, where self matches are not included in calculating the probability. Thus, alower value of SampEn indicates more self-similarity (i.e., high order) of the studied time series whereas ahigher value of SampEn points to higher complexity of the time series.The calculation of SampEn requiresa priori determination of two unknown parameters, m and r (the length of data N is up to the user). The suggested values of r are located in the range of 0.1-0.2 times the standard deviation of the studied signal x. In this work, we looked for the r value that grants the global SampEn maximum because this maximum value leads to the correct interpretation of signals complexity [22]. So, this r of interest is fixed in our simulations at r=0.1. The value of m can be computed via the estimation of false nearest neighbor [23]. In our particular BWR stability discipline, m that was found through false nearest neighbors for our signals was most of the time equal to 2. Therefore, m is fixed at m=2 for all of our computations. SampEn is, theoretically speaking, a fraction in the interval 0<SampEn<∞. However, the next two formulas can be used to find the lower bound and the upper bound of SampEn for fixed values of m and N.The lower bound is computed as(9)2N-m-1N-m-1The upper bound is computed as(10)ln⁡N-m+ln⁡N-m-1-ln⁡2For all of our simulationsm,r and N are fixed at m=2, r=0.1×std(x) and N=300 data points. Thus the lower bound is practically ≈0 whereas the upper bound is close to 11. However, the upper bound value was never attained for any simulation.


## 2.2. The Noise-Assisted Multivariate Empirical Mode Decomposition

The multivariate empirical mode decomposition [16], commonly referred to as the MEMD, is a nonlinear filter to make the classic EMD [18] suitable for processing of multichannel signals. The behavior of the MEMD was analyzed in the presence of white Gaussian noise [20, 24] and it was found that the MEMD in essence acts as a dyadic filter bank on each channel of the multivariate input signal; such MEMD property is illustrated in Figure 1 and its steps are given below.Figure 1
Averaged spectra of IMFs (1-9) obtained from 50 realizations of eight-channel white Gaussian noise via MEMD. Reproduced from [20].The MEMD has the property of aligning the IMFs ormodesfrom different channels across the same frequency range and the NA-MEMD technique was developed to aid resolve themode mixingproblem observed in the current MEMD, two key issues for real applications. The NA-MEMD technique, which makes use of the quasi-dyadic filter bank properties of MEMD on white noise (Figure 1), is capable of reducing themode mixingproblem for classes of signals where the quasidyadic filter bank structure proves useful (such property applies for BWR signals). The NA-MEMD adds q extra channels of white noise to the input studied signal before being decomposed by the MEMD. The addition of these q extra channels helps to establish a uniformly distributed reference scale which, in turn, results in correspondingmodesexhibiting a quasi-dyadic filter bank structure. The extracted IMFs corresponding to the q extra channels of white noise are then discarded yielding a set of IMFs linked with only the original input signal of interest. The NA-MEMD procedure is summarized in the following next steps [20]:Step 1.
Create an uncorrelated Gaussian white noise time-series (q-channels) of the same length as that of the input.Step 2.
Add the noise channels (q-channels) created in Step 1 to the input multivariate (L-channels) signal, obtaining an (L+q)-channel signal.Step 3.
Process the resulting (L+q)-channel multivariate signal using the MEMD [16] algorithm, to obtain multivariate IMFs.Step 4.
From the resulting (L+q)-variate IMFs, discard theq channels corresponding to the noise, giving a set ofL-channel IMFs corresponding to the original signal.
However, it should be mentioned that the NA-MEMD for mitigating themode mixingbecomes more useful for signals in which the dyadic filter bank decomposition is relevant. This is the case for the studied BWR signals.


## 3. BWR Stability Methodology Based on the SampEn and the NA-MEMD

The methodology to detect stable or unstable states (this methodology works with raw data. So, there is no need to preprocess the analyzed data) in a core of a BWR, based on the Sample entropy and the NA-MEMD, is presented in the next steps:(1)
The studied raw signal (APRM or LPRM) obtained from the BWR is segmented in windows of 60 s of duration.(2)
Each segment (of 60 s of time span) is considered as just an independent channel (L=1) added withq independent channels of white Gaussian noise and decomposed through the NA-MEMD (for all of our computer simulationsq=2).(3)
After decomposition, theq channels corresponding to the noise are discarded, giving a set of IMFs corresponding to each decomposed studied segment.(4)
The Hilbert transform of each IMF (i.e., the HHT) is computed to get the instantaneous frequency of each IMF in the studied segment.(5)
When tracking these frequencies, it is possible to get themode(ormodes) linked to the density wave oscillation. In this way, only themode (ormodes) associated to BWR instability is considered for further processing whereas the remaining extractedmodes are ignored.(6)
The SampEn of the trackedmodeof interest (IMF) is computed for N=300 (the length of the analyzed segment of 60 s), m=2 and r=0.1×std(segment) to attain the global maximum of the SampEn that leads to the correct interpretation of the studied segment complexity [22].(7)
The mean, median, and the variance of the SampEn values are computed and averaged along all the studied segments of 60 s of the original signal.AhighSampEn estimate (>1) points out high irregularity and high unpredictability of the studied BWR signal (LPRM or APRM) and thishighvalue might indicate that the studied signal could be merely broad band non-coherent noise (BBNCN). This BBNCN, from the BWR point of view is a sign of stable behavior. On the other hand, alowSampEn value (<0.8) points to high predictability and regularity of the studied BWR signal. Suchlowvalue might indicate that the studied signal could be a cyclic one. Such waveform is a sign of BWR unstable behavior, due to the presence of a density wave that manifests in the LPRM or APRM signals as a cyclic waveform with a natural frequency that oscillates very close to 0.5 Hz. For forced convection phenomena, in the circulation natural convection the effect of the chimney changes this value [8, 9].


## 4. Results and Discussions

Now we validate the methodology given in Section3 with the next data sets, wherefs is the sampling frequency of each set:(1)
Data set from a typical BWR:
(a)
208 steady state condition signals (i.e. stable signals),fs = 5 Hz.(b)
96 unstable signals,fs = 10 Hz.(2)
Case 4 of the Forsmark stability benchmark [5]:
(a)
1 APRM signal plus 22 LPRM signals,fs = 12.5 Hz.(3)
Case 9 cycle 14 of the Ringhals stability benchmark [6]:
(a)
72 LPRM signals (located within 2levelsof the Ringhals reactor, eachlevelhas a total of 36 LPRMs,suchlevelsare called 2 and 4),fs = 12.5 Hz.


### 4.1. Stable Signals

The aim of this test is to observe the demeanor of the BWR in steady-state conditions. The results obtained are presented in Figures2–8. Figure 2 shows the plot of one APRM stable signal of a typical BWR. Figure 3 shows in red dotted line, a studied segment of 60 s of duration of the APRM signal with its power spectral density (PSD) estimate (bottom of figure). Figure 4 shows the first four extracted IMFs of the studied signal, which happens to be noise that is decomposed in a quasidyadic fashion by the NA-MEMD. Figure 5 shows the power spectral density (PSDs) of some IMFs.Figure 2
Stable signal from a typical BWR.Figure 3
Analyzed segment of the studied APRM signal.Figure 4
A plot of the first 4 extracted IMFs in the analyzed segment.Figure 5
PSD estimate of the first 6 extracted IMFs in the analyzed segment.Figure 6
Tracked IMF (IMF 3) of interest and its associated instantaneous frequency (IF) around 0.5 Hz computed through the Hilbert-Huang transform.Figure 7
SampEn along time for the studied APRM stable signal.Figure 8
Computed mean(SampEn) of the studied 208 stable signals.Figure6 shows the tracked IMF (IMF 3) of interest whose instantaneous frequency (IF) oscillates around 0.5 Hz. The computed SampEn of this IMF is 0.9351. SuchhighSampEn value indicates that the studied signal is BBNCN, so the studied segment is stable and far from an unstable behavior. In this stable case, amodeexist whose IF is close to 0.5 Hz. However, the IMF 3 time series looks distant from a cyclic function we look for when an unstable event triggers in the core. All of the IFs were computed through the Hilbert-Huang transform of the extracted IMFs. Figure 7 shows the computed SampEn along time for all of the studied segments of 60 s; most of the SampEn estimates have a value slightly greater than 1 (i.e., the studied signal is stable). Table 1 shows the mean(SampEn), median(SampEn) and std(SampEn) of the time series shown in Figure 7.Table 1
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed APRM.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) APRM 1.1288 1.0994 0.1398Figure8 shows a plot of the computed mean SampEn estimates of the 208 signals (each SampEn estimate is a blue point). A pattern appears, all of the mean SampEn estimates for the studied signals are higher than 1. Thus, the mean SampEn values point to BBNCN. Through SampEn estimation, the proposed stability methodology is able to detect stable signals and classify them in the stable category if their SampEn is greater than 1.


### 4.2. Unstable Signals

The studied signal in this subsection stems from an instability event. This instability case is registered by 96 LPRM signals. For reasons of space, only the analysis of one signal is presented. Figure9 shows the plot of the studied signal (from LPRM 2); this signal looks like a cyclic function. Figure 10 shows a studied segment of 60 s (red dotted line) and its PSD estimate (bottom of figure), the segment looks like a noisy cyclic function. Figure 11 shows the first 4 extracted IMFs of the studied segment; the first 3 IMFs are linked with acquisition noise, whereas IMF 4 looks like a cyclic function (this is in fact the type of waveform commonly associated with density wave unstable events). Figure 12 shows the PSD estimates of the first 6 IMFs, the energetic content of IMF 4 is highly concentrated around 0.5 Hz, and the energetic content of the other IMFs is meaningless next to the energy of IMF 4. Figure 13 shows the time domain series of IMF 4 (which happens to be themodeassociated to BWR instability) and its associated IF that oscillates around 0.5 Hz. The IF was computed via HHT. Now, this IF (IF 4) carries physical meaning [7]. The SampEn estimate of IMF 4 (IMF ormodelinked to instability, see IF in Figure 13 and PSD estimate in Figure 12) is 0.6181; this estimate is a clear indication of signal regularity and low complexity. The SampEn points to a simple function; in this case a cyclic function that in our context is linked with a density wave type of instability. Thus, the studied segment is unstable. Figure 14 shows the SampEn estimates along time for all of the studied segments of the signal of interest; all estimates are smaller than 1.Figure 9
Unstable signal from a typical BWR.Figure 10
Analyzed segment of the studied LPRM 2 signal.Figure 11
A plot of the first extracted 4 IMFs in the analyzed segment.Figure 12
PSD of the first 6 extracted IMFs in the analyzed segment. The PSD of the IMF 4 is highly concentrated around 0.5 Hz.Figure 13
Tracked IMF (IMF 4) of interest linked to BWR instability and its associated instantaneous frequency (IF) around 0.5 Hz computed through the HHT.Figure 14
SampEn estimates of the IMFs ormodes of interest along time.Table2 shows the mean(SampEn), median(SampEn), and std(SampEn) of the time series shown in Figure 14. Now, the tracked IF (IF 4) along time carries physical meaning. Thus, the mean(IF) and std(IF) values of all of the studied segments are now provided in this table.Table 2
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed LPRM 2. Now the mean(IF) and std(IF) are included.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) Mean(IF) [Hz] Std(IF) [Hz] LPRM 2 0.5425 0.5258 0.0544 0.5374 0.0030Figure15 shows a plot of the computed mean SampEn estimates of the 96 unstable signals (each SampEn estimate is a red point). A pattern appears, all of the mean SampEn estimates for the studied signals are smaller than 0.8. Thus, the mean SampEn values point to visually regular signals such as cyclic functions. Through SampEn estimation, we are able to build a binary classifier that separates BBNCN signals (linked to stability) from cyclic ones (linked to density wave instability). One core feature of SampEn is that it accommodates for complex nonlinear and nonstationary data whereas conventional Decay Ratio estimates must assume beforehand that the studied BWR signals behave linearly (an assumption that is false in real life).Figure 15
Computed mean(SampEn) of the studied 96 unstable signals.


### 4.3. Comparison between Stable and Unstable States

To complete the analysis of the stable/unstable signals through SampEn, Figure16 shows a comparison of the mean SampEn values of the studied sets (stable set of 208 signals plus a set of 96 unstable signals). The computed SampEn mean values of the stable signals fluctuate very close to 1.2 whereas the SampEn mean values of the unstable signals oscillate very close to 0.6. So the SampEn detects with success stable BWR signals from unstable ones via complexity analysis (of course with the aid of the NA-MEMD to denoise these signals and isolate themodelinked to instability). There is a distance of ≈0.6 between the two states (stable from unstable). So, by fixing a threshold value around 0.9, it is possible to differentiate one state from the other properly. Any SampEn estimate of a BWR signal segment higher than 0.9 is stable whereas any lower SampEn value from this threshold is unstable.Figure 16
Comparison of the mean(SampEn) values for both studied sets (a stable set of 208 APRM signals versus an unstable set of 96 LPRM signals).


### 4.4. Case 4 of the Forsmark Stability Benchmark

Case 4 of the Forsmark stability benchmark [5]. This event is considered a difficult case to study by the complexity of the phenomena. For reasons of space, only the analysis of one case will be presented in detail. This case contains a mixture between a global oscillation and a regional oscillation. This particular event corresponds to a situation where the neutronic power reactor suffers unusual unstable problems, i.e., presenting a mix of oscillation modes. The C4_APRM and C4_LPRM_x signals correspond to average power monitors (APRM) and local range monitors LPRM. The entire case 4 consists of a total of 23 signals, 22 LPRMs and 1 APRM.Nonetheless, only the analysis of the APRM is presented in this work. Figure17 shows the studied APRM signal of interest. Figure 18 shows a studied 60 s seconds segment (the dotted red line) that is decomposed through NA-MEMD and its PSD estimate (bottom of figure). Figure 19 shows the first 5 extracted IMFs of the studied segment (as a reminder, the IMFs of the 2 noise extra channels are discarded). Figure 20 shows the PSD estimates of the first 6 IMFs; we highlight that there is a bit of mixing between IMF 4 and IMF 5. Nonetheless IMF 5 is slightly closer to 0.5 Hz than the IMF 4 and has more energy. Figure 21 shows a plot of the IMF 5 of interest that is linked to BWR density wave instability and its associated instantaneous frequency (IF) that oscillates around to 0.5 Hz (we highlight that IMF 5 is slightly mixed with IMF 4 but this is a result due to the complexity of the Forsmark Case 4. Nevertheless, IMF 5 has more energy than IMF 4). The SampEn of this IMF 5 is 0.5443. Thus this segment is unstable (SampEn points to high regularity and low complexity). Figure 22 shows the SampEn estimates of the IMFs ormodes of interest along time for all of the studied segments of 60 s. All of the SampEn estimates are smaller than 0.7. So, the studied APRM signal is unstable.Figure 17
Forsmark Case 4 APRM signal.Figure 18
Studied segment of the Forsmark Case 4 APRM.Figure 19
A plot of the first 5 IMFs in a segment of the studied Forsmark APRM.Figure 20
PSD estimate of the first 6 IMFs extracted in the analyzed segment.Figure 21
IMF 5 linked to BWR instability and its associated IF (oscillating close to 0.5 Hz).Figure 22
SampEn estimates of the IMFs ormodes of interest along time.Table3 shows the mean(SampEn), median(SampEn), and std(SampEn) of the time series shown in Figure 22 (SampEn estimates of IMF 4). The mean(IF) and std(IF) are also given in this table. Figure 23 shows a plot of the computed mean SampEn estimates of the 23 unstable signals (each SampEn estimate is a red point). A pattern appears; all of the mean SampEn estimates for the studied signals are smaller than 0.8. Thus, the mean SampEn values point to visually regular signals such as cyclic functions. Through the proposed methodology, we are able to classify these 23 signals in the unstable category with success because their SampEn is smaller than 0.8.Table 3
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed APRM. Now the mean(IF) and std(IF) are included.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) Mean(IF)[Hz] std(IF)[Hz] Forsmark APRM 0.5369 0.5295 0.0411 0.5017 0.0234Figure 23
Computed mean(SampEn) of the studied 23 unstable signals (SampEn of the IMFs linked to instability) of the studied Forsmark Case 4.


### 4.5. Case 9 Cycle 14 of the Ringhals Stability Benchmark

The data of this case comes from measurements of the Swedish BWR reactor Ringhals 1 [6]. This case consists of a total of 72 LPRMs distributed in twofloorsorlevelsof 36 LPRMs each. As stated before, only the analysis of one signal (LPRM 1) is detailed in this work. Figure 24 shows a plot of the studied LPRM 1 signal (of level 4). Figure 25 shows a plot of the studied segment (the red dotted line) that is decomposed through NA-MEMD and its PSD estimate (bottom of figure). Figure 26 shows the first extracted 5 IMFs (by the way, IMF 4 looks like a cyclic function) of the segment. Figure 27 shows a plot of the PSD estimates of the first 6 IMFs. The PSD estimate of IMF 4 is highly concentrated around 0.5 Hz. Figure 28 shows the tracked IMF 4 of interest linked to BWR instability, this IMF 4 looks like a cyclic function and its associated IF is almost a line centered around 0.5 Hz. The computed SampEn of this IMF 4 is 0.4544 (the studied IMF 4 is unstable. Thus, the studied segment is unstable). Figure 29 shows the computations of the SampEn along time for all of the studied segments of 60 s. All of the estimates are smaller than 0.8 and located between 0.4-0.6. So, the studied signal is unstable and highly regular like a cyclic function.Figure 24
Ringhals Case 9 cycle 14 LPRM1.Figure 25
Studied segment of the Ringhals Case 9 cycle 14 LPRM 1.Figure 26
A plot of the first 5 IMFs extracted from the analyzed segment of the studied Ringhals case.Figure 27
PSD estimate of the first 6 IMFs extracted from the analyzed segment.Figure 28
IMF 4 linked to BWR instability and its associated IF (very close to 0.5 Hz).Figure 29
SampEn estimates of the IMFs ormodes (IMFs 4) of interest along time.Table4 shows the mean(SampEn), median(SampEn), and std(SampEn) of the time series shown in Figure 29. The mean(IF) and std(IF) of all of the studied segments of the targeted LPRM are also given in this table. Figure 30 shows the computed mean(SampEn) values of the 36 LPRMs of level 2 of the Ringhals reactor. It is observed, that all of the SampEn estimates are smaller than 0.8 (The percentage of classification of the studied signals in the unstable category is of 100%). The entire floor 2 is unstable and the SampEn estimates are located in the range 0.5-0.7. Figure 31 shows the computed mean(SampEn) values of the 36 LPRMs of level 4 of the Ringhals reactor. It is observed that all of the SampEn estimates are smaller than 0.8 (The percentage of classification of the studied signals in the unstable category is of 100% as well). The entire floor 4 is also unstable and the SampEn estimates are located in the range 0.45-0.7. The entire case 9 is unstable and SampEn was perfectly able to detect that behavior.Table 4
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed LPRM 1. Now the mean(IF) and std(IF) are included.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) Mean(IF)[Hz] Std(IF)[Hz] Ringhals LPRM 1 0.4668 0.4730 0.0497 0.5165 0.0328Figure 30
Mean(SampEn) values for the 36 LPRMs of level 2 of the Ringhals BWR unit.Figure 31
Mean(SampEn) values for the 36 LPRMs of level 4 of the Ringhals BWR unit.


## 4.1. Stable Signals

The aim of this test is to observe the demeanor of the BWR in steady-state conditions. The results obtained are presented in Figures2–8. Figure 2 shows the plot of one APRM stable signal of a typical BWR. Figure 3 shows in red dotted line, a studied segment of 60 s of duration of the APRM signal with its power spectral density (PSD) estimate (bottom of figure). Figure 4 shows the first four extracted IMFs of the studied signal, which happens to be noise that is decomposed in a quasidyadic fashion by the NA-MEMD. Figure 5 shows the power spectral density (PSDs) of some IMFs.Figure 2
Stable signal from a typical BWR.Figure 3
Analyzed segment of the studied APRM signal.Figure 4
A plot of the first 4 extracted IMFs in the analyzed segment.Figure 5
PSD estimate of the first 6 extracted IMFs in the analyzed segment.Figure 6
Tracked IMF (IMF 3) of interest and its associated instantaneous frequency (IF) around 0.5 Hz computed through the Hilbert-Huang transform.Figure 7
SampEn along time for the studied APRM stable signal.Figure 8
Computed mean(SampEn) of the studied 208 stable signals.Figure6 shows the tracked IMF (IMF 3) of interest whose instantaneous frequency (IF) oscillates around 0.5 Hz. The computed SampEn of this IMF is 0.9351. SuchhighSampEn value indicates that the studied signal is BBNCN, so the studied segment is stable and far from an unstable behavior. In this stable case, amodeexist whose IF is close to 0.5 Hz. However, the IMF 3 time series looks distant from a cyclic function we look for when an unstable event triggers in the core. All of the IFs were computed through the Hilbert-Huang transform of the extracted IMFs. Figure 7 shows the computed SampEn along time for all of the studied segments of 60 s; most of the SampEn estimates have a value slightly greater than 1 (i.e., the studied signal is stable). Table 1 shows the mean(SampEn), median(SampEn) and std(SampEn) of the time series shown in Figure 7.Table 1
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed APRM.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) APRM 1.1288 1.0994 0.1398Figure8 shows a plot of the computed mean SampEn estimates of the 208 signals (each SampEn estimate is a blue point). A pattern appears, all of the mean SampEn estimates for the studied signals are higher than 1. Thus, the mean SampEn values point to BBNCN. Through SampEn estimation, the proposed stability methodology is able to detect stable signals and classify them in the stable category if their SampEn is greater than 1.


## 4.2. Unstable Signals

The studied signal in this subsection stems from an instability event. This instability case is registered by 96 LPRM signals. For reasons of space, only the analysis of one signal is presented. Figure9 shows the plot of the studied signal (from LPRM 2); this signal looks like a cyclic function. Figure 10 shows a studied segment of 60 s (red dotted line) and its PSD estimate (bottom of figure), the segment looks like a noisy cyclic function. Figure 11 shows the first 4 extracted IMFs of the studied segment; the first 3 IMFs are linked with acquisition noise, whereas IMF 4 looks like a cyclic function (this is in fact the type of waveform commonly associated with density wave unstable events). Figure 12 shows the PSD estimates of the first 6 IMFs, the energetic content of IMF 4 is highly concentrated around 0.5 Hz, and the energetic content of the other IMFs is meaningless next to the energy of IMF 4. Figure 13 shows the time domain series of IMF 4 (which happens to be themodeassociated to BWR instability) and its associated IF that oscillates around 0.5 Hz. The IF was computed via HHT. Now, this IF (IF 4) carries physical meaning [7]. The SampEn estimate of IMF 4 (IMF ormodelinked to instability, see IF in Figure 13 and PSD estimate in Figure 12) is 0.6181; this estimate is a clear indication of signal regularity and low complexity. The SampEn points to a simple function; in this case a cyclic function that in our context is linked with a density wave type of instability. Thus, the studied segment is unstable. Figure 14 shows the SampEn estimates along time for all of the studied segments of the signal of interest; all estimates are smaller than 1.Figure 9
Unstable signal from a typical BWR.Figure 10
Analyzed segment of the studied LPRM 2 signal.Figure 11
A plot of the first extracted 4 IMFs in the analyzed segment.Figure 12
PSD of the first 6 extracted IMFs in the analyzed segment. The PSD of the IMF 4 is highly concentrated around 0.5 Hz.Figure 13
Tracked IMF (IMF 4) of interest linked to BWR instability and its associated instantaneous frequency (IF) around 0.5 Hz computed through the HHT.Figure 14
SampEn estimates of the IMFs ormodes of interest along time.Table2 shows the mean(SampEn), median(SampEn), and std(SampEn) of the time series shown in Figure 14. Now, the tracked IF (IF 4) along time carries physical meaning. Thus, the mean(IF) and std(IF) values of all of the studied segments are now provided in this table.Table 2
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed LPRM 2. Now the mean(IF) and std(IF) are included.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) Mean(IF) [Hz] Std(IF) [Hz] LPRM 2 0.5425 0.5258 0.0544 0.5374 0.0030Figure15 shows a plot of the computed mean SampEn estimates of the 96 unstable signals (each SampEn estimate is a red point). A pattern appears, all of the mean SampEn estimates for the studied signals are smaller than 0.8. Thus, the mean SampEn values point to visually regular signals such as cyclic functions. Through SampEn estimation, we are able to build a binary classifier that separates BBNCN signals (linked to stability) from cyclic ones (linked to density wave instability). One core feature of SampEn is that it accommodates for complex nonlinear and nonstationary data whereas conventional Decay Ratio estimates must assume beforehand that the studied BWR signals behave linearly (an assumption that is false in real life).Figure 15
Computed mean(SampEn) of the studied 96 unstable signals.


## 4.3. Comparison between Stable and Unstable States

To complete the analysis of the stable/unstable signals through SampEn, Figure16 shows a comparison of the mean SampEn values of the studied sets (stable set of 208 signals plus a set of 96 unstable signals). The computed SampEn mean values of the stable signals fluctuate very close to 1.2 whereas the SampEn mean values of the unstable signals oscillate very close to 0.6. So the SampEn detects with success stable BWR signals from unstable ones via complexity analysis (of course with the aid of the NA-MEMD to denoise these signals and isolate themodelinked to instability). There is a distance of ≈0.6 between the two states (stable from unstable). So, by fixing a threshold value around 0.9, it is possible to differentiate one state from the other properly. Any SampEn estimate of a BWR signal segment higher than 0.9 is stable whereas any lower SampEn value from this threshold is unstable.Figure 16
Comparison of the mean(SampEn) values for both studied sets (a stable set of 208 APRM signals versus an unstable set of 96 LPRM signals).


## 4.4. Case 4 of the Forsmark Stability Benchmark

Case 4 of the Forsmark stability benchmark [5]. This event is considered a difficult case to study by the complexity of the phenomena. For reasons of space, only the analysis of one case will be presented in detail. This case contains a mixture between a global oscillation and a regional oscillation. This particular event corresponds to a situation where the neutronic power reactor suffers unusual unstable problems, i.e., presenting a mix of oscillation modes. The C4_APRM and C4_LPRM_x signals correspond to average power monitors (APRM) and local range monitors LPRM. The entire case 4 consists of a total of 23 signals, 22 LPRMs and 1 APRM.Nonetheless, only the analysis of the APRM is presented in this work. Figure17 shows the studied APRM signal of interest. Figure 18 shows a studied 60 s seconds segment (the dotted red line) that is decomposed through NA-MEMD and its PSD estimate (bottom of figure). Figure 19 shows the first 5 extracted IMFs of the studied segment (as a reminder, the IMFs of the 2 noise extra channels are discarded). Figure 20 shows the PSD estimates of the first 6 IMFs; we highlight that there is a bit of mixing between IMF 4 and IMF 5. Nonetheless IMF 5 is slightly closer to 0.5 Hz than the IMF 4 and has more energy. Figure 21 shows a plot of the IMF 5 of interest that is linked to BWR density wave instability and its associated instantaneous frequency (IF) that oscillates around to 0.5 Hz (we highlight that IMF 5 is slightly mixed with IMF 4 but this is a result due to the complexity of the Forsmark Case 4. Nevertheless, IMF 5 has more energy than IMF 4). The SampEn of this IMF 5 is 0.5443. Thus this segment is unstable (SampEn points to high regularity and low complexity). Figure 22 shows the SampEn estimates of the IMFs ormodes of interest along time for all of the studied segments of 60 s. All of the SampEn estimates are smaller than 0.7. So, the studied APRM signal is unstable.Figure 17
Forsmark Case 4 APRM signal.Figure 18
Studied segment of the Forsmark Case 4 APRM.Figure 19
A plot of the first 5 IMFs in a segment of the studied Forsmark APRM.Figure 20
PSD estimate of the first 6 IMFs extracted in the analyzed segment.Figure 21
IMF 5 linked to BWR instability and its associated IF (oscillating close to 0.5 Hz).Figure 22
SampEn estimates of the IMFs ormodes of interest along time.Table3 shows the mean(SampEn), median(SampEn), and std(SampEn) of the time series shown in Figure 22 (SampEn estimates of IMF 4). The mean(IF) and std(IF) are also given in this table. Figure 23 shows a plot of the computed mean SampEn estimates of the 23 unstable signals (each SampEn estimate is a red point). A pattern appears; all of the mean SampEn estimates for the studied signals are smaller than 0.8. Thus, the mean SampEn values point to visually regular signals such as cyclic functions. Through the proposed methodology, we are able to classify these 23 signals in the unstable category with success because their SampEn is smaller than 0.8.Table 3
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed APRM. Now the mean(IF) and std(IF) are included.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) Mean(IF)[Hz] std(IF)[Hz] Forsmark APRM 0.5369 0.5295 0.0411 0.5017 0.0234Figure 23
Computed mean(SampEn) of the studied 23 unstable signals (SampEn of the IMFs linked to instability) of the studied Forsmark Case 4.


## 4.5. Case 9 Cycle 14 of the Ringhals Stability Benchmark

The data of this case comes from measurements of the Swedish BWR reactor Ringhals 1 [6]. This case consists of a total of 72 LPRMs distributed in twofloorsorlevelsof 36 LPRMs each. As stated before, only the analysis of one signal (LPRM 1) is detailed in this work. Figure 24 shows a plot of the studied LPRM 1 signal (of level 4). Figure 25 shows a plot of the studied segment (the red dotted line) that is decomposed through NA-MEMD and its PSD estimate (bottom of figure). Figure 26 shows the first extracted 5 IMFs (by the way, IMF 4 looks like a cyclic function) of the segment. Figure 27 shows a plot of the PSD estimates of the first 6 IMFs. The PSD estimate of IMF 4 is highly concentrated around 0.5 Hz. Figure 28 shows the tracked IMF 4 of interest linked to BWR instability, this IMF 4 looks like a cyclic function and its associated IF is almost a line centered around 0.5 Hz. The computed SampEn of this IMF 4 is 0.4544 (the studied IMF 4 is unstable. Thus, the studied segment is unstable). Figure 29 shows the computations of the SampEn along time for all of the studied segments of 60 s. All of the estimates are smaller than 0.8 and located between 0.4-0.6. So, the studied signal is unstable and highly regular like a cyclic function.Figure 24
Ringhals Case 9 cycle 14 LPRM1.Figure 25
Studied segment of the Ringhals Case 9 cycle 14 LPRM 1.Figure 26
A plot of the first 5 IMFs extracted from the analyzed segment of the studied Ringhals case.Figure 27
PSD estimate of the first 6 IMFs extracted from the analyzed segment.Figure 28
IMF 4 linked to BWR instability and its associated IF (very close to 0.5 Hz).Figure 29
SampEn estimates of the IMFs ormodes (IMFs 4) of interest along time.Table4 shows the mean(SampEn), median(SampEn), and std(SampEn) of the time series shown in Figure 29. The mean(IF) and std(IF) of all of the studied segments of the targeted LPRM are also given in this table. Figure 30 shows the computed mean(SampEn) values of the 36 LPRMs of level 2 of the Ringhals reactor. It is observed, that all of the SampEn estimates are smaller than 0.8 (The percentage of classification of the studied signals in the unstable category is of 100%). The entire floor 2 is unstable and the SampEn estimates are located in the range 0.5-0.7. Figure 31 shows the computed mean(SampEn) values of the 36 LPRMs of level 4 of the Ringhals reactor. It is observed that all of the SampEn estimates are smaller than 0.8 (The percentage of classification of the studied signals in the unstable category is of 100% as well). The entire floor 4 is also unstable and the SampEn estimates are located in the range 0.45-0.7. The entire case 9 is unstable and SampEn was perfectly able to detect that behavior.Table 4
Mean(SampEn), Median(SampEn), and std(Sampen) of the analyzed LPRM 1. Now the mean(IF) and std(IF) are included.
Signal Mean(SampEn) Median(SampEn) Std(SampEn) Mean(IF)[Hz] Std(IF)[Hz] Ringhals LPRM 1 0.4668 0.4730 0.0497 0.5165 0.0328Figure 30
Mean(SampEn) values for the 36 LPRMs of level 2 of the Ringhals BWR unit.Figure 31
Mean(SampEn) values for the 36 LPRMs of level 4 of the Ringhals BWR unit.


## 5. Conclusions

In this work the Sample Entropy, a measure that provides an index of signal complexity or irregularity of a time series, was used to estimatecomplexity of a BWR time series for stability analysis. The SampEn was tested in conjunction with the noise-assisted multivariate empirical mode decomposition, a nonlinear filter that decomposes nonstationary data from non-linear sources. The NA-MEMD is used in this application to extract correctly the associated mode to the density wave oscillation, before to compute the SampEn.The proposed SampEn + NA-MEMD rolling window methodology allowed categorize stable BWR signals from unstable ones. The computed SampEn values of thestablesignals oscillated around 1.2 (ahighSampEN value is linked to a broad band noncoherent noise (BBNCN) associated with stable behavior). The SampEn values ofunstable signals fluctuate close to 0.6 (they might be lower than that, but overall such SampEn estimates are below the threshold SampEn value of 0.8). Alowvalue of SampEn points in this case to a simple cyclic function that is linked to density wave instability. 100% of the typical BWR studied cases were classified with success into stable or unstable categories. 100% of the Forsmakr Case 4 studied LPRMs confirm BWR instability through our proposed SampEn + NA-MEMD methodology. 100% of the studied LPRMs of level 2 and 4 of the Ringhals Case 9 cycle 14 confirm instability. According to these experiments, the SampEn paired with the NA-MEMD is a suitable candidate to be used as a BWR stability indicator.In the future, the SampEn might be tested in conjunction with other suitable nonlinear indicators (possibly fractal dimension estimators [25]) to build a robust BWR stability monitor based on 2 nonlinear measurements, where one of them is going to be the SampEn. In general, the SampEn is a suitable and easy to implement BWR nonlinear stability indicator for real time series analysis.



---

*Source: 9852925-2018-11-01.xml*