Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Dynamic Simulation and Control of a Supercritical CO2 Power Conversion System for Small Light Water Reactor Applications Shih-Ping Kao and Jonathan Gibbs Massachusetts Institute of Technology 77 Massachusetts Avenue, 24-215 Cambridge, MA 02139-4307 USA [email protected], [email protected] Pavel Hejzlar Massachusetts Institute of Technology 77 Massachusetts Avenue, 24-215 Cambridge, MA 02139-4307 USA [email protected] Abstract - Strategies to control a supercritical CO2 (S-CO2) power conversion system (PCS) have been investigated for light water reactor applications. Due to sharp changes in density and thermal properties near the pseudo-critical region, controlling the S-CO2 Brayton cycle presents unique challenges in nuclear power applications. A dynamic simulation code (SCPS), based on real-gas and momentum integral models, has been developed at MIT to evaluate control strategies for a small pressurized light water reactor equipped with a compact S-CO2 PCS. The NIST REFPROP subroutines and tables were integrated with the code to evaluate the CO2 properties at each timestep. The single-shaft, S-CO2 PCS consists of a radial compressor and turbine power train and three HEATRICTM Printed Circuit Heat Exchangers (PCHE) used as intermediate, recuperator, and pre-cooling heat exchangers. The primary system is that of a scaled-down commercial PWR. The simulation results show an oscillatory heat transfer behavior near the pseudocritical region inside the precooler, similar to that observed in an experiment conducted by the Seoul National University. Consequently, the control programs were developed to maintain compressor inlet conditions away from and above the pseudo-critical region, thus avoiding heat transfer oscillations and compressor surge instabilities. An integrated control system, which couples the reactor power control with turbine throttle and bypass control, has been developed and shown to successfully control a small pressurized light water reactor coupled to a compact SCO2 PCS for a wide range of operational transients.

ABSTRACT Strategies to control a supercritical CO2 (S-CO2) power conversion system (PCS) have been investigated for light water reactor applications. Due to sharp changes in density and thermal properties near the pseudo-critical region, controlling the S-CO2 Brayton cycle presents unique challenges in nuclear power applications. A dynamic simulation code (SCPS), based on real-gas and momentum integral models, has been developed at MIT to evaluate control strategies for a small pressurized light water reactor equipped with a compact S-CO2 PCS. The NIST REFPROP subroutines and tables were integrated with the code to evaluate the CO2 properties at each timestep. The single-shaft, S-CO2 PCS consists of a radial compressor and turbine power train and three HEATRICTM Printed Circuit Heat Exchangers (PCHE) used as intermediate, recuperator, and pre-cooling heat exchangers. The primary system is that of a scaled-down commercial PWR. The simulation results show an oscillatory heat transfer

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

behavior near the pseudo-critical region inside the precooler, similar to that observed in an experiment conducted by the Seoul National University. Consequently, the control programs were developed to maintain compressor inlet conditions away from and above the pseudocritical region, thus avoiding heat transfer oscillations and compressor surge instabilities. An integrated control system, which couples the reactor power control with turbine throttle and bypass control, has been developed and shown to successfully control a small pressurized light water reactor coupled to a compact S-CO2 PCS for a wide range of operational transients. KEYWORDS supercritical CO2, power conversion system, light water reactor

1.

INTRODUCTION

The supercritical CO2 (S-CO2) Brayton cycle has been the focus of Gen-IV applications at MIT for several years [1]. Work has focused on a multi-application power conversion system (PCS) rated at 300 MW electric. This project, however, is focused on a PCS rated between 5 and 30 MW electric [2]. Scaling to these sizes may have significant effects on key design parameters, including turbomachinery design, component efficiency and system pressure. For example, radial (centrifugal) compressors are often better suited than their axial flow counterparts for lower throughput applications. Also, radial compressors typically have a much wider off-normal operating range between stall and surge, which is beneficial (if not required) for applications in which load-follow may be more important than commercial baseload-only performance. Another important difference between the large and medium power systems are high-pressure piping size limitations. CO2 has a relatively small specific heat (about four times less than water) and requires relatively large mass flow rates, and thus large high pressure pipes to maintain reasonable pressure drops. Since compressor work is proportional to pressure drops, the optimal power rating is limited to 200 to 300 MW electric, by optimizing cycle efficiency while minimizing the compressor work. PCHE sizing and pressure drops are provided by the manufacturer, thus pipe size becomes a design parameter. The systems in the tens of megawatt range are not constrained by this limit and provide much more freedom for a compact high-performance design. The objective of this study was to determine the dynamic performance and establish the fundamental design sensitivities and limitations for a S-CO2 closed Brayton cycle PCS in the range of 5 to 30 MW electric. This paper presents the results of the following tasks of this study:

2.

•

Development of a dynamic simulation code (SCPS) of the compact S-CO2 PCS cycle in a configuration where both the turbine and the compressor are mounted on a single shaft.

•

Evaluation of various control strategies, part load operation, and cycle response to various transients using the dynamic code.

PWR/S-CO2 PLANT DYNAMIC MODEL

An integrated power plant dynamic simulator using the real-gas and momentum integral numerical method has been developed to simulate the responses of a 10 MW electric PWR

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

