Transient analysis and control of a heat to power conversion unit based on a simple regenerative supercritical CO2 Joule-Brayton cycle

Abstract Supercritical carbon dioxide (sCO2) heat to power systems are a promising technology thanks to their potential for high efficiency and operational flexibility. However, their dynamic behaviour during part-load and transient operation is still not well understood and further research is needed. Additionally, there is not enough literature addressing suitable control approaches when the objective is to follow the dynamics of heat load supplied by the topping process to maximise the power recovery. The current research aims to fill these gaps by proposing a one-dimensional transient modelling formulation calibrated against the major components of a 50 kWe sCO2 test facility available at Brunel University London. The dynamic analysis showed that the system quickly adapts to a 2800s transient heat load profile, proving the flexible nature of the sCO2 system investigated. The turbine by-pass, during startup and shutdown modes of operation, enabled gradual and safe build-up/decline of the pressures and temperatures throughout the loop. The regulation of the inventory in the range 20–60 kg allowed a 30% variation of the turbine inlet temperature with lower penalties on system performance than the turbomachinery speed control. The designed proportional-integral inventory controller showed a rapid response in the control of the turbine inlet temperature around the set point of 773 K during large variations of the heat load.


Introduction
Supercritical carbon dioxide (sCO 2 ) power cycles are an attractive energy conversion technology. One of their main advantages over conventional power generation technologies is the potential improvement in the overall power conversion efficiency due to the attractive properties of the working fluid in the supercritical state. Slightly above the critical point (73.9 bar and 304.25 K), CO 2 has a liquid-like density and gas-like viscosity. This enables a reduction in the compressor power consumption and the design of highly efficient and compact energy conversion systems [1]. Other beneficial properties of CO 2 are its nonflammability and high thermal stability, which make it also suitable for heat to power conversion from high temperature waste heat sources [2]. In the power generation sector, sCO 2 systems have been investigated for nuclear, Concentrated Solar Power (CSP), geothermal and Waste Heat Recovery (WHR) applications [3]. Even in fossil fuelled power plants, sCO 2 systems could offer advantages over conventional steam and gas turbine technologies by facilitating the use of carbon capture and storage systems [4,5] and the ability to adapt to and accommodate large variations of operating conditions imposed by the increasing penetration of renewables in the global energy mix [6].
Because of the wide variation in thermophysical properties and the high pressures in the cycle, the dynamic behaviour of sCO 2 power systems is still not fully understood. Additional research is also needed on the design of suitable control strategies that are paramount for system performance optimisation but also represent a challenge due to the lack of experimental data on compressor and turbine behaviours during transient operating regimes [7][8][9]. Furthermore, the safe operation of the power block must be ensured by avoiding undesired phenomena such as compressor stall/choking, turbine choking, extreme CO 2 pressures and temperatures, extreme shaft rotational speeds etc. [10].
To date, dynamic models of sCO 2 systems have been developed using tools such as ANL's Plant Dynamics Code [11][12][13], MIT's SCPS and TSCYCO [14] or commercial codes such as Modelica/Dymola® [15][16][17][18][19][20], MATLAB® [21][22][23][24], GT-SUITE TM [25,26] ASPEN PLUS® [27,28], Apros® [29]. Control systems have been mainly designed for the recompression and simple regenerative cycle configurations due to the cost and thermal efficiency advantages of these layouts compared to the various alternatives. The majority of the investigated controllers has been the Proportional-Integral type (PI) due to the ability to achieve a null steady-state error and the reduced sensitivity to high-frequency disturbances [10,30]. Some works also employed PID controllers that were set up and tuned through the Cohen-Coon technique [31]. Most of the published research has focussed on control strategies to guarantee the operational stability and maximisation of power output by acting on the compressor and turbine inlet conditions.
Wright et al. in [32] proposed to control the compressor inlet temperature by acting on the heat sink fluid. The combined control of the heat sink flow rate and the cooler by-pass resulted in a constant main

Symbols ξ
Pressure loss coefficient [-] ρ  compressor inlet temperature to avoid the two-phase region. Luu et al. in [33] suggested the control of the recompressor flow split ratio between the two compressors to maintain a sufficient compressor surge margin. The regulation of the recompressor flow control showed marginal benefits when used for load regulation compared to the inventory control. On the other hand, this approach resulted in optimal load tracking and with lower complexity of the overall control architecture [33]. A further control strategy consists in the variation of the working fluid mass in the CO 2 circuit, also known as charge, to fulfil given control requirements, e.g. a temperature or a pressure set point. Such strategy, commonly referred to as inventory control [27], relies on two storage tanks, one for feeding and one for bleeding CO 2 to and from the system respectively. The higher pressure tank supplies the additional CO 2 to the main loop while the lower one stores the working fluid extracted from the system. The regulation of the CO 2 mass (or charge) of the unit allows to change the working fluid mass flow rate and pressures in the circuit and, in turn, the achievable temperatures. Moyssetev et al. in [27] proposed this strategy to control the electric power generated, or electric load, from a 400 MWe nuclear power plant using supercritical CO 2 as working fluid. When the system operated between the 50% and the 90% of its nominal power, the inventory control maximised the cycle thermal efficiency. However, for fast variations of the electric load imposed by the grid, in an operating range above the 90% of the nominal point, the turbine flow by-pass resulted in a more suitable control strategy. At lower power ranges, between the 20% and the 50% of the nominal power, the turbine flow throttling proved to be a better control option despite causing a significant reduction in the system performance. Heifetz and Vilim in [34] developed a mixed control strategy to follow the generator load of a sCO 2 recompression cycle for nuclear applications using both the turbine by-pass valve and the system inventory as control variables. A simple model was used to validate the controller.
Despite these valuable works, the research to date has not reached an agreement on the most suitable control variable(s) for an optimal regulation of sCO 2 system performance. The inventory control seems a promising strategy, even though its setup poses technical challenges yet to be addressed. Among the major drawbacks and limiting factors of inventory control, the main ones are: stability implications (e.g. pressure gradients in the loop) due to the withdrawals/additions of CO 2 , slow response rate compared to throttling and by-pass controls, and the finite capacity of the storage tanks [35].
To fill the aforementioned research gaps, this work aims to investigate the potential of several operating parameters for the control of a 50 kW e sCO 2 system currently being commissioned at Brunel University London. The research methodology considers a transient onedimensional approach developed in GT-SUITE™, a commercial modelling platform. After the assessment of the most suitable control variables and the investigation of the system dynamics (during startup, shutdown and large variations in operating conditions), an inventory control system has been designed. The controller performance has been investigated by simulating different transient heat load profiles.
Novel aspects of the research include: the development of submodels calibrated from actual heat exchanger, turbomachinery and generator equipment; a heat recovery application with direct sCO 2 heating from the industrial flue gas; testing of the designed control strategy on a complex system accounting for the geometric features and materials of the heat exchangers, turbomachines and pipes; simulations with reference to transient waste heat load profiles typical of energy intensive industries. The results of the study are expected to outline and evaluate procedures to start, operate, control and shutdown of experimental rigs and form the basis for the development of controls for fullscale systems.

