Natural Circulation Characteristics at Low-Pressure Conditions through PANDA Experiments and ATHLET Simulations

Natural circulation characteristics at low pressure/low power have been studied by performing experimental investigations and numerical simulations. The PANDA large-scale facility was used to provide valuable, high quality data on natural circulation characteristics as a function of several parameters and for a wide range of operating conditions. The new experimental data allow for testing and improving the capabilities of the thermal-hydraulic computer codes to be used for treating natural circulation loops in a range with increased attention. This paper presents a synthesis of a part of the results obtained within the EU-Project NACUSP “natural circulation and stability performance of boiling water reactors.” It does so by using the experimental results produced in PANDA and by showing some examples of numerical simulations performed with the thermal-hydraulic code ATHLET.


INTRODUCTION
In the framework of the EU-Project NACUSP, experiments on thermo-hydraulic characteristics of natural-circulationcooled boiling water reactors (BWR) have been performed in four sophisticated thermo-hydraulic test facilities which complement each other, ranging from small-scale to largescale, and from low-pressure/low-power operating conditions to the nominal operating conditions for BWRs.The ultimate goal of the project was to improve the economics of operating and future plants through improved operational flexibility, enhanced availability, and increased confidence levels on the safety margins regarding the stability issues in BWRs [1].
This paper presents a synthesis of the results obtained from the experiments performed in the large-scale PANDA facility [2] and from the numerical simulations performed with the thermal-hydraulic code ATHLET.The range of parameters covers a large spectrum of conditions for the low-pressure/low-power range, therefore, being of interest especially with regard to the start-up procedures for naturalcirculation-cooled BWRs (e.g., for the European simplified boiling water reactor (ESBWR)).

PANDA facility
The multipurpose thermal-hydraulic test facility PANDA is located at Paul Scherrer Institute (PSI), Villigen, Switzerland.The facility is designed and used for investigating at large-scale system behavior and phenomena for different light water reactor (LWR) designs.PANDA has a modular structure, which is based on six cylindrical pressure vessels with a total volume of 460 m 3 and four open pools with a total capacity of 60 m 3 (Figure 1).Four of these pressure vessels are arranged in two vertical columns.The two lower vessels simulating the wetwell or suppression chamber are interconnected with two large pipes.The wetwell vessels support the two vessels, simulating the drywell.They are interconnected at about midplane elevation with one large pipe.The fifth vessel simulates the core flooding (GDCS) pool.The sixth vessel is configured to simulate the reactor pressure vessel (RPV).Its dimension allows for performing natural circulation tests at roughly 1 : 1 scale in height.The vessel has an inner diameter of 1.23 m and is 19.2 m in height.The RPV includes a heated section, a riser, and a downcomer (DC).The heated section has a height of 1.3 m and does not reproduce the geometry of a real core.The power is generated by 115 electrical heater elements divided into six individually controlled groups.This allows for conducting tests with any specified power up to 1.5 MW.The heaters are distributed over the cross-section to obtain a uniform radial distribution over the whole power range.The riser has a height of 9.5 m (starting at the top of the heated section).Four rectangular open pools are placed on top of the drywell vessels and are equipped with heat exchangers/condensers.The facility is well insulated in order to minimize heat losses: an insulation layer of 300 mm is applied on the RPV and 100 mm on system lines and pools.The overall height of the facility is about 25 m.The maximum operating conditions are 10 bar and 200 • C.
The basic instrumentation of the facility includes about 1000 sensors to measure temperatures, pressures, pressure differences, water levels, flow rates, gas concentrations, fluid phases, and electrical heater power.The data acquisition system maximum scanning rate is 5 Hz.The standard operating scanning rate is 0.5 Hz.The facility is remotely controlled from a graphical-display man-machine interface with process visualization.

PANDA facility configuration for natural circulation tests
For the NACUSP project only parts of the PANDA facility were used.probe is used to measure the level swelling due to the void production in the RPV.

