Improving the Gordon-Ng Model and Analyzing Thermodynamic Parameters to Evaluate Performance in a Water-Cooled Centrifugal Chiller

: The Gordon-Ng models are tools that have been used to estimate and evaluate the performance of various types of chillers for several years. A 550 TR centrifugal chiller plant facility was available to collect data from July and September 2018. The authors propose rearranging variables of the traditional (GNU) model based on average electric consumption and through a thermodynamic analysis comparable to the original model. Furthermore, assumptions are validated. Then, by estimation of the parameters of the new model using least square ﬁtting with ﬁeld training data and comparing to the GNU model and Braun model (based on consumption), it was shown that the proposed model provides a better prediction in order to evaluate consumption of a centrifugal chiller in regular operation, by improving the coe ﬃ cient of variation (CV), CV = 3.24% and R 2 = 92.52% for a ﬁltered sub-data. Through an algorithm built from steady-state cycle analysis, physical parameters ( S gen , Q leak,eq , R ) were estimated to compare with the same parameters obtained by regression to check the inﬂuence of the interception term in the model. It was found that without an interception term, the estimated parameters achieve relative errors (ER) below 20%. Additional comparison between external and internal power prediction is shown, with CV = 3.57 % and mean relative error (MRE) of 2.7%, achieving better accuracy than GNU and Braun model. validation, B.F.; formal analysis B.F. and A.B.; investigation, B.F.; resources, A.B.; writing—original draft preparation, B.F.; writing—review and B.F., A.B. and P.C.; visualization, B.F.; supervision, A.B. and P.C.; project administration, A.B. and P.C.; funding acquisition, B.F., A.B. and


Introduction
Nowadays, there is a great effort dedicated to improving heating, ventilation and air conditioning (HVAC) equipment efficiency, especially in chillers, which are equipment used extensively in buildings and industrial air conditioning installations. Therefore, it is crucial to evaluate their performance continually for optimization processes and failure prediction. This situation creates a necessity to generate more accurate models that allow simulation to predict their behavior. Traditionally, engineering models are classified into three categories: black-box, grey-box, and physical or mechanistic models [1][2][3][4][5]. These can be steady-state or dynamic models. The former has the disadvantage that only applies to the data range used to fit model parameters. The latter method needs exact information about the equipment geometry, materials, and data that usually only manufacturers can provide. That makes the analysis of a process very hard when it is not possible to obtain manufacturer data.
Many authors have made several thermodynamic models for chillers. One of the first research in the performance evaluation of a centrifugal chiller was carried out by Braun [6], who proposed a new empirical model based on consumption as a quadratic function. Subsequently, Gordon et al. developed 1.
The accuracy of the parameters obtained by the regression model.

2.
The internal and external power prediction of the model.
Both results can be valuable tools to evaluate the performance of the centrifugal chiller and eventually fault detection.
The paper is organized as follows. The first part contains mathematical and thermodynamic derivation of the models, description of the facility equipment and datasets, and thermodynamic analysis to estimate physical parameters through an algorithm and assumptions of the GNU model. From the training data, physical parameters were obtained and compared, using the regression method (with and without an intercept term) and the algorithm. Also, the accuracy of the model was investigated and contrasted to fundamental GNU and Braun models through statistical metrics.
The assumption of independence of residuals was checked. From the testing data, the accuracy and prediction capacity of the three models were compared over different sub-data chosen according to its setpoint supply temperature. Finally, some results and conclusions are discussed briefly.