System description
The modelled sCO 2 system refers to a plug and play 50 kW e sCO 2 unit for Waste Heat Recovery (WHR) applications available at Brunel University London. The sCO 2 facility in Fig. 1 is based on a simple regenerative Joule-Brayton cycle (Fig. 1a) and employs three heat exchangers as follows: a micro-tube gas/sCO 2 heat exchanger, which performs the actual and direct heat recovery from the heat source; a Printed Circuit Heat Exchanger (PCHE), used as recuperator; and a Plate Heat Exchanger (PHE) which is employed as gas cooler. PCHEs have been widely used in the oil and gas sector because of their reliability to operate with working fluids at high pressures and temperatures [36], and for this reason have been also adopted in sCO 2 power applications [37]. Even though PCHEs have been used also as gas coolers in other sCO 2 test rigs [38], in this test facility a more conventional technology, as PHEs, has been selected as gas coolers in order to substantially reduce costs [39], which is of paramount importance in WHR and any industrial application. Indeed, in the gas cooler, the working fluid achieves the lower pressure levels in the cycle (between 75 bar and 85 bar), which today are common in current CO 2 refrigeration applications, where PHEs are successfully employed. An overview of the facility is shown in Fig. 1b. A gas-fired Process Air Heater simulates the industrial waste heat source and provides flue gases up to 1.0 kg/s and 1073.15 K. The heat sink of the facility is a 500 kW dry-cooler system. The dry-cooler rejects the heat from a water/glycol mixture which, in turn, removes the residual heat from the CO 2 stream in the gas cooler. Variable speed drives for the water pump and fans allow the variation of the heat removal and rejection rate. A series of electric heaters can be eventually used to warm up the water and then indirectly heat up the CO 2 to achieve supercritical conditions at the system startup. This scenario is however not simulated in this research.
The sCO 2 packaged unit is located in a standard 20ft container (Fig. 1c) and gathers all the system components except for the microtube heat exchanger, also known as primary heater, which is located along the exhaust line of the Process Air Heater (Fig. 1b). The primary heater has been specifically designed to enhance the gas/sCO 2 heat transfer without causing an excessive pressure drop on the flue gas side, and it is connected to the sCO 2 loop in the container through two 10meter pipes.
On the right-hand side of Fig. 1c it can also be seen the Compressor-Generator-Turbine (CGT) assembly. The assembly considers two radial machines, a compressor and a turbine, mounted on the same shaft of a synchronous brushless generator. A complex ancillary system is used to ensure lubrication of the bearings and to reduce windage losses. Two motorised compressor and turbine by-pass globe valves (CBV and TBV respectively), shown in Fig. 1a and c, have been foreseen to control the system at nominal, startup, shutdown and emergency operations. The CBV valve is located between the compressor outlet and gas cooler inlet while TBV is placed between the turbine inlet and outlet (Fig. 1c). A data acquisition system has been also installed in the container (left-hand side of Fig. 1c) to allow the remote control and monitoring of the unit.
The sCO 2 test facility has been designed for a 50 kWe power output, which is a very small scale for sCO 2 technology. In fact, at these operating conditions, the turbomachinery impellers become very small and, in turn, their efficiency is likely to be affected by severe losses. Another challenge for small-scale sCO 2 systems is the unitary cost per kilowatt installed. In fact, the capital expenditure is not linear with the power size especially for what concerns the heat exchangers, which typically absorb most of the equipment budget [40]. Therefore, unless costs are significantly reduced, sCO 2 technology is primarily competitive for applications beyond 1 MWe power output and heat source temperatures above 620 K (~350 • C). Despite the power mismatching, this research at pilot scale is expected to advance the knowledge in the field of sCO 2 power research and further support the market uptake of sCO 2 technology. In fact, the cycle layout of a full-scale system will stay as close as possible to a simple regenerative layout to achieve commercial appeal. Furthermore, as reported in [41] sCO 2 turbomachinery will still be radial as long as the system power output stays below 10 MWe. For these reasons, the analysis that will be reported in the following sections are likely to be transposable to develop and operate full-scale sCO 2 power systems.
Further details on the Brunel's high-temperature heat to power conversion facility are available in [42].