equipped with a compact, single-shaft S-CO2 PCS. In this plant simulator, a scaled-down commercial PWR model has been incorporated to simulate the feedback response from the primary system. An integrated control program coupling a reactor power control system with a PCS turbo-compressor control system has been developed to test the controllability of a 10 MW electric PWR/S-CO2 power plant. Simulation transient results covering a wide range of operating conditions have been analyzed in order to achieve optimal performance for a given set of operating conditions. Figure 1 illustrates the conceptual layout of a 10 MW electric, single-shaft compact S-CO2 PCS (the man in the figure is 1.83 m/6 ft tall). Table 1 lists the key performance data of the design.

Figure 1. Single-Shaft Compact S-CO2 PCS Layout [2] Table 1. S-CO2 PCS Performance Characteristics [2] Net cycle efficiency (%) Thermodynamic efficiency (%) CO2 mass flow rate (kg/s) Turbine specific work (kJ/kg) Compressor specific work (kJ/kg) Cooling water inlet temperature (˚C) Cooling water outlet temperature (˚C) Cooling water mass flow rate (kg/s) Cooling water pump work (kW) Recuperator effectiveness (%)

21.28 22.88 253.65 60.14 17.81 29.56 39.26 894.48 169.26 95.1

The primary system model is represented by a lumped-parameter thermal-fluid transport model, consisting of four components: a reactor core, hot-leg/vessel exit plenum, coldleg/vessel inlet plenum, and the primary side of an intermediate heat exchanger (IHX). Reactor power is calculated by a six-group point-kinetics model and a seven-group decay power model. Fuel heat transfer is calculated by a lumped fuel/cladding model for a cylindrical fuel rod with Zircaloy cladding. The primary system pressure and coolant mass flow rate are assumed to be constant. The primary system model represents a first-order approximation with respect to fuel and IHX conduction heat transfer and the thermal fluid heat transport process. The computed dynamic thermal-hydraulic parameters provide the

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

necessary inputs to calculate the reactivity feedbacks found in a typical pressurized water reactor. The secondary PCS model was developed using a more sophisticated approach. A modified integral momentum model was used to calculate the real-gas, thermal-fluid transport process of S-CO2. The model developed represents a single-shaft recuperated PCS, consisting of a generator, a radial turbine and compressor, IHX, recuperator, and precooler heat exchanger of the PCHE type, as shown in Figure 2.

(BYP) 

(GEN) 

Figure 2. PWR/S-CO2 Plant Control Schematic Diagram 2.1. PCHE Heat Transfer Model The highly compact design of HEATRICTM printed circuit heat exchangers (PCHEs) makes them favorable compared to the conventional tube/shell heat exchanger designs for the S-CO2 PCS, due to their significantly higher heat transfer surface area to volume ratio. PCHEs are constructed by diffusion bonding a series of thin metal plates (typically 1.5 mm thick) together, with each plate having zigzag or straight flow channels (typically semicircles of 2 mm diameter) chemically etched on the surface. The hot and cold fluids flow through alternating plates in counter-current flow; however, the exact plate configuration, channel type, and flow arrangement can be customized based on heat transfer and pressure drop constraints. A nodal, lumped-parameter approach [3] has been developed to calculate the heat transfer rate inside a PCHE. The average nodal temperature differences are evaluated along the flow path of the heat exchanger. Figure 3 illustrates the nodalization scheme of a row of single hot and cold channel pairs (unit cell) in the PCHE. The single arrows represent the direction of fluid flow, whereas the double arrows represent conduction and convection heat transfer through the dividing wall.

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

H h ,in

H h,1

H h,1

Tw,1

H c ,1

H c, 2

Tw, N −1

Qc,j

H c , 2 H c , j −1

H c, j

H h, N −1 H h , N −1 H h, N Qh,N-1

Tw, j

Qc,2

H c ,1

H h, j H h, N −2

Qh, j

Tw, 2

Qc,1 c ,0

H h, j

Qh,2

Qh,1

H

H h , 2 H h, j −1

H h, 2

Qc,N-1

H c, j H c, N −2

H h, N

Qh,N

Tw, N Qc,N

H c , N −1 H c , N −1 H c, N

H c ,in

Figure 3. PCHE Model Nodalization Scheme The discretized time-dependent energy conservation equations for the hot and cold side fluid, using a semi-implicit, forward difference approximation, are given by Eqs. (1) and (2):

(ρAc ∆z )h

H hi +, j1 − H hi , j

(ρAc ∆z )c

H ci ,+j1 − H ci , j

∆t

∆t

= m& h (H hi +, j1−1 − H hi +, j1 ) − Qhi , j

(1)

= m& c (H ci ,+j1 − H ci ,+j1−1 ) + Qci , j

(2)

Since for transient conditions, Qh,j ≠ Qc,j , the following lumped-parameter transient wall heat conduction equation is included in the model: ρC p Aw ∆z

Twi,+j1 − Twi, j ∆t

= Qhi +, j1 − Qci ,+j1

(3)

where Q hi = hh′ Ph ∆z (Thi − Twi )

(4)

Qci = hc′ Pc ∆z (Twi − Tci )

(5)

and

The wall temperature is the volume-averaged structure temperature and the heat transfer coefficients are given by: ⎛ 1 c ⎞ hh′ = ⎜⎜ + L ⎟⎟ ⎝ hh 2 k w ⎠

−1

⎛1 c ⎞ hc′ = ⎜⎜ + L ⎟⎟ ⎝ hc 2k w ⎠

−1

(6)

and (7)

The hot and cold side convection coefficients are determined from the Gnielinski correlations for turbulent flow as: fc (Re − 1000 ) Pr 8 h= ⎛ 2 ⎞ f 1 + 12.7⎜⎜ Pr 3 − 1⎟⎟ c ⎠ 8 ⎝