PANDA test matrix
The test matrix allowed the RPV power and pressure to be varied, as well as other parameters influencing the natural circulation behavior such as core-inlet hydraulic resistance and RPV water level.The tests were performed at constant power, balanced by a corresponding condenser heat removal capacity.The specified RPV power was selected to match the desired steady state pressure in the RPV.The basic test matrix is shown in Table 1.For a given RPV pressure level, the maximum energy is removed if the condenser is fully submerged in the pool.The corresponding power for the four specified pressure levels is given in Table 1 (column "high").For running tests at lower power, the condenser heat removal capacity is adjusted by lowering the water level in the pool, thereby reducing the active condenser area.Three tests with different power ("low," "middle," and "high") were performed for each of the four different pressure levels.
The core-inlet hydraulic resistance coefficient k has been varied as follows: for the basic tests a BWR typical value was used (k basic = 30).Selected tests were performed with a low inlet resistance (k low = 7) and with high inlet resistance (k high = 500).The RPV water level has also been varied: most of the tests were performed with the collapsed water level at 12.8 m above bottom of RPV.This level represents about the nominal value for the ESBWR.For few tests the water level was reduced close to top of riser which is at 11.0 m above RPV bottom.The selected values were 11.1 m and 11.4 m.
The following three series of tests with totally 25 experiments have been performed:  An overview of all PANDA tests is included in Table 2.

