Transient numerical analysis of charging a heat accumulator in a combined heat and power plant

This paper presents a technical description and characteristics of heat accumulation system operation in one of the Polish combined heat and power (CHP) plants. The total thermal capacity of the nonpressure hot water accumulator is 75 MWh. An important issue in the case of hot water tanks is thermal stratification caused by the difference in density between higher-temperature water, located in the top of the hot water storage tank and colder water in its bottom. The article presents the results of numerical analysis based on axial-symmetric model of hot water tank using Ansys Fluent software. Simulations with the Computational Fluid Dynamics (CFD) method were carried out for the process of charging the heat accumulator lasting 10 hours. The results were compared with the measurement data from the CHP plant.


Introduction
In recent years, interest in using hot water tanks as energy storage has increased significantly. This is mainly due to the changes taking place on the electricity market and in particular the increase in production from renewable energy sources (RES). Research carried out all over the world mainly concerns issues related to storage and management of energy in the form of heat for the needs of single-family, multi-family buildings and large hot water tanks used in power plants, heating plants and combined heat and power plants.
In the literature, you can find many works on stratification studies. These works mostly concern the influence of geometrical parameters such as the shape and location of the inlet nozzle, the shape of the tank, and the use of additional plates improving the separation, at the height of the transition zone. The paper [1] presents the results of studies on the influence of the inlet nozzle shape on the thermal stratification. Kurşun and Ökten in [2] examined the dependence of stratification on the angle of inclination of a rectangular tank for different slenderness of the tank. The article [3] presents the results of numerical tests of the new design of the tank inlet. Similar tests have been presented in the paper [4].
Nowadays, there has been a growing interest in large hot water storage tanks to improve heat supply to consumers. The article [5] presents numerical analyses of two-dimensional temperature and velocity fields in vertical cylindrical large-volume heat accumulators used in district heating.
The article [6] presents an analysis of the possibility of using hot water pressure tanks to improve the * Corresponding author: pdzierwa@pk.edu.pl dynamics of a 200 MW power unit. Thanks to the installation of pressure tanks, the flow of energy supplied to the boiler together with fuel are used not only to produce electricity but also to heat water in pressure tanks.
This article uses an axial-symmetrical model to solve the task of charging the hot water tank. The use of a two-dimensional computing domain allowed to significantly shorten the simulation time in relation to three-dimensional models. The correctness of the applied calculation method has been confirmed by the good consistency of the results with the measurement data from the real object.
2 Technical description of the heat storage system In Figure 1 a schematic diagram of a heat storage system with a hot water tank is shown.  The thermal energy storage system is an indirect system and the water circulating between the heat exchanger and the tank is heated by network water, which is used to collect heat from the cooling process of gas engines. The role of the heat accumulator is performed by a vertical, cylindrical steel tank with a height of 20 meters and an internal tank diameter of 11 meters. The tank filled with water is connected in a closed circuit with a heat exchanger, which is a hydraulic separation between network water and water stored in the storage. This heat exchanger is a counterflow heat exchanger of 15 MW.
The flow of water between the heat exchanger and the heat accumulator is forced by a pump, which is mounted on a pipeline supplying heat storage medium to the battery exchanger. The maximum permissible volume of water flowing through the pump is 210 m 3 /h.
To control the direction of water flow in the heat accumulator circuit, automatic control valves marked with "AV" in the drawing are used. In case of heat absorption by cold water taken from the lower part of the tank, valves "AV1" and "AV4" are opened. The other two valves remain closed until the system operation mode changes. The mains water is directed to the battery exchanger from the discharges located in the supply pipelines of three main heat consumers.
After cooling down in the heat exchanger, the circulating water in the heating system is directed to a collective collector for all return pipelines from the heating networks of the mentioned heat consumers. From there, the water is directed to heat exchangers, where the network water is heated and directed to a distribution manifold on the supply.
The recovery of accumulated thermal energy takes place by means of heat flow in the opposite direction than in the case of charging. Hot water stored in the tank is taken from the upper part and flowing through the battery exchanger is cooled down and then returned to the lower part of the heat accumulator. The system is set to discharge mode by closing valves "AV1" and "AV4" and opening the other two valves on the secondary side of the heat accumulator exchanger.