Empirical and Semiempirical Models
Following the Ng-Gordon method, a chiller cycle is analyzed [7]. Figure 1 shows a simple-vaporcompression cycle for a conventional chiller plant in which all the measurable variables are shown. Some assumptions are required:  A total heat loss term, Q LT (kW), is defined in Figure 1. This term is calculated from Equation (1): where Q lev (kW) is the heat losses at the evaporator, Q lc (kW) is the heat rejected at the condenser, and Q lcomp (kW) is the heat losses at the compressor. Using the First Law of Thermodynamics for the cycle: Entropy net change for the chiller expresses as: where S gen (kW/K) is the irreversibility due to the entropy generated in the chiller, friction losses in the compressor, and the expansion device, the cooling load, Q ev (kW), it is the heat or heat transfer rate removed from the building through water flow. It can be calculated in several ways. From the heat exchanger where inlet and outlet refrigerant temperatures and evaporator geometry must be known, and from thermal capacitance (Equation (4)) approach: where m ch c ch is the thermal capacitance of cooling water (kW/K). With . m ch as mass flow rate (kg/s).
Similarly, for the removed heat in the condenser: It is assumed, in Equation (3), that heat exchange through boundaries occurs at constant temperatures. It is, therefore, necessary that both temperatures are in terms of measurable variables.
To estimate those temperature, according to the ε-NTU method: where Q max defined as: where C ch < C re f Due to overheating and subcooling being negligible, an isothermal process takes place in the evaporator and condenser. For the evaporator, the temperature can be expressed according to Equation (8): Similarly, for the condenser: With r c = ε c C cw and Replacing Equations (2), (8) and (9) in Equation (3): Furthermore, rearranging this equation and expressed in terms of A, B, C, D, and E coefficients: where, A = T cho r ev r c Q lev − Q lcomp + T cho r ev S gen Q lev − Q lc − Q lcomp + T cho T ci r ev r c S gen − Q lev r ev Q lev − Q lc − Q lcomp − Q lev T ci r ev r c (13) B = r ev Q lev − Q lc − Q lcomp + T ci r ev r c + r c Q lev − Q lcomp + T ci r c S gen + S gen Q lev − Q lc − Q lcomp (14) C = Q lev r c − T cho r ev r c − T cho r ev S gen (15) D = r ev + r c + S gen (16) According to Ng [8], these equations can be approximated based upon actual measured values for reciprocating chiller. There is no evidence in centrifugal chillers; for this reason, these terms will be estimated from real data. Due to evaporator and condenser are well insulated (pipe and shell heat exchangers), Q lc is negligible, and some terms related to the resistances (r c , r ev ) are neglected because of their low values.
A −T cho r ev r c Q lev − Q lcomp + T cho T ci r ev r c S gen − Q lev T ci r ev r c (17) B T ci r ev r c (18) Replacing Equation (1) into Equation (17), an expression is obtained for the average heat losses from the evaporator and condenser to the ambient. Ng and Gordon [8] call this term Q leak,eq Q leak,eq = Q lev + Q lcomp T cho T ci − T cho (21) A W comp expression is obtained by simplifying Equation (11): where, Gordon's model original parameters, in Equation (22), are maintained. But now the output variable is W comp . A fourth term is added. The expression can be rearranged to have the same parameter with physical meaning, S gen , Q leak,eq and R: Then, Equation (24) becomes: where,

Fundamental Gordon-Ng Model
Gordon-Ng (GNU) is a simple thermodynamic model derived from the First and Second Laws of Thermodynamics. This model applies to all chiller types [9] and is based on the coefficient of performance (COP) as a linear function of refrigeration load and temperatures that easily measurables in a cooling plant. The model takes the following form: where b 1 , b 2 and b 3 have physical meaning: b 1 = S gen = Internal entropy generation. b 2 = Q leak,eq = Rate of heat losses or gains from the chiller. b 3 = R = Total heat exchanger thermal resistance.