⎛ kf ⎜ ⎜d ⎝ eq

⎞ ⎟ ⎟ ⎠

⎛ d eq ⎞ ⎜⎜1 + ⎟K L ⎟⎠ ⎝

(8)

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

where K is a correction factor to account for wall heating effects and is defined as: ⎛ Pr K = ⎜⎜ ⎝ Prw

⎞ ⎟⎟ ⎠

0.11

; for Pr >1.5,

(9)

and ⎛ Tj K =⎜ ⎜T ⎝ w

⎞ ⎟ ⎟ ⎠

0.45

; for 0.5 < Pr < 1.2

(10)

The subscript w in Eqs. (6) through (10) implies that the fluid properties are evaluated at the wall temperature where the Moody friction factor, fc, is defined as: ⎞ ⎛ 1 ⎟⎟ f c = ⎜⎜ ⎝ 1.8 log Re − 1.5 ⎠

2

(11)

and the hydraulic diameter, deq, is defined as: d eq =

π d c2 ⎞ ⎛ d 2⎜ π c + d c ⎟ ⎠ ⎝ 2

(12)

Note that the Gnielinski correlation represented by Eq. (8) holds true only for straight channels. However, the recuperator (REC) has zig-zag channels. The models for the heat transfer coefficient and pressure drop in zig-zag channels are the same as described in Ref. 1, based on the measurements of HEATRICTM PCHE performance at Tokyo Institute of Technology (TIT) [4]. 2.2. PCS Dynamic Simulation Model The modeling methodology of the PCS simulation code is based on the Momentum Integral (MI) method [7], in which the mass, energy and momentum conservation equations are integrated over the flow channel, eliminating the computational limitations of the sonic effect in the fine-mesh scheme. In this approach only the channel averaged mass flow rate is calculated, and the fluid is assumed to be sectionally incompressible but thermally expandable in each control volume, while the fluid is compressible over the integral control volume. The assumptions applied in this MI model are summarized as follows: •

• • • •

A single system pressure is assumed for the hot and cold side of the system. Thus, two system pressures, one for the hot side (or the turbine discharge side) and the other for the cold side (or the compressor discharge side), are modeled as state variables in the conservation equations for the mass, energy and momentum. Rate of energy change due to viscous dissipation and kinetic energy are negligible Fluid is assumed incompressible but thermally expandable in each control volume. A linear enthalpy profile is assumed inside each control volume. The specific volume of CO2 is assumed to vary linearly with enthalpy (this assumption is very close to reality, even for pressures near the critical point, as shown in Figure 4).

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Figure 4. Specific Volume of CO2 Versus Enthalpy (7.5 – 19.5 MPa) [8] Using the above assumptions and enthalpy as a state variable instead of temperature, overcomes severe numerical difficulties due to the nonlinearity of density near the critical point, as shown in Figure 5. This approach is expected to be more robust and provide a faster solution than the traditional engineering codes.

Figure 5. CO2 Density Versus Temperature Near the Critical Point (7.5 – 11 MPa) [8] The MI mass and energy conservation equations can be expressed as follows [5]:

Where,

dM i = m& i −1 − m& i + m& injection − m& leak dt

(13)

dE i = (m& H ) i −1 − (m& H ) i + (m& H ) injection − (m& H ) leak + Qi ,hx − Qi , wall dt

(14)

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Qi , hx =

heat transfer rate from the heat exchanger;

Qi , wall = conduction heat transfer rate to the structure.

where the total mass and energy content of control volume i, Mi and E, can be defined in terms of volume averaged quantities, mi and ei, as follows: Mi = mi Vi

(15)

Ei = ei Vi - P Vi

(16)

and

and the volume-averaged density, mi, and the volume-averaged specific energy, ei are defined in terms of enthalpy, H, and pressure, P, as: mi =

1 ρ ( P , H ) dV Vi V∫i

(17)

ei =

1 ρ ( P, H ) HdV Vi V∫i

(18)

and

The integrals in Eqs. (17) and (18) can be evaluated analytically if the spatial profile for H within the cell and the relation between density, ρ, and H are known. In this model, it is assumed that the transient is slow enough that a linear profile always holds for enthalpy, H. Therefore, by a linear transformation, mi and ei can be expressed as functions of the inlet and exit enthalpies and pressure: ∫H i −1 ρ ( P , H )dH Hi

mi =

H i − H i −1

(19)

and ∫H i −1 ρ ( P, H ) HdH Hi

ei =

H i − H i −1

(20)

Applying the approximation that the specific volume of CO2, v i , varies linearly with enthalpy, Hi, mi and ei can be solved analytically as follows:

v ⎞ ln⎛⎜ i vi −1 ⎟⎠ ⎝ mi = vi − vi −1 and

⎡ H − H i −1 ⎤ ⎡ (H i −1vi − H i vi −1 )⎤ ⎛ vi ⎞ ⎟⎟ ei = ⎢ i ⎥ln⎜⎜ ⎥+⎢ 2 ⎣ vi − vi −1 ⎦ ⎣ (vi − vi −1 ) ⎦ ⎝ vi −1 ⎠

(21)

(22)

The closed-loop momentum conservation equation, integrated over the flow components of turbine (TRB), main compressor (CMP), and turbine bypass valve (BYP) in terms of the volume-averaged mass flow rate, < m& i >, leads to the following:

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

⎛L⎞

∑ ⎜⎝ A ⎟⎠

CMP