Technical data of the heat accumulator
This hot water storage tank can be classified as an open heat accumulator according to the classification given in [7] and [8]. The heat storage medium is water. Hot water storage tanks are qualified as medium-term energy storage, but in this case, it is a short-term heat accumulator because it operates on a daily cycle. Due to the operating parameters of 98 o C/50 o C, the tank can be qualified as a medium-temperature thermal energy storage.
The working capacity of the tank is 1,800 m 3 . Due to the volume of the accumulated medium, the tank can be classified as an average heat accumulator. This heat accumulator can store 75 MWh of thermal energy. The cylindrical part of the tank has a height of 20 meters and an internal diameter of 11 meters. It is a welded construction made of structural steel S355. Table1 shows the basic technical data of the analysed heat accumulator.

Numerical simulation CFD
This chapter describes the process of preparing a CFD simulation for charging the hot water tank. The individual subchapters relate to the geometric model, physical, boundary and initial conditions, the calculation mesh, and the selection of time steps.

Geometric model
The heat accumulator system is a symmetrical system so that a flat design domain with axisymmetric boundary conditions is used ( Figure 3). A fixed height of the water in the tank equal to the height of the cylindrical part of the tank, i.e., 20 m, and a fixed internal diameter of the tank equal to 11 m are assumed. The analysis included thermal insulation of 200 mm thickness. For the sake of simplicity, the tank wall was omitted due to the insignificant value of thermal resistance. For thermal stratification to take place, plates with a diameter of 1 m and a thickness of 125 mm were placed 1 m away from the tank inlet and outlet. The inlet and outlet diameters are the same and are 250 mm. The length of the tank connection stubs placed in the lower and upper part of the tank in the axis of symmetry is equal to 5 times their diameter.

Mathematical model
To model the charging process, it is necessary to solve the transient problem of heat and momentum transfer. Due to conducted observations of momentum and heat transfer processes in the tank, it was assumed that the flow is laminar. The following simplifications were adopted in the calculation model:  the fluid is incompressible except for the density in the term representing the mass forces in the equation of motion, written for the x-axis direction,  the effect of changes in fluid viscosity is negligible,  axial symmetry of flow was assumed due to the installation of plates on the water inlet and outlet. The two-dimensional axisymmetric flow-heat problem is defined by the following mass, momentum, and heat balance equations

Physical conditions
The change of water density as a function of temperature, for normal conditions, has been approximated by the function: ( ) = −0.002767 � -0.17963 + 1,003.914 (5) where: ( ) -temperature dependent water density function, -fluid temperature.
The thermo-physical properties of mineral wool, which is the external insulation, adopted in the calculations, are given in Table 2. It is worth noting that the written dependence of water density (5) on temperature was considered in numerical calculations of the last term in equation (2) by means of the UDF (User Defined Function) procedure.
The variable ambient temperature is based on measurement data. A constant value of the heat transfer coefficient on the surface of the external insulation equal to 10 W/m 2 /K was assumed.

Initial condition
Based on measurement data from 20 temperature sensors installed in the tank, a 5th-order polynomial function describing the initial water temperature distribution was determined using the mean square approximation method. In water layers 1 meter high, located directly at the inlet and outlet stub, fixed temperature values equal to the temperature value from the first and last measurement point, respectively, were assumed. Figure 4 shows the initial temperature distribution of the tank water and the temperature values that were recorded by the measuring system on the actual object at the initial moment. The parameters of water at the inlet to the tank were established according to the measurement data taken from the SCADA system with a sampling rate of 1 minute. The values for each time step were approximated linearly between individual measurement points.