Biquadratic Braun Model
The model proposed by Braun is an empirical or black-box one, also called the biquadratic model. Predicts a standardized power consumption based on only two inputs with five coefficients and interception term (Equation (27): where X = Q/Q des and Q des is designed cooling capacity of the chiller. Z 0 is the designed input power of the chiller. As such,

Materials and Methods
The chiller plant operates at the Engineering Building at Universidad del Norte, Colombia, and manufactured by York ® (Optiview TM Software C.OPT.01.08A.300 version 10.0, Johnson Controls, PA, USA). The plant comprises one centrifugal chiller unit of 550 TR, four water-cooling towers, and a cold-water loop with constant primary water flow rate through the evaporator. The condenser circuit has constant-flowrate. The centrifugal chiller is a single-stage compressor type with prerotation vanes (PRV) for capacity control and a variable speed driver. A proprietary system, Metasys, manage inputs, and measured outputs. The chiller working fluid is R134a. Technical data can be seen in Table 1. The OptiViewTM York control center provides temperature and pressure data. Pressure transducers and thermistor with DC voltage as output mounted in various locations have different ranges. The uncertainties of the measured values are listed in Table 2. The data used for the model correspond to a stable state, with the restriction that the temperature values do not change above 1% during a time interval of 10 minutes. Table 2 shows the data intervals collected during July and September 2018. Data filtering was done as follows: • Wrong data due to sensor malfunction or system errors.

•
Unusual data is rejected: negatives and zeros values and large power values.

•
Values related to all aspects of the dynamic system: Powering up and down of equipment. • Data, according to Standard AHRI 550/590 [27].  Original data was collected through de building automation system (BAS), corresponding to 4464 data points for July and 4320 for September. Both months have highly cooling load dynamics. September will be used to internal prediction and July for external prediction. Both data have been filtered for two ranges of outlet water chilled temperatures relatives to Standard rating Air conditioning, heating and refrigeration institute (AHRI) conditions between 4.4 • C and 8.9 • C. July has been averaged to build a load and consumption profile. Table 3 shows the number of data points for every group of data after filtering.

Estimation of Thermodynamic Parameters
Through a thermodynamic analysis, the b and c values for Equations (1) and (27) are calculated. The estimation of parameters is made as closely as possible of steady-state operation with data collected. Some assumptions must be included, even a superheating value (not allowable from data). This analysis requires next assumptions: Figure 2 is saturated liquid at evaporator pressure.   From First Law of Thermodynamics applied to evaporator, condenser, and compressor, the following expressions have been obtained: Also, Q c from Equation (6): Based on next assumptions: h 1 = h sh Enthalpy at superheating temperature at evaporator pressure. h 3 = h f @Tc,sat Enthalpy of saturated liquid at condenser saturated temperature. h 4 = h 3 Assuming isenthalpic-throttling device. h 2 = h dis Enthalpy of overheated vapor at discharge temperature and condenser pressure.

As
. m re f is unknown, an algorithm to determinate Q lev , Q lc , Q lcomp , and . m re f converging from a seed value is developed. Parameters Q leak,eq , R, and S gen can be calculated from Equations (21), (23), and (3). The flow chart of the algorithm is shown in Figure 3. This program fastening calculates properties in different states from temperature and pressure data. Then, from the energy balances compute the heat leaks, loads, and consumption, through an iterative cycle assuming an initial value of refrigerant flow mass for residual convergence of 10 −5 . Finally, estimates parameters set by the below Equations and Formulas (31) and (32).
Effectiveness from Equation (6) can be expressed for both heat exchangers, as:

Regressor Parameters Estimation
Regressor parameters from Equations (24) and (26) can be estimated by regression analysis. Although some authors have studied the use of Time series or Non-linear regression [15,25], in the first case, there is no enough data by year to analyze a seasonal behavior; in the second case, the author concluded that there is not a significant practical improvement by using a Non-linear procedure. It should have been noted that both Equations ((24), (26)) are implicit, so the response variables (W comp , COP) are not corresponding to the response model variable y. In order to estimate physical parameters with moderate computational effort, an Ordinary Least Squares method has been used.
On the other hand, there is a discussion about the inclusion of the interception term. Reddy recommends using an interception purely for prediction purposes and not to use it for estimate parameters [4]. In order to investigate it, GNU and GNU rearranged models have been fitted in two manners, one corresponding to Equation (25), and another one including the interception term: Equation (25) is regressed for both models to estimate physical parameters and compare them with parameters obtained by thermodynamic analysis. Then, Equation (33) is regressed for both models to evaluate internal and external predictions. Next metrics have been used to compare: The square correlation coefficient R 2 and R 2 adjusted: where y mi is the modeling (or predicted) value, and y ei is the experimental value of the i-th data, n is the number of data points and p the number of parameters to be regressed. The mean relative error (MRE) and the relative error (RE) for parameters: The predictive accuracy can be assessed considering the coefficient of variation (CV) of the root mean squared error (RMSE), which indicates how well the model satisfies the predictions.

Results and Discussion
From the September data, a profile can be built to assess the chiller behavior and contrast it with thermodynamic analysis. Three variables are evaluated: cooling load (calculated with Equation (5)), coefficient of performance (COP), and the consumption per cooling load (In TR). The September profile is depicted in Figure 4. Note that it seems to be two regimes of load, close to 0.2 and close to 0.7. A histogram ( Figure 5) shows this behavior and essential information about the distribution of data for model assumptions. There is a normal distribution biased to the left relative to points of low cooling load.  In order to assess the behavior of performance and efficiency according to cooling load regimes, it should be added an efficiency profile (in kW per TR) to load profile for a sub-data of September. Figure 6 shows a recurrent profile per day; at the beginning and the end of the day, the chiller operation is at low load, but most of the daily operation is at partial load, approximately 80% of the design load which is a typical centrifugal chiller performance design [28]. At its maximum and partial load above 60%, efficiency is good (close to manufacturer value: 0.55), but it drops at 20% of the cooling load. Although the rearranged model is still an implicit equation, this structure can be useful to analyze further the effect of some input variables that have been collected as data, which indeed could be integrated as sub-models for building loads, pumps, and cooling towers. Cooling load has been calculated from Equation (4). Temperature values are sensed, mass flow is given by manufacturer, and the pressure drop for water is not allowable. However, further work needs to include variable mass flow, and a model for heat pump can be integrated to Equation (22) through Q ev and T cho .
Similarly, T ci is a crucial variable linked to a cooling tower. Further work is dedicated to analyzing the effect of the condenser water inlet temperature through a cooling tower model integrated to Equation (22). Some research could be done in the direction of estimate the consumption of a whole chiller plant.
It is crucial because, from Figure 5, it is evident that the chiller operates much time at partial load and is necessary to analyze equipment irreversibilities in order to improve the performance of this centrifugal chiller. To evaluate COP and consumption, scatter charts have been plotted (see Figure 7, Figure 8).   From plots 7 and 8, it is observed that according to theory [7,8] at low cooling rates, it seems to appear linear behavior between cooling load and COP dominated by losses as fluid friction and throttling. At high cooling rates, the curve is dominated by losses from heat transfer. Input power is proportional to the cooling load above 50%, but it is not clear at low cooling rates. Chiller performance is related to the three physical parameters. Those can be estimated by analysis and regression. Table 4 shows the values for every model regression applied to data and sub-data filtered with their corresponding ER%. For all sub-data, there is a good accuracy of the parameter R with ER% below 10%. Q leak,eq and S gen show poor accuracy; this is because of the strong dependency of the cooling rate that makes it hard to assume a constant term. Centrifugal chiller control by variable speed driver (VSD) and prerrotation vanes (PRV) influence S gen generating significant friction at low cooling load. The better accuracy has been achieved close to T cho set point (6.7 • C) when it is expected a more stable operation. See Figure 9 to observe the behavior of the physical parameters for the sub-data of September. Q leak,eq and S gen are highly dynamic. Values listed in Table 4 are averaged values; in the same way, the regressor parameters are assumed as constant, but they are not. However, the heat exchanger's thermal resistance is very stable, with few critical peaks which correspond to cooling load peaks.  Centrifugal chillers have been used to provide large cooling in buildings, and their efficiency is high, as observed in the present work. Even though they seem to have stable behavior, they are exposed to dynamical changes in loads and temperatures so that any steady-state model will be inaccurate.
For the 550 TR chiller plant facility, not all the data was available due to limitations in the control system, and for that reason, many assumptions were made to develop a simple thermodynamic model. From Figures 6-8, it is clearly noticed that there is a segmented range of operation with different load regimes. Because of that, different behavior could be expected for each segment. However, it is possible to find a model with enough robustness for both: low and partial load.
The fundamental model GNU depends on assumptions made to its development and final terms. Original authors [8] validate the approximations for A, B, C, and based on experimental measures for reciprocating chillers. The rearranged model proposed in this paper applied to centrifugal chillers should be validated, too, for the same approximations. Table 5 shows the values of the coefficients denoting the suitability of the model.

Validation with Training Data
Training data extracted from September data has been used to fit three models in order to compare their internal prediction power: GNU, GNU rearranged, and Braun. Metrics R 2 , R 2 adj, MRE, and CV are shown in Figures 10 and 11. Braun model has a CV of less than 5% except for the larger data. MRE is less than 10% for all the sub-data, which demonstrates good accuracy for Braun empirical model. However, parameters of this model do not have physical significance for thermodynamic performance purposes; they only would be useful for external prediction or fault detection and diagnosis. In the case of semiempirical models, MRE is less than 10% for both GNU and GNU rearranged. Also, both of them achieved CV less than 5% for the sub-data close to the setpoint, but GNU does not have a good square correlation R 2 (41%), while GNU rearranged has an excellent R 2 (97%). In all cases, the Braun model has better accuracy than the rearranged model. This result fits with the good performance of the model in the literature. However, the Braun model is empirical, and its coefficients have no practical meaning, which makes the Braun model an excellent predictive tool but less useful to evaluate performance through physical parameters. (see Figure 11).  From Figures 9-11 is possible to analyze the relationship between the physical parameters and accuracy of the GN models. As the subdata is close to the set point the more the model is accurate. Values close to setpoint is reflected in a stabilized behavior of the parameters S gen , Q leak,eq and R. another consequence is the assumption of constancy of those parameters linked to the steady-state framework of the model.
Other works over centrifugal chillers [10] show better metrics for the GNU model (CV below 5%) applied to the same capacities (550 TR). The water flow rates are also constant, but they are not explicit about the type of control used for partial loads. GNU fundamental model was not designed over centrifugal chillers with variable speed drives, and inlet guide vanes control used at the same time, and it is not accurate enough when dealing with dynamic changes caused by these types of control. Variation in pressure and flow rate must be reflected in COP, which is an instantaneous parameter. By modifying the GNU model, consumption and cooling load can be obtained as an average of cumulative values, which are measured in a period and are not instantaneous.
In addition, changing variables help to stabilize unbiased data. An example of this can be seen in Figure 12, where the independence of residuals for both models is depicted. Based on metrics and predicted values (Figure 13), GNU rearranged model achieves good internal accuracy to estimate energy consumption and to analyze performance.  The three physical parameters: internal entropy generation, rate of heat losses or gain from the chiller, and total heat exchanger thermal resistance, were compared with the values obtained by an algorithm in order to examine the influence of including the intercept term. Results are listed in Table 6 next to ER%. It should be noted that the GNU rearranged model as well as GNU model do not estimate the physical parameters when the interception term is included. From Table 6, it is evident that ER% is larger than 100%. Furthermore, there is a negative value for the entropy generated by the GNU model. Summarizing, models with interception terms are not useful to estimate parameters, nevertheless, they can be used for prediction.

Validation with Testing Data
Testing data for July 2018 was obtained from the BAS system, it is not used to fit models, but it is collected from the same equipment. For that reason, it is used to validate the prediction capacity outside the training data. Every sub-data was averaged and tested. The results are listed in Table 7.
It would be expected to achieve a better external prediction accuracy where the internal prediction was the best. Results show that a GNU customized model estimates consumption with better accuracy than other models for testing data, being CV% improved for extensive averaged testing data. MRE is less than 10% for both GNU models and below 5% for GNU modified. Braun model does not have good accuracy for testing data. Additional research must be done in this area. Figure 14 plots the power prediction of each model.  The two GNU models show similar behavior at low cooling rates (reflected as power consumption) except the Braun model that overestimates values. Dynamic changes are not well captured in the whole operation. GNU rearranged overestimates values of consumption at partial and high loads, especially peak values.

Toward Validation with Any External Data
Empirical based models are applicable only for the same range of data from which they were fitted. However, it is valuable to validate the proposed model for the same centrifugal chiller under another set of conditions. Future work is in the direction of getting a new data set from another month of the year 2018. There is an interest in investigating the behavior of the model after and before maintenance routines. As has been set the scope of better accuracy, the proposed model should include the parameters regressed listed in the first row of Table 6. The proceedings are as follows:

Conclusions
In this paper, a thermodynamic modeling analysis was proposed to analyze, predict, and optimize the performance of a centrifugal chiller. According to a simple thermodynamic analysis (several assumptions were made), a steady-state model was built based on easy measured variables and three physical parameters which can be obtained by regression. The three parameters were also estimated through an algorithm built from a thermodynamic analysis of the steady-state cycle. After comparation between the fundamental GNU and the proposed model was found that the proposed model showed less ER% for the thermal resistance of heat exchangers (less than 5%) and ER% below 20% for entropy generation and leak energy, for a sub-data close to the setpoint value, 6.7%. Two sets of data were used, a training data corresponding to September and testing data (split in two) from July in order to validate internal and external prediction, respectively. It was also proved that the appropriate inclusion of an intercept term is valuable if the aim of the model is only for prediction purposes. The validation of assumptions made by GNU model along centrifugal chillers, was made. however, further experimental work should be made and contrasted with the results obtained here.
Analysis indicates that (for internal data), the accuracy of the proposed model is better than the fundamental GNU and close to the Braun model. With R 2 and CV of 92%, 40%, 97% and 3.24%, 3.50%, 2.14%, respectively. The model achieved better accuracy as the sub-data is filtered around the setpoint value, where the operation is stable. External validation through an average testing data denoted good accuracy, with CV% less than 5% comparing to GNU and Braun models, but it tends to overestimate high cooling loads. It is essential to keep in mind that models are developed under steady-state and ideal conditions, without considering control routines/devices. Further work will be done in this direction. Also, it is necessary to develop dynamic modeling for vapor compression cycles considering the simultaneous effects of inlet guide vanes and the variable speed driver. Funding: This research was funded by Colciencias through Convocatoria Doctorados Nacionales 727 of Fondo Nacional de Financiamiento para la Ciencia, la Tecnología e Innovación FCTeI, by Sistema General de Regalías SGR and Universidad del Norte.

Acknowledgments:
The authors would like to thank Universidad del Norte and Colciencias through Convocatoria Doctorados Nacionales 727 of Fondo Nacional de Financiamiento para la Ciencia, la Tecnología e Innovación FCTeI, by Sistema General de Regalías SGR for supporting this research.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study, in the collection, analyses, or interpretation of the data.