Modelling methodology
The model of the sCO 2 heat to power conversion system has been developed in the commercial software GT-SUITE™. This tool is based on a one-dimensional (1D) formulation of Navier-Stokes equations and on a staggered grid spatial discretization [43]. Each component can be independently modelled through accurate input data, whose availability is crucial for the overall model development. The required inputs relate to the geometrical features of the components as well as their performance data. This information can result either from an experimental campaign or, as in the current case, from more detailed or complex models (e.g. 3D RANS CFD).
The components modelled as equivalent 1D objects are heat exchangers and pipes while the turbomachines, valves and the receiver are treated with a lumped approach. The 1D models discretise the components into a series of capacities such that manifolds are represented by single volumes while pipes are divided into one or more volumes. These volumes are connected by boundaries. The scalar variables (pressure, temperature, density, internal energy, enthalpy, etc.) are assumed to be uniform in each volume. On the other hand, vector variables (mass flux, velocity, mass fraction fluxes, etc.) are calculated for each boundary [43].
Each capacity considers the algebraic sum of all the incoming and outgoing mass flow rate contributions from the neighbouring capacities that occur at the boundaries, as per the continuity equation (1).
The pressure dynamics in the system is calculated through the momentum equation (2), which neglects body forces and takes into account the algebraic sum of momentums, pressure forces and dissipations due to friction and pressure drops through the boundaries [43]. In pipes, the latter two terms are respectively related to distributed (i.e. due to surface roughness) or concentrated (i.e. due to bends) pressure losses.
The energy equation (3) is expressed in terms of total enthalpy. This formulation is required for the further implicit integration scheme employed by the solver for the analysis of energy systems whereas resolving fast dynamics (e.g. indicating pressure in positive displacement machines) is not the end goal [43]. Neglecting variations of potential energy, for a given capacity, the rate of change of total enthalpy depends on the volume capacity variations, the enthalpy fluxes and the heat transfer phenomena. The solution of the energy equation requires the computation of the local heat transfer coefficient through calibrated heat transfer correlations.

Turbomachines
Compressor and turbine have been modelled through performance maps since it is assumed that turbomachines have faster dynamics compared to heat exchangers. Therefore, performance maps have no intrinsic dynamics and merely act as look-up tables to provide the boundary conditions to the equipment upstream and downstream the machines. Table 1 summarizes some geometrical details about their aerodynamic design. Additional information on the turbomachinery design is available in [45].
The turbomachinery performance maps result from 3D RANS CFD simulations whose modelling methodology has been discussed in [25,46]. The 3D modelling approach has been validated through experimental data available from the Sandia National Laboratories [32]. The inlet boundary conditions of the 3D model are the total pressure and temperature as well as the flow direction, which is considered normal to the boundary. Outlet average static pressure has been chosen as outlet boundary condition. Fig. 3a and b show the turbine maps in a typical format typically adopted in turbomachinery [47]: the operating curves and the total-tostatic isentropic efficiency of the machines are reported as a function of the pressure ratio, the reduced mass flow rate and the reduced revolution speed. For the pressure and temperature ranges in which the turbine operates, it is reasonable to assume an ideal gas behaviour for the CO 2 and a unique reference state for the normalisation, namely 650 K and 145 bar. As shown by the turbine operating curves reported in Fig. 3a, the maximum mass flow rate that the machine can deliver before choking occurs at 3.4 kg/s. From the isentropic efficiency map shown in Fig. 3b instead, it is possible to notice that the optimal operating range for the designed turbine occurs at high revolution speeds (over 80,000 RPM) and for a pressure ratio between 1.7 and 2.3.
Assuming a single reference state is not valid for the compressor, which instead operates close to the critical point. In this region, the real gas properties of the fluid must be considered and the use of reduced quantities via normalization can lead to errors in the compressor performance predictions [48]. To overcome this issue, the model considers multiple compressor maps at four reference states that span the whole CO 2 critical region (308.15 K at 70, 75, 80 and 85 bar). Each of these maps has been generated by maintaining constant the inlet conditions of the working fluid (pressure and temperature) and changing the outlet static pressure at different revolution speeds. At least five pressure ratios for each revolution speed are required to ensure an accurate interpolation of the operating curves. Beyond the speed range available from the calibration data, a linear extrapolation method is used to predict the performance of the turbomachines. Fig. 4a and b show the compressor performance maps, operating curves and total-to-static isentropic efficiency for the reference state at 75 bar and 308.15 K. The small distortion located near the surge line at high revolution speeds is due to the supercritical CO 2 thermophysical properties, which affect the pressure changes inside the compressor at different speeds.
The design point for the compressor used in the current study approaches the critical point. This thermodynamic condition has a positive effect on the choke line. In fact, its slope suggests that an operation of the compressor close to the critical point would be beneficial in terms of choke margin compared to the one typical of conventional machines using ideal gases. This is in line with the findings in [49]. In particular, for a pressure ratio of 2, the compressor can process a minimum mass flow rate of 2.1 kg/s before choking. For lower pressure ratios, for instance 1.3, the minimum mass flow rate of working fluid processable is 1.1 kg/s, which shows the advantage of using sCO 2 as working fluid in terms of operational flexibility (Fig. 4a). The compressor optimal operating region is located close to the choke line, as shown in Fig. 4b, where the isentropic efficiency assumes a value of 0.8 in a range of revolution speeds and pressure ratios between 50000-90000 RPM and