d < m& > CMP d < m& > TRB d < m& > BYP ⎛L⎞ ⎛L⎞ +∑⎜ ⎟ +∑⎜ ⎟ dt dt dt ⎝ A ⎠ TRB ⎝ A ⎠ BYP

(23)

= − ∑ Fi (< m& > i ) + ∑ ρ i g sin θ i + J CMP − J TRB i

i

Where, ⎡ f ( L / Dh ) + K f ⎤ ⎥ < m& > i < m& > i ; 2 ρA 2 ⎣ ⎦i

∑ F (< m& > ) = ⎢ i

i

< m& > i =

i

1 m& dL; Li ∫i

f = friction factor; K f = form loss coefficient; < m& >i = < m& >TRB, < m& >CMP, or < m& >BYP; and

J TRB = characteristic head of turbine. J CMP = characteristic head of compressor. The characteristic performance map is a function of normalized volumetric flow rate (or mass flow rate), rotor speed, and fluid density ratios, such as: ⎛^ ^ ^⎞ J TRB , J CMP = F (m& m& ref , Ω Ω ref , ρ ρ ref ) = F ⎜ m& , Ω, ρ ⎟ ⎠ ⎝

(24)

The assumptions implied by the momentum conservation equation are: • •

The pressure loss in each flow component is small compared to the system pressure of the hot or cold side; and Fluid is incompressible within the hot and cold sides of the system.

The system dynamic model for a single-shaft, variable-speed turbine-compressor can be derived from the angular momentum conservation equation given by the following equation: ⎛ dΩ ⎞ I shaft ⎜ ⎟ = TTRB − TCMP − TGEN ⎝ dt ⎠

(25)

Where, Ishaft = moment of inertia of the shaft; TTRB = torque generated by the turbine; TCMP = torque applied by the main compressor; and TGEN = torque applied by the electric generator load. The torques for the turbine and main compressor are given by the characteristic functions expressed in terms of normalized shaft speed, Ω volumetric flow rate (or mass flow rate, m& ), and fluid density, ρ, ratios, as follows:

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

⎛^ ^ ^⎞ T = F (m& m& ref , Ω Ω ref , ρ ρ ref ) = F ⎜ m& , Ω, ρ ⎟ ⎝ ⎠

(26)

2.3. Turbomachinery Performance Maps Since the turbomachinery characteristic maps for the S-CO2 PCS were not available, a set of performance maps developed for a generic centrifugal compressor [9] were incorporated into the SCPS code. Assuming the fluid is relatively incompressible, compressor pressure rise can be scaled as a function of reference pressure rise, ∆Pref , density, ρ, and shaft speed, Ω, as follows: 2 ρ exit Ω ∆P = ∆Pref ρ exit ,ref Ω 2ref

(27)

The compressor mass flow rate can be scaled as a function of mass flow rate, m& ref, density, ρ, and shaft speed, Ω, as follows:

m& = m& ref

ρexit , ref Ω ref ρexit Ω

(28)

A polynomial fit of the generic centrifugal compressor maps yields the following equation for the compressor pressure rise [9]: Y = −4.5942( X ) + 10.487( X ) − 7.9392( X ) + 3.0484 3

2

Where, m& X =

ρ exit ,ref Ω ref ρ exit Ω m& ref

(29)

(30)

and Y =

∆P ρ exit ,ref Ω 2ref ∆Pref

ρ exit

Ω2

(31)

Similarly, the radial compressor efficiency is estimated by the following polynomial function: Yη = −1.4386( X ) + 3.2804( X ) − 2.2997( X ) + 1.4599 3

2

Where, Yη =

η η ref

(32)

(33)

Similarly for the turbine, because the performance maps for the axial turbine were not available during the project, a pressure head versus mass flow rate relationship is used to model the turbine head given by:

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

∆PTRB = KTRB m& 2 / ATRB

(34)

The compressor work is given by the following relation:

Q& CMP =< m& CMP > J CMP / ρ CMP / η CMP

(35)

And the turbine work is given by:

Q& TRB =< m& TRB > ( H TRB ,in − H TRB ,out ,isentropic )ηTRB

(36)

2.4. PCS State-space Model Applying a form of the donor-cell differencing method that assumes any change in a cell is propagated uniformly between the inlet and exit [10], the energy conservation equation of Eq. (14) can be expressed as follows, (assuming m& injection << m& i ): ' ⎛ ∂E * ⎞ dH i ⎛ ∂E ⎞ dPi ⎟⎟ + ⎜⎜ = m& i ( H i−1 − H i ) + m& injection ( H injection − H i* ) + Q& i ,hx − Q& i ,wall ⎜ ⎟ ∂ ∂ H dt P dt ⎝ ⎠i ,i−1 ⎝ ⎠i

(37)

where, H i* = 1 / 2( H i + H i −1 );

(38)

⎡ ρ i H i − ρ i −1 H i −1 ⎤ ⎛ ∂E * ⎞ ⎟ = Vi ⎢ ⎜ ⎥ ∂ H ⎠ i ,i −1 ⎝ ⎣ H i − H i −1 ⎦

(39)

Similarly, the mass conservation equation of Eq. (13) can be expressed as follows: ⎛ ∂M * ⎞ dH i ⎛ ∂M ⎞ dPi ⎜⎜ ⎟⎟ +⎜ = m& i −1 − m& i ⎟ ⎝ ∂H ⎠ i ,i −1 dt ⎝ ∂P ⎠ i dt

(40)