Boundary conditions
A time-varying water temperature and a varying mass flux of water entering the accumulator were assumed at the inlet to the accumulator. Heat losses through the cylindrical side surface were considered (boundary condition of the third kind). It was assumed that no heat transfer occurs through the top and bottom surfaces of the accumulator. No slippage was assumed on all walls, i.e., the water velocity on the wall is zero. An acceleration due to the gravity of 9.81 m/s 2 was assumed. The assumed boundary conditions are shown in Figure 3.
Changes of water parameters at the inlet to the tank and changes of ambient temperature as a function of time are shown in Figure 5. At the initial time, mass flow equal to 0 kg/s and inlet water temperature 357 K were assumed. After one minute, the inlet water temperature changed slightly from 363 K to 369 K. In the case of mass flow, the value changed from 10 to 40 kg/s. The highest heat output was delivered to the tank between 4 and 6 hours of the charging process. The flow then decreased to a value of about 20÷25 kg/s and slight changes were observed between successive time steps. The ambient temperature during the whole ±1 K.

Calculation mesh
The effect of mesh density on the solution was checked. The calculations were carried out for an element size equal to 125 mm. Then the calculation mesh was thickened twice in the radial direction and after the height, the number of elements increased by four times. Table 3 shows the dimensions and number of control volumes. A view of parts of the applied calculation meshes in the water domain in the middle of the tank is shown in Figure 6.

Time step
To determine the time step, the maximum water speed at the tank inlet has been calculated. This value was determined for the maximum volumetric flow. In the general case of two-dimensional flows, the Courant-Friedrich-Leve number (CFL) is determined by a formula using both speed components and both control volume dimensions. However, in a simplification, as shown in [2], a time step can be determined as for a onedimensional flow, because the maximum water velocity value is reached at the inlet stub of a tank into which the flow only takes place in one direction. Therefore, according to the formula for CFL, the time step is expressed by the formula [10]: where: ∆ ��� -maximum allowable time step, ∆height of control volume. The analysis of the above formula shows that for the same maximum velocity in the case under examination the time step is directly proportional to the length of the control volume. Four times the density of the grid results in a twofold reduction of the time step value. For each mesh, the time step was selected individually and the number of time steps that must be taken to simulate the entire battery charging process was determined using the formula below: (8) where: ∆� -number of time steps, -charging time of the heat accumulator. Table 4 summarises the results of the time step selection and the calculated number of time steps.        The effect of mesh density on the calculation was compared in Figures 13 to 15. The comparison was carried out based on an analysis of the average temperature of the tank water (Fig. 15), the density of the heat flux that is released to the environment through the cylindrical part of the accumulator (Fig. 16) and the thickness of the transition layer separating hot water from cold water (Fig. 17). The average heat loss stream is related to the external surface of the tank. To determine the thickness of the thermocline, it is assumed that if the water temperature is below 62 o C, then the water is in the cold water zone. The heat accumulator, on the other hand, is in the hot water zone when its temperature is higher than 92 o C.   Based on the analysis of the results for both meshes, very good compliance with the measurement results was obtained, which proves the correct selection of the calculation meshes.

Conclusions
The article presents the results of numerical calculations of the transient flow and heat problem, the aim of which was to model the process of charging the heat accumulator. The results of the CFD simulation were compared with the measurement data from a real hot water tank, which is in one of the Polish heat and power plants. In the comparative analysis, mainly graphic methods were used.
The assumption of the axial-symmetrical issue and the placement of simplified geometry of the inlet stub could have resulted in an increased degree of mixing of water layers with a small temperature difference in hot and cold water zones. Differences in water temperature distribution, especially in hot and cold water layers, do not cause heat stratification disturbances and do not have a significant impact on the thickness of the transition layer.
The chosen size of the control volumes has a direct impact on the results and on the value of the time step, which is also important for the actual calculation time. For both cases considered, good compliance with the measurements was achieved. This shows the correct selection of control volume meshes.