Heat exchangers
The primary heat exchanger employs a micro-tube technology. Similarly to a shell and tube heat exchanger, the high pressure sCO 2 flows in the micro-tubes while the flue gas path is in cross flow along the shell side. As shown in the 1D schematic of Fig. 5a, the heat exchanger is divided in four modules composed of a series of micro-tubes with an outer diameter of 2 mm. Fig. 5b also reports a cross-sectional view of the modules with the typical micro-tube arrangement. The tube spacing in the transverse (X p ) and in the longitudinal (X n ) directions are specifically designed to enhance heat transfer without an excessive pressure drop increase on the flue gas side. Table 2 summarizes the detailed geometrical features of the micro-tube heat exchanger. Fig. 5c instead displays the heat exchanger model block diagram developed in GT-SUITE™. Each module is represented by a predefined template requiring the following geometrical features as inputs: tube diameter, thickness, wetted perimeter, cross sectional area and number of parallel tubes. Because of the 1D modelling approach adopted, the tube staggered configuration cannot be reproduced in the model. Therefore, the use of calibration coefficients was considered.
The recuperator of the sCO 2 system is a Printed Circuit Heat Exchanger (PCHE) made of Stainless Steel 316L (SS316L) plates chemically etched with semi-circular channels which are bonded together through a thermal diffusion process. The hot and cold side channels are stacked alternately on top of each other, and nozzles distribute the fluid stream into several channels. These channels are equally spaced and with identical geometrical features. As such, the model of this heat exchanger can be substantially simplified by considering a series of elementary heat transfer units comprised of a pair of channels surrounded by solid substrates and periodic boundary conditions. A more detailed explanation of the modelling procedure is available in [50]. Table 3 summarises the geometrical features of the PCHE recuperator.
The gas cooler is a Plate Heat Exchanger (PHE). This technology was selected given the less harsh operating conditions (lower pressure and temperature) and the availability of PHE from the CO 2 refrigeration sector. Similarly to the PCHE, the PHE modelling has been carried out discretising the heat exchanger with a sequence of elementary heat transfer units. The coolant fluid has been assumed as pure water. The complete modelling procedure is reported in [30]. Table 4 summarises the geometrical features and material of the PHE gas cooler.
In each heat exchanger, the geometrical inputs are used to define the properties of the equivalent 1D channels. Data provided by the manufacturers or from more complex models (e.g. 3D CFD) are then employed to calculate the best fitting coefficients of the Nusselt-Reynolds (Nu-Re) correlations for the equivalent 1D networks that the heat exchangers are approximated with [43]. The performance data provided by the manufacturer or obtained from more complex models are related to several operating conditions of the heat exchangers.
The pressure drops are computed using a modified version of the Colebrook relationship in Equation (4). In this expression, the Fanning factor is calculated using the explicit approximation of the Colebrook equation proposed by Serghides [51], which is valid for the turbulent regime (Re D greater than 2100). The quantities C 2 and C 3 , which can be calculated using Equations (5) and (6), account for the roughness of the heat exchanger channels Ra. The term C 1 is the calibration coefficient used to adapt the simulation results to the performance data provided by the heat exchanger manufacturer.
Even though this modelling methodology is common to all the three heat exchangers considered, the gas cooler requires an additional correlation to account for possible condensation of CO 2 . In this case, to predict the phase change, the formation of vapor bubbles or liquid droplets is addressed by evaluating the fluid density in each sub volume, while the the two-phase area is computed using the vapour Rayleigh-Plesset formulation in Equation (7) [52].
The results of the regression analyses carried out to calibrate the several heat exchangers are detailed in Table 5, which compares the Re-Nu curve interpolation of the different data provided by the manufacturer against the ones obtained by using the Gnieliski [53] and Dittus-Boelter [54] heat transfer correlations. It can be observed that the Gnieliski correlation provides better predictions of the manufacturer data. For this reason, in the absence of data on heat exchanger performance from experimental tests, the Gnieliski correlation for the calculation of the heat transfer coefficients has been employed in this study [50]. Table 6 lists the number of sub-volumes in which the different heat exchangers have been discretized and their time constant, calculated as the ratio between the heat exchanger mass multiplied by the specific heat capacity (m*cp) of the material and conductance (UA) of the heat exchanger. It can be seen that the primary heater has the lowest thermal inertia and the recuperator the highest due to its much higher thermal mass of the material used for its manufacture.

Valves and other equipment
The valves have been modelled as orifices. The valve template is in the form of a look-up table that provides the flow coefficient as a function of the pressure differential across the device. The globe valves have been designed by the manufacturer to follow an equal percentage characteristic curve. A series of forward and reverse discharge coefficients have been used as a function of the lift position of the valve actuator. These discharge coefficients are used by the software to calculate the effective flow area at the throat, while the pressure ratio across the valves is used to compute the velocity at the throat and, consequently, the mass flow rate through the valve. Equation (8) correlates the valve discharge coefficient to ratio between the lift L and the valve diameter D [55].
Straight pipes and bends are simulated using the 1D modelling approach described in Section 3. In the current study, pipes are considered smooth and insulated, thus friction and thermal losses are neglected. This is assumption can be justifies by the fact that the pipe diameters used in the Brunel sCO 2 test facility have an oversized diameter (2 ′′ ) to reduce pressure losses and are insulated with a ceramic wool layer wrapped between an inner layer of silica wash treated glass cloth and an outer layer of grey PTFE coated glass cloth to minimise thermal losses.
The receiver has been modelled as a container (capacity) with fixed volume. The receiver is situated downstream the gas cooler to absorb the thermal expansion of the fluid in the circuit, decoupling the high from the low-pressure side of the system (Fig. 2). Its volume is 0.165 m 3 and accounts for almost 50% of the overall system capacity.
The inertia of the generator shaft has been taken into consideration in the simulation but the electrical power generation process has not been considered in detail in this study. Parasitic losses of the system ancillaries (water cooling pump and fans, oil pump for turbomachinery lubrication, CO 2 drainage compressors) have been also discarded.