⎛ ρ − ρ i −1 ⎞ ⎛ ∂M * ⎞ ⎛ ∂M ⎞ ⎛⎜ ∂M i ⎞⎟ ⎟⎟ = Vi ⎜⎜ i =⎜ ⎟ ⎟ +⎜ ⎜ ⎟ ⎝ ∂H ⎠ i ,i −1 ⎝ ∂H ⎠ i ⎝ ∂H i −1 ⎠ ⎝ H i − H i −1 ⎠

(41)

where,

Applying Eulerian discretization on Eqs. (37) and (40) yields a system of equations consisting of six energy and two system mass conservation equations as follows: Gn ( yn+1 - yn ) = dn where, y= [H TRB , H CMP , H IHX −CS , H REC − HS , H REC −CS , H PRE −CS , PHS , PCS ] G= the Jacobian of y dn= right-hand sides of Eqs. (37) and (40).

(42)

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

The matrix-vector system of equations (42) is solved using a Gaussian elimination process with a fixed time-step size. The closed-loop momentum conservation equation given by Eq. (23) is linearized and solved explicitly for the average mass flow rate: .

⎧⎛ L ⎞ ⎛L⎞ ⎫ ⎛L⎞ ⎛L⎞ ⎛L⎞ ⎛L⎞ + ⎜ ⎟ + ⎜ ⎟ ⎬∆ m& +⎜ ⎟ +⎜ ⎟ +⎜ ⎟ ⎨⎜ ⎟ ⎩⎝ A ⎠ IHX − cs ⎝ A ⎠ PRE − hs ⎝ A ⎠ REC − hs ⎝ A ⎠ REC − cs ⎝ A ⎠TRB ⎝ A ⎠CMP ⎭ ⎡ ⎤ = ∆t ⎢− ∑ Fi ( m& ) + ∑ ρi g sin θi + J CMP − JTRB ⎥ i ⎣ i ⎦

(43)

The linearization of the angular momentum equation results in the following equation: I shaft ∆Ω = ∆t [CTRB ∆TTRB − CCMP ∆TCMP − CGEN ∆TGEN ]

(44)

where, CTRB =

∂T ∂T ∂TTRB ; C CMP = CMP ; C GEN = GEN . ∂Ω ∂Ω ∂Ω

(45)

Equations (43) and (44) are solved explicitly for the system average mass flow rate and shaft speed. 3. CONTROL MODELS FOR PWR/S-CO2 PLANT Control models for a PWR/S-CO2 plant have been developed based on a load-following scheme similar to that of a commercial PWR. Upon receipt of an electric-load demand signal, the integrated control systems attempt to control the reactor power to produce the generator power matching the demand. This strategy is based on the fact that the primary reactor system has a much larger thermal inertia than the S-CO2 PCS. Therefore, the primary control system is a lead or master controller for the integrated system, while the PCS control system is the slave following the primary control actions. 3.1. Reactor Power Control Model The reactor power control system (RCS) is similar to that of the control rod system designed for a commercial PWR. That is, the reactor power is raised or lowered by inserting or withdrawing the control rods. The control program is simplified by assuming that the integral control rod reactivity is linearly proportional to the rod position. Based on the Doppler and moderator temperature feedbacks in reactor kinetics, the controller program accounts for power defects due to changing fuel and coolant temperatures. In this investigation, the control-rod power control program is designed to respond to fast transients, such as sudden loss of load or power increase. Long-term power operations such as startup and shutdown, where decay power and Xenon and Samarium poisoning must also be accounted for, are not taken into consideration. Figure 6 illustrates the control block diagram for the reactor power control system. The reactor fission power (Qfission) is compared against the power demand (Qdemand), where the power error is filtered through an impulse compensator. The power error is summed with a

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

primary system temperature error compared at the demand power level to yield a total power error signal. Based on this total error, the control rod movement rate is determined from a linear rate program.

Figure 6. Reactor Power Control Block Diagram 3.2. S-CO2 PCS Control Model In order to make PCS power follow reactor power, a turbine power controller is used to regulate the PCS flow through the turbine by changing the control turbine throttle valve position. Figure 7 shows the control block diagram for the turbine power controller. The controller is a classical PI controller, using demand versus turbine power error as the input.

Figure 7. PCS Turbine Throttle Controller Block Diagram

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

In order to prevent shaft overspeed in a sudden loss of load event, a turbine bypass controller is used to open a flow path bypassing the turbine, reducing the turbine torque applied on the shaft, thereby slowing down the shaft speed. A classical PI controller is used to cycle open a turbine bypass valve. Figure 8 illustrates the block diagram for the turbine bypass controller. Due to limited capacity of the bypass flow, designed at 20% of the total flow capacity of the PCS (20% was found to provide sufficient bypass flow, combined with turbine throttle control, to accommodate 100% loss of load without exceeding the turbine speed limit), the turbine bypass controller alone cannot be used to prevent shaft overspeed in loss of load events greater than 20% of full power.

Figure 8. Turbine Bypass Controller Block Diagram The turbine bypass controller has a total of five possible bypass paths or configurations: • • • • •

Turbine bypass: IHX cold side outlet to recuperator hot side inlet. IHX bypass: Recuperator cold side outlet to recuperator hot side inlet. Precooler dump: IHX cold side outlet to precooler hot side inlet. Recuperator bypass: Recuperator cold side outlet to precooler hot side inlet. Cycle bypass: Compressor outlet to precooler hot side inlet.