Natural circulation modes
The coolant temperature increases in the heated section (core), but it may not reach saturation conditions under the given conditions.A rough estimation of the position of the boiling boundary in the core by means of a heat balance was performed.The temperature measured at the inlet of the core was used as the reference temperature, and the coarsest assumptions were made (no heat losses, uniform radial and vertical distribution of power in the core, no subcooled boiling).The nondimensional value z * boil.b reported in Table 2 corresponds to the distance (measured from the inlet of the core, and divided by the length of the core) at which this boiling boundary can be expected (z * boil.b≤ 1: boiling in the core; z * boil.b> 1: no boiling in the core).The result is that boiling in the core did not occur in any test except may be in tests with high core-inlet resistance and relatively high power for which z * boil.b∼ = 1 is found.In the same way, a coarse estimation of the height in the riser at which flashing may be expected to occur has been made.Using the same inlet temperature, calculating the temperature at the outlet of the core, and assuming it is   conserved all along the riser (no heat losses), a distance z * flash (measured from inlet of riser, and divided by length of riser) is estimated and reported in Table 2.For most of the tests, it is expected that flashing should only occur in the upper region or above the riser.Of course, in tests at high core-inlet resistance and relatively high power, and also in tests with low RPV water level (and hence lower pressure head), this basic calculation predicts flashing in lower sections of the riser.
These estimations do not take into account the 3D and local effects that might affect the physical phenomena and the flow regime actually occurring in the RPV.Some interesting information can be retrieved by plotting some of the temperature profiles measured in the central axis of the RPV.Two examples are given in Figure 3, for which, a curve corresponding to the saturation temperature is also plotted on each graph.The saturation temperature has been estimated by using the time-averaged RPV pressure and calculating the pressure head at each height, without considering the possible presence of void.
From Figure 3, it is clear that the validity of the comparison with the saturation temperature may be limited due to the error in temperature measurements.Due to the data reduction procedure, this error slightly increases with temperature (Figure 3).The saturation temperature decreases along the vertical axis as the system operates at low pressure.Flashing should occur when the coolant reaches saturation conditions in the riser.Following this approach, in test H2.1 (shown in Figure 3(a)), the average height of the flashing boundary can be roughly estimated to be at about the level corresponding to the top of the riser.In this case, flashing should be responsible for the shape of the profile above this height, where measured temperatures seem to lie in the vicinity of the saturation curve.This basically confirms the previous coarse prediction (Table 2, z * flash = 0.95).In test H8.1, a different situation is shown (Figure 3 Rough estimations of time-averaged void fraction values processed from differential pressure measurements also indicate that the highest void fraction was assessed at the top of the riser, and especially for cases with high coreinlet resistance/high power (Class 3).These cases showed the largest difference between the swell level (measured by the TDR probe) and the collapsed water level retrieved from differential pressure measurements.For example, in tests HL5.3 and H8.3, a void fraction of the order of 10% was estimated in the upper region of the riser.In the lower sections, values lower than a few percent were found.

Analysis of flow velocity measurements in RPV downcomer
The natural circulation flow rate is of main interest.Hence, this short analysis concerns the velocity measurements from 3 ultrasonic flow meters in the RPV DC.The sensors are located at elevation 5.00 m above RPV bottom with an angle of 120 • azimuthally between each other.
The velocity signals were sampled at a frequency of 0.5 Hz for the whole duration of the test (5 hours).It should be noted that the time constant of the sensors was set to 6 seconds in order to avoid aliasing problems.The autopower spectral density (APSD) of the signals was calculated to identify resonance frequencies of the system.The autocorrelation function (ACF) was used to calculate the decay ratio (DR) of the system.The DR is a widely used parameter to quantify the stability of the system: if DR < 1, the system is linearly stable; if DR > 1, it is unstable.In practice, the DR values presented in Table 2 were estimated by simply calculating the ratio between the second and the first maxima of the ACF (Figure 4(c)).Cross-power spectral densities, coherence, and cross-correlation functions were also calculated to examine coherence, possible phase shifts, and so forth between signals from different flow meters.
An overview of the results of the analysis performed can be found in Table 2.The basic parameters of each test are given in the first columns of this table.Mean values of the velocities calculated over the test period (V 1 , V 2 , and V 3 ) are also reported.V 123 is the average of the three mean velocities.The quantity I 123 should provide an estimate of the intensity of the velocity fluctuations, and was calculated as follows: where Var 123 represents the average of the 3 variances of the signals.It should be realized that this relation could result in very high values for the tests with very low natural circulation velocities, for which the signal-to-noise ratio is expected to be higher (e.g., test H2.3).
Concerning the spectral analysis, the frequency and the corresponding period of the major peak are reported when possible.The value "0" in the column "Δϕ" of Table 2 means that no significant phase shift between the signals was observed, that is, the three velocity signals were oscillating in phase at one given frequency.When possible, a DR was estimated from the ACF and is also reported in Table 2.One column of the table refers to the analysis of signals from the sensor measuring the total pressure in the RPV.These additional observations may help to characterize some of the tests.Looking at the spectra and at the results of this analysis based on the ACF and the estimation of DR values, the PANDA tests can be categorized in two different groups as follows.

1st group
The rows in Table 2 corresponding to these tests are greyed.An illustration of this case is made using test L8.3 results.An excerpt of the raw velocity signals recorded during this test is shown in Figure 4(a).A major peak can be observed on the processed spectra of the 3 sensors reported in Figure 4(b).The time traces of the velocity signals show that the three sensors oscillate "together," that is, with no significant phase shift between each other.This was confirmed by looking at the phase of the complex cross-power spectral densities between the different signals.Just simply looking at a number of raw signals recorded over the test period, it can also be seen that the other measured variables (temperatures, pressure in the RPV, condenser feed, and drain flows, etc.) show fluctuations with the same frequency.The ACF obtained from velocity signal 1 is given in Figure 4(c).Clearly, a DR value indicating stable behavior (DR = 0.64) can be estimated from this graph.

2nd group
An illustration of this group is given using test B3.3 results.Observations made on the raw signals (Figure 5(a)) in this case show rather random behavior.In some other tests belonging to the same group, it is sometimes possible to observe significant phase shifts between the sensors.Figure 5(b) shows that there is no longer a major peak present in the power spectrum.Using a logarithmic scale, the low-frequency part of the spectrum appears to be flatter.Moreover, the ACF does not really allow a DR value to be derived (Figure 5(c)), this DR appearing to be much smaller (actually, very close to 0) compared to the tests of the first group.
It was also observed that these two characteristic modes of behavior, indicated by the ACF plots, did not vary as a function of time during the course of the tests.For tests L8.3 and B3.3, this was done by calculating ACFs over 5 successive (1-hour long) time periods and it clearly showed that, in both cases, all the ACFs for one test at different times were very similar.
Some data from other sensors (condenser feed and drain flow rates, RPV pressure, etc.) have also been analyzed using the same approach as for the velocity data.For the tests from the 1st group, a resonance peak can be observed at the same oscillation frequency.No peak can be found for the tests belonging to the 2nd group.However, in some of these tests (e.g., test B3.3), the condenser feed flow and RPV pressure signals may exhibit a peak in their spectrum whereas no such peak was observed for the velocity signals (see also the second last column of Table 2, concerning the analysis of RPV pressure signals).This indicates that flashing possibly takes place above the top of the riser, not influencing the main circulation flow.These cases fall consistently into Class 1 (see Figure 3).
From the expected mechanism of the flow oscillations, the fluid travelling time in the heated section, and the transit time of the enthalpy perturbations in the adiabatic section, should play a role in determining the period of the oscillations.From Figure 11, it is clear that the oscillation period decreases with increasing average flow rate.The total transit time can be evaluated as the sum of the transit time in the riser section plus half of the transit time in the heated section [3].In PANDA, no velocity measurement is available in the riser.An estimation of the riser velocity was made by simply using the time-averaged value of the velocity in the DC reported in Table 2, and just considering the flow area ratio between the riser and the DC.
Following this simplified approach, the oscillation period was found to be in the order of 1.6 times the transit time of the fluid (Figure 6).The simple way in which the transit time was estimated did not consider the occurrence of twophase flow in the riser.Hence, it is likely that our assessment overestimated this time.Thus, these observations seem to be quite consistent with those presented for the CIRCUS facility [1], for which the period of oscillation was found to be twice the transit time.

Discussion
The differentiation between the two groups of tests is based on DR values derived from the velocity signals.However, the validity and the accuracy of this method can be questioned.Many other methods have been used to extract stability parameters and can be compared to the ACF-based method.For this reason, some of our experimental velocity time traces were processed using a more sophisticated method, namely, the autoregressive method.This method properly applied to our signals yielded about the same DR values as those presented in Table 2.This confirms the main result from this analysis that all tests are stable.
However, by having a closer look to the tests, the different flow regimes that can be assumed, and the three classes defined, can be related to the previous findings regarding the stability of the system.Information reported in Table 2 was used for this purpose.The plots in Figure 3 partly illustrate the following classification.Class 4 (Figure 3(a)).High inlet subcooling and relatively low (or even very low) flow rates characterize this class.It is likely that flashing occurred only above the riser, the natural circulation flow still being in a single-phase regime.This would explain why, in some of these tests, contrary to the velocity signals in the DC, some RPV pressure signals show a peak in their spectrum and allow the derivation of DR values (see, e.g., test B5.2).It is also possible that in some of these tests (notably those performed at the lowest powers), no flashing occurred at all in the whole system.
Class 5.These tests, performed at relatively high power and low inlet resistance, have been characterized by a wellidentified oscillation period and DR values could be derived from the velocity signals.The highest natural circulation flow rates were measured in these tests, which indicate the combined effect of the flow enhancement by steam production and of a low inlet resistance.According to our estimations and analyses, no boiling occurred in the core but flashing did in the upper half of the riser.
Class 6 (Figure 3(b)).Very stable behavior was found for these tests performed with an increased inlet resistance.Under these conditions, the flashing front fluctuates close to the bottom of the riser, which stabilizes very much the two-phase natural circulation.The measured flow rates are among the lowest of the test series, as a result of the high flow resistance.Looking at the parameters (notably the expected position of the boiling boundary) reported in Table 2 and at the characteristic temperature profiles (Figure 3(b)), it seems that stable two-phase natural circulation flow can be assumed in these cases.
The presented classification gives a rather coherent picture of the PANDA test series.It was possible to link the experimental observations with the expected phenomenology.No unstable tests were recorded.However, the characterization of the tests (from single-phase to two-phase flow conditions) is satisfying and quite consistent with what is reported in the literature.Test BLL5.3, which shows a DR of about 0.85, should be very close to unstable conditions.
The test matrix allowed the variation of the RPV water level.However, the actual influence of this parameter on the flow behavior cannot be clearly assessed.Comparing tests B5.3 and BL5.3, it seems that a lower level would have a stabilizing effect.

Thermal-hydraulic model
The thermal-hydraulic code ATHLET [4], which has been developed by GRS (Gesellschaft für Anlagen-und Reaktorsicherheit mbH), was used for the calculation of selected experiments from the PANDA test matrix (see Section 2).The ATHLET input dataset models all main parts of the PANDA configuration used for the NACUSP natural circulation experiments (Figure 7).
The ATHLET model consists of the lower plenum (P1-LP-1), the core section (P1-CORE) with 115 electrical heater elements, the riser (P1-RIS1), upper plenum (P1-UP-1), downcomer (P1-DC), and the upper part of the RPV (P1-RIS2).The isolation condenser (IC) is not modeled.Therefore, a bypass (CIRC-ENT, CIRC-EIN) is connected to the upper part of the RPV and to the downcomer, modeling the IC-feed and drain lines by fills with constant mass flows and given enthalpy (only drain line).In all control volumes, the 5-equation model (separate conservation equations for liquid and vapour mass and energy, mixture momentum equation) and the full-range drift-flux model of ATHLET are used.A valve at the lower end of the downcomer (VLV, see Figure 7) is used to model different k-factors of the core inlet.The cross-section of this valve can be changed to adjust the core inlet flow resistance.

Steady state calculation
At first, a steady state calculation with constant boundary conditions was performed for each PANDA test.With help of these calculations, the pressure losses, heat losses, the RPV water level, and the temperature distribution were adjusted.
The steady state calculation starts with zero power.After a few seconds, the core power is switched on and the power is increased with time.For the calculations, it was not possible to specify a certain k-factor for the core inlet, because the form-loss coefficients are changed within the algorithm of the ATHLET code during the steady state calculation.Therefore, a valve at the lower end of the downcomer with changing cross-section was used to adjust the core-inlet flow resistance.The cross-section of this valve was reduced to a value, which leads to DC velocities corresponding to the measured data.The steady state calculation is stopped if stable conditions are reached.As an example, the results of the steady state calculation for Test B 8.3 and the corresponding experimental data are shown in Figure 8.

Transient calculations
From the PANDA test matrix, 5 experiments were selected for the ATHLET simulations.Table 3 shows the initial conditions calculated with ATHLET in comparison to the measured data.For the experiments the data were recorded over a period of 5 hours.The experimental data in Table 3 represent average values.The initial conditions for the ATHLET simulations correspond to the end of the steady state calculations.Each transient calculation starts from a steady state calculation as described in Section 3.2.To simulate the system behavior with respect to natural circulation stability a short disturbance of the drain line mass flow was used in order to stimulate oscillations in the loop (see Figure 9).The response of the system can be seen in the calculated riser mass flow.
As in the experiments, the DC velocities were used to calculate the spectrum and the decay ratio.Figure 10 shows the downcomer velocity, the spectrum, and the autocorrelation function for two ATHLET calculations (Test B 8.3 and Test B 3.3).Both calculations predict stable behavior and no limit-cycle oscillations occur.The disturbance of the drain line mass flow leads to flow oscillations with decreasing   amplitude.In both cases, the oscillations were damped due to the higher core-inlet loss coefficient (k = 30) used in the B-series experiments.Although the oscillations were damped in the ATHLET simulations, an oscillation period and also a decay ratio can be calculated from the DC mass flow.Table 4 gives a comparison between measured and calculated period and decay ratio for the selected PANDA tests.The experiments and also the ATHLET calculations show that the oscillation period decreases with increasing DC velocity.This behavior is illustrated in Figure 11.The results of the ATHLET calculations show a good agreement with the experimental data.
To demonstrate the ability of ATHLET to calculate an unstable behavior with the PANDA model, an additional calculation with a core-inlet loss coefficient of k = 7, a higher core power, and a lower RPV level was performed.The initial conditions for this calculation are given in Table 3  strong oscillations caused by flashing in the riser section.
Figure 12 shows the calculated DC mass flow, the spectrum, and the autocorrelation function.The oscillations have a period of 62 seconds with a decay ratio >0.9.
In Figure 13, the calculated temperature distribution over the height of the PANDA model is illustrated for t = 5000 seconds.At this time, the dynamic behavior becomes unstable.The subcooling in the core is high and no void production takes place.Only above the middle of the riser section, the temperature reaches saturation conditions and flashing-induced oscillations occur.

CONCLUSIONS
In the PANDA facility, 25 tests have been carried out in order to investigate natural circulation flow behavior and stability, under low-pressure/low-power conditions.Some parameters influencing the natural circulation flow were varied, such as the RPV power and pressure, the core-inlet hydraulic resistance, and also the water level in the RPV.
An analysis of the data from three ultrasonic sensors installed at three circumferential locations in the downcomer was presented.The power spectra show, in a few tests, a major and unique resonance peak.The period of oscillation of the flow rate seems to be of the order of twice the transit time of the fluid in the riser.Analyses of the time series show that these cases are stable.The tests for which more flat spectra were obtained show a behavior that is even more stable.
A phenomenological classification has been applied according to three circulation modes that could be assumed: single-phase circulation (with possibly flashing above the riser), two-phase circulation with flashing in the upper half of the riser, or two-phase flow along the complete riser length.The new experimental data about the behavior of the natural circulation flow under low-pressure conditions represent a valuable extension of the available database, notably because the height of the riser in PANDA is approximately preserved with respect to a real reactor.The data have been used for assessing the capabilities and limitations of thermal-hydraulics codes in predicting natural circulation characteristics at low pressure/low power.The ATHLET simulations show stable behavior for the selected PANDA tests.In all simulations, the oscillations were damped and no limit-cycle oscillations occur.The calculated periods and decay ratios show a good agreement with the experimental results.With help of the additional calculation, it could be demonstrated that ATHLET can also simulate an unstable behavior.In this case, strong oscillations occur, caused by flashing in the riser section.

Figure 2 :
Figure 2: PANDA configuration for natural circulation tests.

Figure 3 :
Figure 3: Measured temperature profiles in the vertical axis of the RPV.

Figure 6 :
Figure 6: Period of oscillation versus fluid transit time in core and riser.
(b)): the fluid temperature clearly following the saturation curve indicates that two-phase flow should have occurred along the complete riser length.To summarize, these different plots help to qualitatively distinguish different cases and define three classes: Class 1. flashing above the riser, or perhaps no flashing at all (Figure 3(a)); Class 2. flashing at a lower elevation in the riser; Class 3. two-phase flow all along the riser (Figure 3(b)).

Figure 7 :
Figure 7: Schematic representation (a) and ATHLET model (b) of PANDA test facility.

Figure 10 :
Figure 10: DC velocity, spectrum, and autocorrelation function for the ATHLET calculations of tests B 8.3 and B 3.3.

Figure 11 :
Figure 11: Oscillation period as a function of the DC velocity for the PANDA tests and the ATHLET simulations.

Figure 12 :Figure 13 :
Figure12: DC velocity, spectrum, and autocorrelation function for the additional ATHLET calculation with higher core power and lower RPV level.

Table 1 :
Basic test matrix (variation of pressure and power).

Table 2 :
Overview of the PANDA tests and analysis results.

Table 4 :
Oscillation period and decay ratio calculated from the DC mass flow (Exp.= experiment, Cal.= ATHLET calculation).
. The results of the simulation show an unstable behavior with