Results and discussion
After the modelling stage, a series of simulations have been carried out to assess the steady-state and transient behaviour of the sCO 2 system. The steady-state analysis has been focused on the assessment of the most suitable control variables for the regulation of the system under different waste heat loads. In particular, the revolution speed of the turbomachines and CO 2 mass charge in the circuit, namely system inventory, have been considered. The mass charge (or inventory) in the sCO 2 loop has been varied from 20 kg to 60 kg while the revolution speed of the turbomachines ranged between 50,000 RPM and 90,000 RPM. During these simulations, the inlet conditions of the heat sink and the heat source (e.g. mass flow rates and temperature) have been kept constant and equal to their nominal values. These quantities are reported in Table 7.
After the steady-state analysis, to further investigate the dynamic response of the sCO 2 system, two transient heat load profiles have been simulated. In particular, the first profile has been obtained by varying the flue gas mass flow rate and maintaining the flue gas inlet temperature equal to the design condition. For the second one, the mass flow rate of the flue gas has been kept constant while the inlet temperature has been varied.

Steady state analysis
It is worth to start by clarifying that the results presented strongly depend on the nature of the turbomachine maps. As such, the same sCO 2 system equipped with different turbomachines might exhibit a different behaviour. Fig. 6a shows the effects of the variations of the    (Fig. 6a). Hence, the TIT is more sensitive to variations of the system inventory rather than the turbomachinery revolution speed. This difference can be explained by referring to Fig. 6b. Increasing the system mass charge at higher revolution speeds leads to higher CO 2 mass flow rates circulating in the system, which allow to lower the temperature of the working fluid at the inlet of the turbine. Operating the system at these conditions is particularly important in order to reduce the turbine inlet temperature in case of high peaks of the waste heat source temperatures. Using the revolution speed as control variable would in fact not allow to keep the TIT in the safe temperature margin. With regards to the CO 2 pressures at the turbine inlet, the variation of both the turbomachinery revolution speed and the CO 2 inventory in their whole respective range does not lead to excessive pressures. For an inventory of 60 kg/s and a revolution speed of 90,000 RPM, the maximum pressure achieved at the turbine inlet is equal to 150 bar (Fig. 7a), which can be considered in the safe operating margin of the system. Fig. 7 shows also the effect of the CO 2 injection in the system on the compressor inlet pressure, which goes from 50 bar to more than 80 bar when the mass charge varies from 20 kg to 60 kg (Fig. 7b).
It is possible to notice that the compressor inlet pressure is only affected by the supercharging of the circuit. On the contrary, the turbine inlet pressure depends on both control variables since the pressure at the inlet of the machine is affected by the lowest cycle pressure (set at the compressor inlet) and also from the cycle pressure ratio, which depends mainly from the turbomachinery revolution speed (Fig. 7c). For this reason, the turbine inlet pressure increases following the rise of both the system inventory and the turbomachinery revolution speed (upper part of Fig. 7a).
The cycle pressure ratio also affects the turbine inlet temperature as can be noticed from comparison of Figs. 6a and 7c. Increased expansion ratios across the turbine (top right-hand side of Fig. 7c) lead to lower temperatures at the turbine outlet and, in turn, to reduced residual heat available at the hot side of the recuperator. This lower amount of residual thermal energy available, leads to lower temperatures at the outlet of the recuperator cold side and, assuming a constant heat supply from the flue gases, at lower temperatures at the turbine inlet.
A further interesting aspect is the effect of turbomachinery revolution speed and the system inventory on the turbine performance and on the system net power output. Fig. 8a shows that an inventory variation in the considered range, for a given fixed revolution speed of 80,000 RPM, does not affect substantially the turbine isentropic efficiency of 73%. On the contrary, the change in the revolution speed has a more pronounced effect on the turbine isentropic efficiency, which drops from the optimal value of 74% at 86,000 RPM to 65% when the speed is decreased to 50,000 RPM (for a system mass charge of 40 kg as shown in Fig. 8a). This trend can be explained by looking at the turbine maps shown in Fig. 3. The map shows that, for a wide operating range of the machine, the isentropic efficiency remains almost constant, while it drops consistently for small cycle pressure ratios and lower CO 2 mass flow rates. These operating conditions are mainly achieved at lower revolution speeds of the compressor-turbine-generator unit.
The losses in the turbine performance directly affect the system net power output, which decreases for the same variation of the revolution speed from 50 kW to 20 kW (Fig. 8b). At the nominal revolution speed of the turbomachines, the net power output of the system strongly depends also by the system inventory, going from 10 kW when the system is filled with 20 kg of CO 2 to 50 kW for a CO 2 mass in the circuit of 40 kg (Fig. 8b). When the system is in fact undercharged, the pressure at the compressor inlet is lower than 74 bar and thus the supercritical conditions of the CO 2 are not achieved, compromising the performance of the unit.