It has been observed that unstable heat transfer oscillations can occur near the pseudo-critical region of CO2, which may result in large density and pressure head changes in the compressor. In order to avoid this unstable region, a flow controller is used to regulate the precooler cooling-water flow rate to maintain the S-CO2 temperature above 31ºC at the precooler exit. Figure 9 shows the block diagram for the precooler temperature controller. As shown, a flow controller is used to regulate the flow control valve position controlled by a signal from the temperature deviation. This control mechanism regulates cooling water flow to control the heat transfer coefficient on the water side of the precooler; therefore, indirectly controlling the CO2 temperature exiting the precooler.

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Figure 9. Precooler Temperature Controller Block Diagram 4. STABILITY ANALYSIS IN THE VICINITY OF THE CO2 CRITICAL POINT Component test transient runs prior to incorporating control mechanisms have been conducted to test the response of an uncontrolled S-CO2 PCS to step changes in boundary conditions. Test results have shown that the PCS responds very rapidly to the introduced changes, with cycle thermal response time constants between 8 to 12 seconds to reach new steady state conditions for a 20 MW electric PWR/S-CO2 plant. However, when the precooler temperature is in close proximity of the pseudo-critical temperatures and CO2 pressure is close to the critical pressure, oscillations occurred. The pseudo-critical ranges are considered as ±0.05 MPa and ±0.2 °C of the critical pressure and temperature, respectively. These oscillations are shown for the case of a 3 °C step decrease of precooler water inlet temperature (at time 0 seconds) as shown in Figure 10 and Figure 11. The oscillations occur near the phase transition between vapor and liquid in close proximity to the critical point. The oscillations are found to be inside the 4th node closest to the pseudo-critical temperature. Using more spatial nodes has been found to increase the magnitude of oscillations since larger variations in properties are captured with finer nodalization. o Pre-cooler Node CO2 Heat Transfer Rate (Case 1: Precooler H2O Inlet Precooler Node CO C) -3C) 2 Heat Transfer Rate (Case 1: Precooler H 2O Inlet ‐3

30

PRE-QH4

Heat Transfer Rate (MW)

PRE-QH3 PRE-QH2 25

PRE-QH1

20

15

10 -1.0

1.0

3.0

5.0

7.0

9.0

11.0

13.0

15.0

Time (s)

Figure 10. Precooler Nodal CO2 Side Heat Transfer Rate

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

o Pre-cooler Node Pressure (Case 1: Precooler H2O Inlet Precooler Node Pressure (Case 1: Precooler H C) -3C) 2O Inlet ‐3

7.60

PRE-PRS-H4 PRE-PRS-H3

Pressure (MPa)

7.55

PRE-PRS-H2 PRE-PRS-H1

7.50 7.45 7.40 7.35 7.30 -1.0

1.0

3.0

5.0

7.0

9.0

11.0

13.0

15.0

Time (s)

Figure 11. Precooler Nodal CO2 Side Pressure The sharp pseudo-critical peak in Cp does not occur only above the critical point, but also below the critical point in liquid and vapor phase regions as long as the state of the fluid is in close proximity of the critical point. The CO2 thermal capacity in this region is extremely sensitive to pressure and temperature changes as shown by the NIST plots in Figure 12. This property anomaly presents a numerical challenge, that is, an isobaric temperature change of 0.1˚C can change the thermal capacity by an order of magnitude. Similarly, a small change in pressure on the order of 0.01 MPa can also change the thermal capacity by an order of magnitude at constant temperature.

Figure 12. Isobaric CO2 Cp vs. Temp. (7.35 – 7.4 MPa at 5KPa intervals) [8] It has been calculated that a small reduction in CO2 temperature from 30.7538°C to 30.7527°C (with pressure of 7.339 MPa just below the critical pressure of 7.377 MPa) results in a significant reduction of CO2 heat transfer coefficient (by 40%), thus reducing the heat transfer rate to the H2O side. This affects CO2 temperature and density and thus system pressure due to CO2 expansion or contraction. Because of the very large change in thermal capacity of CO2 near the critical point, a small change in temperature and pressure near the critical point significantly affects the CO2 properties and thus heat transfer rate, providing a positive feedback. A similar type of wall temperature oscillation has been reported in the

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

laboratory experiments on a S-CO2 loop at Seoul National University [6]. Because of these oscillations and the difficulties with compressor design near the critical point, it has been decided to increase compressor inlet temperature/pressure to 35.5°C/9.24 MPa. The PCS was re-optimized for these new conditions and all subsequent control studies were performed at these new reference conditions. There is a small penalty on cycle efficiency when compressing at a higher temperature, but this can be easily compensated with a small increase of turbine inlet temperature. 5. INTEGRATED TRANSIENT ANALYSIS A series of ten transient tests were performed to test the robustness of the SCPS code simulating the dynamic responses of a coupled plant consisting of a PWR primary system and S-CO2 PCS under a wide range of operating conditions. The transient tests were performed for three categories of transients. In the first category, partial load rejection tests designed to demonstrate the load-leading and load-following capabilities of a 10 MW electric plant were executed. Figure 13, Figure 14, and Figure 15 show the plant responses for shaft speed, compressor pressure ratio, and reactor core and generator power output versus demand, respectively, following a 40% loss of generator load with turbine bypass and turbine throttle control. The test results have demonstrated that the controller actions were able to bring the plant to a new steady-state in about 80 seconds, with 7% shaft overspeed. The second test category assessed the potential of the system to accommodate power changes in the load-leading mode of operation where the generator load follows the reactor and turbine power driven by the reactor and turbine power controllers. To demonstrate the load-leading control characteristics, a 20% step power demand swing with turbine bypass control was performed, where a 20% step decrease in power demand was followed by a 20% step increase in power demand 30 seconds later. The results (shown in Figure 16, Figure 17, and Figure 18) confirmed that the controllers were able to bring the plant to a new steady state with small deviations from the new power demands, which varied between 1.5 to 7.5%. These deviations could be reduced by further tuning of controller parameters. The third category of tests involved transients with sudden changes in boundary conditions, namely 1) the cooling water temperature step was increased and decreased to test the precooler cooling water flow controller’s ability to prevent the PCS from entering into the pseudo-critical region where instabilities can occur, and 2) discharge operation of the inventory control system at full power to test the plant and control systems response to a reduction in CO2 inventory. Figure 19, Figure 20, and Figure 21 show the results from a transient under sudden 3 ºC decrease in precooler cooling water temperature. The results indicate that the controllers are able to maintain the PCS CO2 temperature above the pseudocritical region thus avoiding the undesirable heat transfer instability mentioned earlier. The last scenario investigated was a 100% loss of load (LOL) to test the controllers’ ability to control shaft overspeed under the most limiting operational conditions. The controller configuration employed turbine bypass and turbine throttle for this test. The transient response for shaft speed is shown in Figure 22, and the compressor exit and inlet pressure ratio is shown in Figure 23 and Figure 24, respectively. The results confirm that turbine throttle runback and turbine bypass control working together can effectively limit the shaft overspeed to under 15%. It was also found that turbine bypass control alone cannot accommodate the 100% load rejection because the compressor starts to choke due to a large flow rate increase.

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Therefore, bypass control in combination with turbine throttle control needs to be used. This method was found to be effective in both full load rejection and partial load rejection events. Moreover, partial load rejection cases can be controlled solely by turbine throttle and do not need bypass control. For the 10 MW electric PWR/S-CO2 power plant, both load-following operations, where the reactor and turbine power controllers react to power mismatch signals as a result of external load change, and load-leading operations, where the generator load follows the reactor and turbine power, were investigated. In the latter mode, the reactor and turbine power controllers adjusted power based on the prescribed power demand. However, the traditional control strategy based on pressure changes upstream of the turbine throttle was not possible because the stored energy in CO2 is very small due to the compactness of the cycle and the small heat capacity of CO2. Finally, it is important to note that all the simulations and analyses of PCS control were performed using approximations of turbomachinery off-design characteristic curves, based on typical ideal gas machines, because the compressor off-design curves for S-CO2 turbomachinery were not available. Normalized shaft speed (40% LOL)