Transient heat load variations
Even for small power outputs, heat to power conversion technologies such as Steam or Organic Rankine cycles have relatively large components, with high thermal mass and working fluid inventories. This leads to a high thermal inertia which influences the system dynamics [30]. On the other hand, for the particular sCO 2 system under consideration, the inertial effects are less prominent and a higher operational flexibility can be obtained due to the high density of the working fluid and the small size of components. Fig. 9 demonstrates this concept, showing the response of the CO 2 temperature at the inlet of the turbine (continuous line) and the compressor (dashed line) to a variation of the heat load at the primary heat exchanger resulting from a change in the flue gases mass flow rate. In particular, the time varying mass flow rate profile consists of a number of ramps and plateaus. The time length of the plateaus is of 400 s while the ramp lasts 200 s. During the plateaus, the mass flow rate of the flue gases assumes the same minimum and No delays are noticeable in Fig. 9 with regards to the dynamic response of the turbine inlet temperature. Hence, for this specific system, the thermal inertia effects are negligible. In particular, the turbine inlet temperature (red line) promptly varies from the value of 743.15 K to 798.15 K, when the gas flow rate goes from 1.00 kg/s to 1.25 kg/s and drops to 683.15 K when the gas flow rate decreases to 0.75 kg/s. The same trend holds for the pressure, which goes from a minimum of 135 bar at 0.75 kg/s of flue gas up to a maximum of 140 bar when the flue gas mass flow rate increases to 1.25 kg/s.
The changes of the working fluid thermodynamic conditions at the compressor inlet change less notably following the heat source mass flow rate variations, with the compressor inlet pressure remaining almost constant and equal to 72.5 bar. This effect is mainly due to the receiver placed downstream the gas cooler, whose large capacity contributes to decouple the hot and the cold sides of the sCO 2 loop.
The inverse response of the compressor inlet temperature, which goes from 312.15 K to 313.15 K (Fig. 9) when the mass flow rate of the flue gases decreases from 1.25 kg/s to 0.75 kg/s, can be explained by the increase of the CO 2 mass flow rate circulating in the system. A higher heat load leads to a decrease of the working fluid density, due to the increase of its temperature, and consequently to a reduced CO 2 mass flow rate (Fig. 9). However, the cooling load supplied by the heat sink remains fixed, and therefore a lower temperature at the compressor inlet is achieved. The opposite occurs for a decrease of the hot source mass flow rate (Fig. 9). Fig. 10 shows the response of the system to a variation of the flue gas inlet temperature (black line) from 723.15 K to 1123.15 K with same time duration of ramps and plateaus imposed for Fig. 9. The trends of the thermodynamic conditions (pressure and temperature) of the working fluid at the inlet of the compressor and the turbine resemble the ones shown in Fig. 9, but more accentuated given the deeper effect of the hot source inlet temperature on the cycle parameters. In particular, the turbine inlet temperature goes from 623.15 K for the minimum temperature value of the flue gases (723.15 K), up to 923.15 K for its maximum one (1123.15 K), while the compressor one goes from 312.15 K to 314.15 K for the same temperature variation of the flue gases (Fig. 10).

Startup and shutdown
After a pre-heating phase in which the working fluid is brought above the supercritical state to prevent any liquid pocket in the loop, the turbomachines are span to their nominal revolution speed. To investigate the system dynamics during these conditions, a rotational speed profile has been imposed to the single-shaft turbomachinery unit as boundary condition as shown in Fig. 11. In particular, the turbine and compressor speed increases from 0 to 86,000 RPM in 900 s during the startup phase, it is then held to the nominal value of 86,000 RPM for 1200 s and then is decreased from 86,000 RPM to 0 RPM in 700 s. The ramp rates have been set according to the maximum acceleration allowed for the generator rotor, which is equal to an increase of 100 RPM every second for the startup stage and to a 130 RPM decrease every second for the shutdown (black dashed-line in Fig. 11). During the whole transient simulation, the heat source mass flow rate and inlet temperature have been kept constant and equal to their nominal values of 1 kg/s and 923.15 K respectively.
When the turbine and compressor bypass valves (Fig. 2) are both closed, for low values of the turbomachinery speed, the temperature at the turbine inlet would be equal to 923.15 K due to the low mass flow rate in the system and the low pressures achieved in the circuit (cyan line, Fig. 11). To overcome this issue, two different strategies for the turbine by-pass globe valve (TBV) and the compressor by-pass globe valve (CBV) have been simulated. In the first case, both the TBV and the CBV are fully open during the acceleration and deceleration decrease ramps, and closed when the nominal value of 86,000 RPM is achieved. Fig. 7. Contour maps of the turbine (a) and compressor (b) inlet pressure, and the cycle pressure ratio (c) as a function of the turbomachinery revolution speed and the system inventory.   In the second strategy, the CBV is kept fully closed and only the TBV is operated. Similarly to the first strategy, the TBV is opened during the startup and shutdown stages (when the revolution speed of the turbomachines is increased and decreased respectively) and closed when the nominal value of 86,000 RPM is achieved. In both strategies, the valve opening profiles are the same: the second strategy differs from the first one only because the CBV is not operated.
The standalone opening of the TBV (magenta line, Fig. 11) achieves better results in terms of turbine inlet temperature reduction, since it allows the temperature to decrease after a limited time period of 225 s, bringing it down from 923.15 K to 753.15 K when the speed increases from 25,000 RPM to 86,000 RPM. A much less pronounced effect is produced by the simultaneous opening of both valves (grey line, Fig. 11). In this case, the turbine inlet temperature is reduced only from 923.15 K to 903.15 K when the turbomachinery speed exceeds 55,000 RPM. Therefore, the beneficial effect due to the opening of the turbine by-pass valve is counter-balanced by the opening of the compressor bypass, which decreases the mass flow rate going through the primary heater and therefore leads to a temperature increase at the turbine inlet.
To better understand how the turbine by-pass valve opening affects the dynamics of the working fluid temperature at the turbine inlet, reference can be made to Fig. 12. This chart displays the startup and shutdown turbine mass flow rate and pressure ratio evolution for the three strategies presented in Fig. 11 superimposed to the turbine operating map. When TBV and CBV are both closed (cyan circles) the working fluid flows entirely through the turbine and this prevents the achievement of the suitable pressure ratio across the machine during the revolution speed increasing and decreasing ramps (as shown in the low part of the left-hand side of Fig. 12, where the pressure ratio starts to increase only for revolution speeds higher than 50,000 RPM during the unit startup and shutdown). As a consequence, a lower CO 2 mass flow rate circulates in the loop and higher temperatures are achieved at the heater outlet despite the same heat load. When TBV and CBV are fully opened (grey circles), the mass flow rate flowing into the turbine decreases. During the startup and shutdown, this leads to higher pressure ratios than the former case (lower left-hand side part of Fig. 12). However, since the compressor by-pass valve is also open, the compressor is not able to fully pressurise the fluid due to the lower mass flow rate. On the contrary, the sole opening of the turbine by-pass valve (magenta circles), allows to achieve a more gradual build-up of the pressure ratio across the machine, which leads to higher mass flow rates circulating throughout the loop and a faster achievement of the nominal working fluid temperature at the turbine inlet (Fig. 12). From this analysis it can be concluded that during the startup and shutdown phase of the unit, reducing the mass flow rate throughout the turbine and normally operating the compressor, leads to a faster and more gradual increase of the working fluid mass flow rate circulating in the loop and thus the achievement of the nominal operating conditions of the system limiting the exposure of the components to extreme temperatures. Fig. 12 further allows to explain the small fluctuations of the turbine inlet temperature depending on the strategy considered in Fig. 11. This is mainly due to the different build-up of pressure in the machine for each different case. Indeed, it is possible to observe that, during the startup and the shutdown of the turbomachinery, the pressure ratio across the turbine changes with time following a different path for each of the considered strategies. These pressure ratio variations lead to small fluctuations of the working fluid mass flow rate processed by the turbine and, in turn, in the whole circuit. Since the heat load is assumed as constant in these simulations and the energy balance must hold at the primary heater, these variations of working fluid mass flow rate result in temperature fluctuations at the turbine inlet.

Control system design
The results from the steady-state analysis showed that, for the system considered, the regulation of the inventory represents a more suitable option to control the turbomachinery inlet conditions without excessively compromising the system performance and operational stability. In particular, the control strategy here proposed acts on the CO 2 mass in the circuit to regulate the turbine inlet temperature, which is an important parameter both from a cycle performance perspective as well as to prevent unsafe operating conditions due to unforeseen peaks in the waste heat load.
The inventory control envisaged works as follows: for a given change of the heat load, if the turbine inlet temperature becomes higher than the set point value, the system is charged with additional mass of working fluid, which flows from the high pressure inventory control tank to the circuit. The increased CO 2 mass in the system (or system inventory) leads to a rise of the minimum cycle pressure at the compressor inlet and the CO 2 mass flow rate processed until the control error is eliminated. On the other hand, if the temperature at the turbine inlet drops, the working fluid is withdrawn from the main loop and flows into the lowpressure inventory control tank, lowering the minimum cycle pressure  and mass flow rate.
The addition and the withdrawal of the CO 2 is thus possible thanks to two inventory control tanks, one at higher and one at lower pressure compared to the circuit charging point. The tank used to inject the CO 2 in the system is at higher pressure than the circuit charging location, while the tank used to store the extracted CO 2 from the system is at a lower one. By modulating these pressure differentials between the high/ low pressure tanks and the circuit charging point (which can be done by using two pressure regulators located downstream the tanks) the control system is able to set the right amount of CO 2 to be injected/withdrawn into/from the circuit.
The hardware implementation of this control approach would require, in addition to the two tanks at high and low pressure connected to the same charging point, a three-way valve to exclude the addition or withdrawal line depending on the desired control action. Two pressure regulators should be placed upstream the tanks to be able to change the pressure level between the two extreme values of the high-and lowpressure tanks. From a modelling perspective, this hardware architecture can be simplified by considering just one inventory tank able to assume a certain level of pressure between the maximum and minimum values of 110 bar and 55 bar, as shown in Fig. 2. The tank is also assumed to have an infinite capacity while the dynamics of the pressure regulators and of the three-way valve have been neglected in the current study.
The two inventory tanks are connected to a circuit feeding point between the gas cooler and the receiver, as shown in Fig. 2. This design choice allows to avoid detrimental effects on the compressor operation due to a possible fast charge or discharge of the system. In fact, the large volume of the receiver smooths the sudden pressure variations following the increase or decrease of the working fluid mass charge in the system. The working fluid temperature injected in the main loop is assumed to be at 308.15 K.
The modulation of the inventory control tank pressure is commanded by a Proportional-Integral (PI) controller to reduce the error between the actual turbine inlet temperature and the set point. The derivative term has been set to zero because of the possible high control actions resulting from noisy signals. To prevent the saturation of the control action caused by the integral term, an anti-wind up scheme has also been adopted. The proportional and integral coefficients have been calculated via a lambda tuning procedure [56]. In particular, a first order relationship has been considered between the mass flow rate of the fluid injected in the system and the controlled process variable (turbine inlet temperature). Therefore, the control output has been modified in the entire admissible range and the process variable response analysed [43]. From the calculated time constant (τ) and process gain (K), the resulting proportional and integral coefficients of the PI controller are 0.177 and 0.011 respectively.
To firstly assess its performance, the turbine inlet temperature is regulated assuming the transient heat load profile reported in Fig. 10, where the turbine inlet temperature of the waste heat source is varied between the two plateaus of 623.15 K and 1123.15 K while the mass flow rate has been kept constant and equal to the nominal value of 1 kg/ s. Both the uncontrolled and controlled response of the turbine inlet temperature system have been reported in Fig. 13. The temperature set point has been assumed equal to 773.15 K and it is displayed as a dashed grey line. It is possible to notice that the controller is able to regulate the   (Fig. 13). From a performance perspective, it is possible to notice that the control of the CO 2 inventory leads to an improvement of the system thermal efficiency (which is calculated with respect to thermal energy recovered by the waste stream) and power output when the temperature of the flue gas increases. In particular, for the controlled case an increase of the 25% and of the 101% in the system thermal efficiency and net power output respectively can be observed, compared to the uncontrolled case, when the flue gas temperature achieves the maximum value of 1123.15 K. On the contrary, when the flue gas temperature decreases to 723.15 K, the action of the inventory controller is detrimental, leading to a system operating condition where no electric power can be extracted by the heat recovered (net power output of − 10 kW and system thermal efficiency equal to − 0.05). These effects are due to working fluid mass flow rate injection and extraction imposed by the controller. Indeed, following an increase of the flue gas temperature, the loop is charged with additional mass of working fluid, while the system is discharged in the opposite case. Therefore, when the temperature of the flue gases is too low, to achieve the desired temperature set point, the inventory control system imposes an excessive extraction of CO 2 mass from the loop, leading to an operating condition where the turbine cannot supply the required mechanical power required from the compression to pressurise the working fluid. Indeed, for a lower mass flow rate circulating in the circuit, the isentropic efficiency of the compressor decreases, causing the detrimental effect on the system performance. Table 8 shows the different coordinates of the turbine and compressor efficiency maps when the waste heat source assumes the minimum temperature value of 723.15 K for both the uncontrolled and controlled system response.
Despite then the inventory control has been considered as one of the most suitable control strategies for sCO2 systems in the power generation sector [57], where the system thermal efficiency is optimised depending on the electric load set by the grid, this does not apply when waste heat recovery applications are considered. When indeed the objective is to maximise the power recovery depending on the dynamics of heat load supplied by the topping process, namely a heat load following strategy, the use of an inventory control system presents some criticalities since introduces a strong detrimental effect on the system performance for reduced thermal levels of the waste heat source. This limitation could be overcome by the adoption of multi-level control system combining multiple control actions depending on the heat load available.
To test the performance of the controller on a real industrial case, a second transient heat load profile has been simulated. The profile refers to a typical process in an energy intensive industry. Fig. 14 , Fig. 14).
The controlled response shows that the designed controller not only eliminates the error between the actual turbine inlet temperature and the set point, but also shows good dynamic performance. The slight overshoots at 300 s, 860 s and 1240 s due to the control action saturation are indeed limited.
It is possible also to notice that the turbine inlet temperature in the uncontrolled system is not able to achieve the set point value of 773. 15 (Fig. 14). From a performance perspective, for a such waste heat load profile, the inventory control leads to a net improvement of the system power output and thermal efficiency since no substantial decreases of the waste heat source inlet temperature occur.