RPM (norm.)

1.08 1.06 1.04 1.02 1.00 0.98 0

20

40

60

80

100

Time (sec)

Figure 13. Normalized Shaft Speed (40% LOL) Compressor Pressure Ratio (40% LOL)

Pressure Ratio

3.0

2.5

2.0

1.5 0

20

40

60

80

Time (sec)

Figure 14. Compressor Pressure Ratio (40% LOL)

100

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Core & Generator Power vs. Demand (40% LOL)

Normalized Power

1.2 Core Generator Demand

1.0 0.8 0.6 0.4 0.2 0

20

40

60

80

100

Time (sec)

Figure 15. Normalized Reactor and Generator Power Vs. Demand (40% LOL) Normalized Shaft Speed (20% Power Swing)

RPM (norm.)

1.05

1.00

0.95

0.90 0

20

40

60

80

100

Time (sec)

Figure 16. Normalized Shaft Speed (20% Power Swing) Compressor Pressure Ratio (20% Power Swing)

Pressure Ratio

2.5 2.4 2.3 2.2 2.1 2.0 1.9 0

20

40

60

80

100

Time (sec)

Figure 17. Compressor Pressure Ratio (20% Power Swing)

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Core & Generator Power vs. Demand (20% Power Swing)

Normalized Power

1.2 1.1 1.0 0.9 Core Generator Demand

0.8 0.7 0.6 0

20

40

60

80

100

Time (sec)

Figure 18. Normalized Reactor and Generator Power vs. Demand (20% Power Swing) Normalized Shaft Speed (Cooling Water Temp. -3 ºC Change)

RPM (norm.)

1.050 1.025 1.000 0.975 0.950 0

20

40

60

80

100

Time (sec)

Figure 19. Normalized Shaft Speed (Cooling Water Temp. -3 ºC Change) Compressor Pressure Ratio (Cooling Water Temp. -3 ºC Change)

Pressure Ratio

2.6

2.4

2.2

2.0 0

20

40

60

80

100

Time (sec)

Figure 20. Compressor Pressure Ratio (Cooling Water Temp. -3 ºC Change)

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Norm. Core & Gen. Power vs. Demand (Cooling Water Temp. -3 ºC Change)

Normalized Power

1.10 Core Generator Demand

1.05 1.00 0.95 0.90 0

20

40

60

80

100

Time (sec)

Figure 21. Normalized Reactor and Generator Power vs. Demand (Cooling Water Temp. -3 ºC Change) Normalized Shaft Speed (100% LOL, Turbine Bypass) 1.2 RPM (norm.)

1.1 1.0 0.9 0.8 0.7 0.6 0

10

20

30

40

50

60

Time (sec)

Figure 22. Normalized Shaft Speed (100% LOL)

Compressor Exit Pressure (100% LOL, Turbine Bypass)

Pressure (MPa)

20.2 20.0 19.8 19.6 19.4 19.2 19.0 0

10

20

30

40

50

Time (sec)

Figure 23. Compressor Exit Pressure (100% LOL)

60

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Compressor Inlet Pressure (100% LOL, Turbine Bypass)

Pressure (MPa)

9.4 9.2 9.0 8.8 8.6 0

10

20

30

40

50

60

Time (sec)

Figure 24. Compressor Inlet Pressure (100% LOL)

6.

CONCLUSIONS

Supercritical CO2 is a unique working fluid that can be used in a Brayton cycle to offer a more compact PCS footprint as compared to a Rankine cycle or a Brayton cycle using an ideal gas working fluid. However, due to its large variations in thermal properties near the critical point, which can result in undesirable thermal and pressure oscillations, the S-CO2 PCS must be designed in such way that the operation in the pseudo-critical region is avoided. In this research, a set of control programs have been designed and tested for a hypothetical 10 MW electric PWR/S-CO2 PCS plant to demonstrate the feasibility of such applications. Simulation test results from a wide range of operating conditions have shown that an integrated control system employing a reactor power controller coupled with turbine throttle, bypass control, and cooling water temperature control mechanisms, can successfully control the plant to avoid the regions of instabilities. However, feasibility of this application needs to be analyzed again when representative S-CO2 turbomachinery performance data become available. ACKNOWLEDGMENTS The authors would like to thank Lockheed Martin for facilitating this research project. Special thanks to Drs. Nathan Carstens and Yifang Gong, and Mr. Dustin Langewisch for their contributions in this project. Finally, the authors are grateful for the project guidance of Professor Michael J. Driscoll. REFERENCES 1. 2. 3. 4.

Hejzlar P. et al. “Supercritical CO2 Brayton Cycle for Medium Power Applications”, Final Report, MIT-ANP-PR-117, MIT Center for Advanced Nuclear Energy Systems, April 2006. Hejzlar P. et al. “Dynamic Evaluation of Supercritical CO2 Brayton Cycle for Medium Power Applications”, Final Report, MIT-ANP-TR-121, MIT Center for Advanced Nuclear Energy Systems, September 2008. Dustin Langewisch, Shih-Ping Kao, M.J. Driscoll, “HXMOD – A Heat Exchanger Module for Power Cycle Dynamics Modeling”, MIT-GFR-040, August, 2006. Ishizuka T., Y. Kato, Y. Muto, K. Nikitin, N.L. Tri, and H. Hashimoto, “Thermal Hydraulic Characteristics of Printed Circuit Heat Exchanger in a Supercritical CO2 Loop”,

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Proc. Of NURETH 11, Avignon, France, October 2-6, 2005. Kao, S. P., “A Multiple-Loop Primary System Model for Pressurized Water Reactor Sensor Validation”, Ph. D. thesis, MIT, 1984. 6. Kim, J. K., et. al., “Wall Temperature Measurement and Heat Transfer Correlation of Turbulent Supercritical Carbon Dioxide Flow in Vertical Circular/Non-Circular Tubes,” NURETH-11, Avignon, France, October 2–6, 2005, Nuclear Engineering & Design, Vol. 237 (2007). 7. Meyer, J. E., “Some Physical and Numerical Considerations for the SSC-S Code,” BNLNUREG-50913 (1978). 8. NIST Reference Fluid Thermodynamic and Transport Properties Database (REFPROP): Version 7.0, http://www.nist.gov/srd/nist23.htm, 2007. 9. Carstens, N. A., “Control Strategies for Supercritical Carbon Dioxide Power Conversion Systems,” Ph. D. Thesis, Dept. of Nuc. Sci. & Eng., MIT, June 2007. 10. Weaver, W. L., et. al, “A Few Pressure Model for Transient Two-Phase Flows in Networks, “ Trans. ANS, 28, 273-274, 1981. 5.

NOMENCLATURE A Cp cL dc deq e E f g h H H I J k K L m& m P Ph Pc Q t T T V v W

flow area, cross sectional area (m2) specific heat capacity (J/kg-K) conduction thickness (m) channel diameter (m) hydraulic diameter (m) specific energy (J/kg) energy (J) Moody friction factor gravity acceleration constant (kg/m-s2) heat transfer coefficient (W/m2-K) enthalpy (J/kg) node-average enthalpy (J/kg) moment of inertia (kg-m2) characteristic head (Pa) thermal conductivity (W/m-K) correction factor to account for wall heating effects length (m) mass flow rate (kg/s) volume-averaged density (kg/m3) pressure (Pa) heated perimeter (m) cooled perimeter (m) heat rate (W) time (s) temperature (K), or torque (N-m) node-average temperature (K) volume (m3) specific volume (m3/kg) work (J)

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009

Greek letters ρ ρˆ ∆z ∆t Ω ˆ Ω

ε η

density (kg/m3) normalized density node length (m) time step size (s) shaft speed (RPM) normalized shaft speed heat transfer convergence band (%) efficiency (%)

Subscripts c cold, channel cs cold side eq equivalent h hot hs hot side i node number w wall ref reference Superscripts n nth time step Abbreviations BYP turbine bypass CMP compressor GEN generator IHX intermediate heat exchanger PI proportional-integral PCS power conversion system PRE precooler RCS reactor coolant system REC recuperator TRB turbine

Proceedings of S-CO2 Power Cycle Symposium 2009 RPI, Troy, NY, April 29-30, 2009