Conclusions
This research investigated the dynamic behaviour of bottoming supercritical CO 2 (sCO 2 ) heat to power systems and discussed control approaches with reference to typical operating conditions of sCO 2 power blocks, such as startup, shutdown and transient heat load profiles from the topping industrial process. The sCO 2 high temperature heat to power conversion system at Brunel University London, based on the simple regenerative Joule-Brayton cycle with a nominal electrical power output of 50 kW has been considered as test case. The research methodology employed a one-dimensional (1D) modelling approach to perform offdesign, transient and control studies.
A steady-state analysis has been performed to assess the most suitable control variables for the regulation of the inlet conditions at the turbine. The system inventory showed to be a good parameter for the control of the turbine inlet temperature having only a small impact on turbine performance. A reduction in CO 2 from 60kg to 20kg would indeed result to only 1.0% reduction in the turbine isentropic efficiency.
Transient variations of the flue gas flow rate and temperature in a time interval of 2800 s always led to changes in the turbine inlet temperature. The turbine inlet pressure was slightly affected, especially at variable flow conditions but constant temperature of the heat source (max 3% of design point for a variation of the heat source mass flow rate of the 25%). More importantly, unlike ORC systems of similar power output, the turbine inlet temperature variation did not experience any lag with respect to the topping process. This demonstrates the flexible nature of the sCO 2 system investigated. Furthermore, during the startup and shutodown of the unit, the opening of the turbine by-pass globe valve could guarantee a safer operation of the turbine and a quicker build-up of the pressures and temperatures in the sCO 2 loop.
Based on these outcomes, an inventory control strategy has been designed to regulate the turbine inlet temperature during transient heat load profiles. The PI controller showed good accuracy and dynamic performance. The single envisaged location for additions/withdrawals of CO 2 downstream of the gas cooler prevented instability in the compressor operation for rapid changes of the waste heat source temperature.
The analysis further highlighted the possible detrimental effect on system performance introduced by the inventory regulation should a heat load following control strategy be considered. A sudden temperature drop at the waste heat source would action excessive withdrawals of working fluid mass flow rate from the main sCO 2 loop and, in turn, lead to negative power outputs being the compression power higher than the turbine one.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
This section summarises the design specifics for the High Temperature Heat To power Conversion (HT2C) at Brunel University London [42]. Fig. A1 shows the cycle layout and heat exchangers diagrams while Table A1 lists the design conditions for the working fluid (CO 2 ) as well as heat source (flue gas from natural gas combustion) and heat sink (75% water, 25% glycol).

Fig. A1
. Design specifics of the Brunel's HT2C facility: sCO 2 system layout (a) and Temperature vs. heat load diagrams of primary micro-tube heater (b), printed circuit recuperator (c) and plate gas cooler heat exchangers (d).

Table A1
Process variables at design conditions with reference to Fig. A1