A Novel Data Management Methodology and Case Study for Monitoring and Performance Analysis of Large-Scale Ground Source Heat Pump (GSHP) and Borehole Thermal Energy Storage (BTES) System

: Aalto New Campus Complex (ANCC) is a recently inaugurated educational facility at Aalto University, located in Otaniemi (Espoo), Finland. Within over 40,000 m 2 , it comprises two faculties, a shopping center, recreational areas, and a metro station. ANCC is also a large-scale application of Ground Source Heat Pump (GSHP)–Borehole Thermal Energy Storage (BTES) in Finland, comprising an irregular BTES ﬁeld of 74 boreholes with an overall length of roughly 23 km and 4 million m 3 of energy storage. Therefore, accurate monitoring of the GSHP–BTES energy system is crucial for sustainable and efﬁcient long-term operation. Due to the fundamental issues affecting the accuracy of all thermal energy meters, a novel methodology adjusting for consistency of the measured data (in order to accomplish daily energy balance on both sides of the GSHP) is developed. The proposed methodology is used also in conjunction with reconstruction of missing relevant data before April/May 2020 by applying linear regression techniques. The developed data management is considered essential due to its capability to handle measured data with high uncertainty (thermal meters) by using highly accurate data regarding the GSHP power demand. Additionally, operational data and relevant GSHP performance indicators for the 18-month period starting from July 2019 is presented and analyzed.


Introduction
The energy footprint of buildings accounted for 29% of overall primary energy consumption globally in 2018, according to BP Energy Outlook 2020 [1].In Europe, heating and cooling represent roughly half of all energy used in buildings and industry.As a result, making heating and cooling more sustainable and efficient is a priority for the EU [2].In this context, ground source heat pumps (GSHP) in tandem with underground thermal energy storage (UTES) is an attractive technological option for efficient dispatching of heating and cooling loads, the integration of renewable energy sources (RES), and waste heat.They are interesting also when aiming for further decarbonization of the existing heating and cooling networks [3,4] and especially effective when applied in a centralized/shared way in a district level [5].The high potential of GSHP for efficient heating and cooling systems has fomented their exponential deployment and growth during the last decades.This has happened especially between 2015 and 2020 when their utilization almost doubled to 167 TWh/year and the installed capacity increased by 54%, reaching 77.5 GWt [6].
Finland is one of the top leaders in the world in heat pump utilization, with some 900,000 units, of which 140,000 are GSHP [6,7].However, most of these heat pumps are small units installed in individual dwellings or apartment buildings, with still very few medium-/large-scale realizations.In 2019, Kallio [7] estimated that there were only around 20-25 large-scale GSHP applications in Finland (applications with total borehole length over 10 km).
One of these large-scale GSHP-Borehole Thermal Energy Storage (BTES) successfully operating projects in Finland is the Aalto New Campus Complex (ANCC), with a total borehole length of 23 km powered by a GSHP with installed capacity of 0.8 MWt.The GSHP is the heart of the energy system, providing simultaneously heating and cooling all year round.Therefore, accurate monitoring of the GSHP operation and its interaction with the BTES field is crucial for efficient and sustainable operation in the long-term.
Spitler and Gehlin [8] highlighted the difference between the buildings' energy performance indicators in the design phase vs. how buildings perform in real conditions, denominated as the "performance gap".The difficulty in determining the concrete reasons for such a performance gap is an important question raised by Spitler [9].He considers that measurements and continuous monitoring of the operation of the energy system is essential for identifying the potential problems.There are many factors to be considered, with potential impact on the overall energy performance, such as a building's envelope and airtightness; quality and correct operation of heating, ventilation, and cooling systems; varying outdoor and indoor conditions (weather, occupancy, etc.); specific parameters related to GSHP-UTES energy systems (e.g., heat pump Coefficient of Performance (COP), source/sink temperatures, balanced interactions between the GSHP and the BTES field, and long-term impact), etc.Moreover, Spitler's vision of a widespread monitoring network for building performance combined with machine learning and data analytics can improve the operation, can increase the capability of fault detection, and can potentially benefit all implicated stakeholders in the near future [9].Consequently, the enormous network of all smart meters, data sensors, and measuring devices and the massive amount of data they continuously generate would require a multidisciplinary big data processing environment combined with artificial intelligence algorithms [10].
The increasing interest in evaluating the long-term performance of GSHP operation during the last decades has promoted different approaches such as the SEPEMO scheme (Seasonal PErformance factor and MOnitoring for heat pump systems) [11].SEPEMO defines four different system boundaries, e.g., the first boundary takes into account Heat Pump (HP) refrigeration cycle, the second one additionally includes the power demand of the circulation pumps on the source-side, the third also takes into account the auxiliary heating and cooling, while the fourth additionally includes the distribution-side circulation pumps and fans [8].The SEPEMO scheme was primarily developed for noncomplex residential heat pump systems (not specifically GSHPs).Despite the enormous efforts done so far, nowadays, there have been certain limitations of SEPEMO boundaries acknowledged [12], especially when dealing with sophisticated GSHP configurations delivering simultaneously heating and cooling all year round [8,12,13].Spitler and Gehlin [8] additionally highlighted the complexity to allocate the power demand of the air handling units (AHU) separately to heating, cooling, and ventilation, which for example influences the correct implementation of Seasonal Performance Factor (SPF) boundary 4.
The aim of the present research is to develop and implement a novel data management methodology to assess the performance of a complex GSHP-BTES system, the New Campus Complex of Aalto University, under the conditions of high uncertainty related to measured data.The complexity of such an implementation lies in the delivery of continuous and simultaneous cooling and heating on both sides (evaporator and condenser) of the GSHP while using different heat sources and sinks for balancing the system at the same time: the BTES field acting as both source and sink, and dissipation (heat sink) and district heating (DH) as additional peak heat sources (see Figure 1).The additional challenge related to fundamental issues in the accuracy of all thermal energy meters was accepted as a challenge and a starting point for developing an innovative methodology for measured data management.The consistency of correlated measured data of complex energy systems in buildings has been rarely addressed so far in previous research.Several studies have used measured data to calibrate and optimize building energy simulation models [14-16].However, none of them questioned the accuracy and checked for the consistency of the acquired measured data.Data validation and reconciliation (DVR) in industrial thermal processes is a technique addressed in several research studies [17,18].The role of DVR is to extract reliable and accurate information from raw measured data and to produce a consistent dataset resulting in the most probable process operation.DVR belongs to a more general category of data management [17].The DVR procedure normally intends to correct measurement errors and can be expressed mathematically as an optimization problem for optimally correcting the measurement data in a way that the adjusted values are consistent with the laws of conservation.The objective function is based on the least squares method with various sets of constraints, e.g., in thermal processes-constraints normally based on flow, mass, or energy balance [17].
The novelty of this paper is to introduce a specifically developed DVR methodology for data management adjusting the measured data for consistency in order to achieve the necessary energy balances on both sides of the GSHP.Application of the generalized DVR method, which handles raw measured data with unmeasured quantities (GSHP condenser and evaporator loads), tends to decrease the uncertainty after data reconciliation [17,18] and is a certainly indispensable solution to cope with the inherent high uncertainty of ANCC thermal energy meters.The present work also develops a method for reconstruction of some essential missing measurements before April/May 2020 by applying linear regression techniques based on different correlated measured data (Section 2.5).

Materials and Methods
Aalto New Campus Complex and its energy system is introduced in Sections 2.1 and 2.2 The measurement equipment is described in Section 2.3, and the methodology used for data management is presented in Section 2.4.Additionally, the linear regression models utilized for reconstructing essential missing measured data are presented in Section 2.5.

Description of the University Complex
Aalto New Campus Complex (ANCC), which was inaugurated in 2018, is located in Otaniemi (Espoo) at 60.186 • N/24.825 • E, only 10 km from Helsinki.ANCC is a polyfunctional 4/5-storey modern building (Figure 2) and includes the following main parts: the School of Arts, Design, and Architecture (ARTS); the School of Business (BIZ); as well as a shopping center, recreational areas, cafeterias, restaurants, and a metro station (METRO).

GSHP-BTES Tandem as Core of the Energy System
The heart of the ANCC energy system is the GSHP, composed of nine centralized heat pump modules operating in a cascade mode, with a total nominal capacity of 790 kW (at brine temperature 0 • C/output temperature 35 • C).Heat pumps generate simultaneously heating on the condenser side and cooling on the evaporator side (Figure 3).Most of the annual cooling demand is delivered as free-cooling energy directly from the boreholes.During the cooling-dominated season (summer), most of AHU-cooling demand is provided directly by the evaporator side of the heat pumps, and at the same time, all waste heat generated on the condenser side is dissipated.
The ground source heat pumps are connected to an irregular BTES field composed of 74 groundwater-filled boreholes with an average effective depth of 302.5 m and a diameter of 115 mm, containing a plastic PE single U-tube pipe with a diameter of 40 mm.The borehole field includes two 11 kW circulation pumps, each one capable of delivering 21 L/s at 200 kPa.Both pumps work in parallel connection by adjusting their output with variable frequency drives (VFD).
Additionally, the energy system is connected to a DH network, which is utilized as a peak heat source in addition to the GSHP.Domestic hot water (DHW) is directly generated by DH; thus, the GSHP is used only for heating.Overall, heating is delivered within a low temperature network (supply temperature 30-45 • C adjusted depending on outdoor temperature) through three different circuits: AHU heating, space heating and snow melting (via embedded pipes in several pedestrian zones).The maximum return temperature is set to 30 • C, and every time it surpasses this value (in summer), excess heat is dissipated.On the other hand, space cooling is provided all year round through radiant systems (supply/return temperatures 12-16 • C/15-18 • C).AHU-cooling operates only with high cooling demand (summer mode) and is provided mostly by an HP evaporator (average supply/return temperatures 10 • C/16 • C).The switching between summer and winter operation depends on several parameters such as the outdoor dry air temperature (ODA), cooling demand, temperatures of the cooling circuits, etc.It occurs normally with an ODA between 12 • C and 15 • C.  Electrical measurements are usually very accurate (with 0.5-1% accuracy), and this is the case of all power meters within electrical panels, such as the one located in the main Heat, Ventilation, and Air Conditioning (HVAC) mechanical room.All circulation pumps (CP) in the room are controlled by VFD, which also can acquire and measure power demand with an accuracy of 5%.The GSHP power demand is calculated by subtracting all circulation pumps from the overall electricity demand of the main HVAC mechanical room; thus, the measured accuracy of GSHP power demand is around 1-1.5%.

Thermal Energy Meters
Within the monitoring system of ANCC, there are 12 thermal energy meters measuring heating and cooling demand (both AHU-and space-demand) of 3 different zones, namely, the School of Business (BIZ), the Metro station (METRO), and both Schools of Arts and Business (ARTS + BIZ).Additionally, there are 3 thermal energy meters for measuring the BTES field, snow melting, and dissipation energy, respectively.The nomenclature of all thermal energy meters of ANCC is presented in Table 1.All heating and cooling demands of the ANCC have been monitored since July 2019.However, loads of the BTES field and dissipation have been acquired only since April/May 2020.All thermal energy meters are not calibrated yet, and there is a high uncertainty regarding their measurements due to the following reasons:

•
The accuracy of energy meters depends on different factors: errors related to the flow, sensors' delta T, calculator, quality of installation, proximity to disturbances and pumps, dirt in fluid, gas entrainment, sampling interval, etc. [19].

•
ANCC meters record data only once per hour regarding flow, temperature, and energy (instead of doing it more frequently, with reduced sampling interval).This can generate significant errors, which are difficult to estimate without calibration.

•
All ANCC meters are battery-powered; thus, the sampling interval can be an issue, since a much longer sample interval is used to reduce power consumption [19].Meters with slow response can cause significant errors around 20-30% [19].

•
The location of most of the energy meters is far away (sometimes hundreds of meters of pipe connections) from the main HVAC room; thus, heat losses/gains in distribution pipes and heat exchangers are not accounted for and can introduce additional errors.

Seasonal Performance Factors (SPF)
Currently, there is no possibility to separately measure the GSHP compressor and its internal circulation pumps.As mentioned previously, GSHP power demand is monitored as a difference between the overall power demand of the main HVAC room and all circulation pumps.Therefore, considering HP total power demand and the power demand of the BTES circulation pumps, the logical SEPEMO performance metrics is SPF with system boundary 2 [20].The SPF scheme (boundary 2) considers the net heating and cooling loads covered by the GSHP (plus BTES circulation pumps) and takes into account the free cooling from the ground (e.g., space cooling).Additionally, heat pump COP is calculated as a reference to the refrigeration cycle (the ratio between HP condenser output and HP power input).
The SPF (H2) performance calculation takes into account the HP net heat input to the system (the portion of total heating demand covered by the GSHP).It should be noted that this approach is crucial in the summer, since the GSHP is mostly used to deliver AHU-cooling and only a small fraction of HP condenser supplies net heat to the system (low heat demand); thus, most of HP condenser output has to be dissipated as waste heat.The measured performance for heat pump COP as well as heating, cooling, and combined SPF (boundary 2) is conducted according to the following equations: where ∑Q HP,cond , ∑Q cool , ∑Q heat , ∑Q DH,heat , and ∑Q net_heat are, respectively, the sum of thermal loads over the studied period for HP condenser, total cooling demand, total heating demand, DH input for heating, and total net heating demand (MWh); ∑W HP and ∑W CP,BTES are, respectively, the sum of power demands of the GSHP and BTES circulation pumps (MWh); while additional indices indicate the allocation of power demands for heating and cooling.
During the winter period, all heating loads are considered to be generated by the GSHP (its power demand is allocated only to heating), while power demand of the BTES circulation pumps is allocated proportionally to net heating and cooling loads.In the summer, the GSHP evaporator operates directly against AHU-cooling loads and covers a big share of them (thus, HP power demand is allocated proportionally to the covered AHU-cooling and net heating loads), while the rest of the loads are free-cooling (thus, the power demand of BTES circulation pumps is allocated entirely to space-cooling loads).

Calculation of Heat Pump COP
The measured values of GSHP COP determined with Equation (1) are compared with two different COP estimations based on the source and sink temperatures of the heat pump cycle, namely Lorentz COP [21] and Carnot COP.Reinholdt et al. [21] presented a calculation of Lorentz COP Lor , defined as follows: , where where T lm,H and T lm,L are the logarithmic mean temperatures of the heat pump condenser and evaporator, respectively.T HPC and T HPE are condenser and evaporator temperatures, while the notation I/O stands for inlet/outlet temperatures of the evaporator and S/R stands for supply/return temperatures of the condenser (expressed in Kelvin).For Lorentz efficiency η Lor , Reinholdt [21] suggested values between 0.5 and 0.6 (for very efficient large-scale heat pumps).A more conservative value of 0.39 was adopted for the ANCC case study.On the other hand, Carnot COP takes into account only evaporator inlet temperature T HPE,I and condenser supply temperature T HPC,S , and for ANCC case study, Carnot efficiency η Car = 0.43 is adopted:

Measured Data Management through Uncertainty Factors 2.4.1. Uncertainty Factors
A specifically developed DVR procedure with daily resolution is adopted due to the significant thermal inertia of both heating and cooling hydraulic networks, which in practice makes impossible to achieve, e.g., hourly balance on both sides of the GSHP.This resolution is good enough for depicting system daily operation; however, it is not able to capture more precise energy ramping or hourly peak loads.
The proposed particular DVR methodology is used in conjunction with a reconstruction of missing relevant data before April/May 2020 (Section 2.5), allowing for reconstruction of the whole energy system since July 2019.Therefore, the proposed methodology is essential due to its capability to handle measured data with high uncertainty (thermal meters) by using highly accurate data (power meters) as a robust procedure balancing the evaporator and condenser sides of the GSHP.
Four different uncertainty factors were introduced for all measured heating, cooling, BTES, and dissipation loads in order to accomplish a daily balance on both the evaporator and condenser sides of the GSHP.
Relative errors ∆ evap,i and ∆ cond,i are formulated on both sides of GSHP for each day i of operation: where Q cool,i and Q heat,i are, respectively, total measured cooling demand (negative sign) and total heating demand (positive sign) (MWh); Q evap,i and Q cond,i are the GSHP loads on the evaporator/condenser side (MWh) (positive sign, they are not measured and are determined with the DVR procedure); Q DH,i , Q BTES,I , and Q diss,i are the measured loads related to DH (net input for heating, positive sign), BTES (negative for heat extraction from the ground and positive for heat injection), and dissipation (positive sign), respectively, (MWh); and f cool,i , f BTES,i , f heat,i , and f diss,i are the uncertainty factors for cooling, BTES, heating, and dissipation.

Optimization Problem for Data Management
The optimization problem is defined in order to minimize the sum of squared relative errors ∆ evap,i and ∆ cond,i (eventually to 0) by varying the aforementioned uncertainty factors f cool,i , f BTES,i , f heat,i , and f diss,i and the loads on the condenser Q cond,i .It operates with daily values and is formulated on a monthly basis (n is the number of days in a month): Objective function: Subject to constraints: The initial values of all uncertainty factors are set to 1, while the initial values of HP condenser Q cond,i are chosen in order to fulfil the balance on the condenser side (∆ cond,i = 0).The optimization problem is solved using the Generalized Reduced Gradient (GRG) nonlinear method in MS Excel Solver.
The obtained uncertainty factors with the optimization procedure are used to adjust the loads by multiplying them with the corresponding uncertainty factors: total cooling demand Q cool,I composed of AHU-cooling Q AHU-cool,i and space-cooling Q spa-cool,i ; BTES loads Q BTES,i ; total heating demand Q heat,i composed of AHU-heating Q AHU-heat,i , spaceheating Q spa-heat,i , and snow melting Q snow-melt,i ; and dissipation loads Q diss,i .The adjusted data fulfills the necessary energy balances on both sides of the GSHP (evaporator and condenser, which loads are also determined) and therefore is consistent.The detailed procedure applied for data management is depicted in Figure 4.
An illustration of the specific DVR methodology for 12 December 2020 is presented below.The following measured data are obtained (daily loads, in kWh):

Reconstruction of Missing Data Measurements
There are three different types of essential missing data (before April/May 2020) needed to be reconstructed based on different correlated variables: the power demand of all pumps located in the main HVAC room, BTES field interactions (both heat extraction and injection), and dissipation loads (which in the summer have relevant weight in the energy balance of GSHP condenser).

HVAC Pumps
Normally, there are twin pumps in most of the circuits, such as those corresponding to BTES field, heating, and space-and AHU-cooling and dissipation.It was possible to obtain manually the totalized power consumption (kWh) and number of hours of operation for each pump from the automation system within the available period between 13 March 2018 and 16 June 2020 (Table 2).Table 2 presents a good insight of the system's operation, and despite its general and limited information, it is possible to draw some valuable conclusions:

•
Pumps of the BTES field account for more than half of the total consumption (51%).

•
All pumps on the BTES and heating circuits are used very intensively for around 64-70% of the time of the studied period, followed by space cooling (44-45%).

•
Pumps on the AHU-cooling circuit (used only during the summer period) are active for around 12% of the time, indicating that summer mode is very short.

•
Dissipation pumps are used very sporadically, only for 4% of the time.
The data with a 5-min resolution for HVAC pumps' power demand are only available from the middle of April 2020.However, there exists data with a 5-min resolution for the frequency percentage (FP) sent to each pump by the Variable Frequency Drives (VFD), which can be used to restore the missing power demands.The approximate linear correlation between pump power (kW) and the cube of the frequency factor x = (FP/100) 3 are well known.Therefore, they are utilized for linear regression based on the data between mid-April and mid-June 2020.The results are shown in Table 3, and considering the weight of each pump within total power demand, the overall standard error is close to 5%.If the outcome of any regression equation turns out to be negative, a minimum value of zero is taken.Using these regression models, missing HVAC pumps' power demand is reconstructed with an hourly resolution between July 2019 and April 2020.

BTES Field and Dissipation
The energy meter measuring the BTES field started operating properly only since May 21st, 2020.Prior to this, BTES energy data had to be estimated based on the frequency percentage of both pumps (FP1 and FP2) and the difference between flow and return temperature (∆T), as well as by clustering the hourly-based data depending on whether heat is extracted from the ground (negative sign) or injected into to ground (positive sign).The multivariable regression for the BTES field is based on ∆T and the combined frequency factor CFF = 3 FP1 100 3 , since CFF is approximately proportional to the pumping flow rate of the BTES circuit.The resulting regression equations as well as the corresponding R 2 value and the standard error are presented in Table 4. Again, if the outcome is of different sign than expected, a zero value is taken.The high correlations in both regressions (0.94 < R 2 < 0.96) are reflected by low standard errors between 13% and 16%.The energy meter monitoring the dissipated heat started acquiring data only since 23 May 2020.Prior to this, the regression model of dissipation heat rate was based on only one available data point: the entering fluid temperature TE 24 within the dissipation internal circuit (Table 4).The outcome can be only nonnegative (thus, for TE 24 > 31.8 • C).The lower correlation (R 2 = 0.732) implies that the standard error rises to over 42%.

Raw Measured Data
The raw measured data were retrieved from the automation system during 18 months between July 2019 and December 2020.All missing data were reconstructed according to the regression models introduced in Section 2.5 (estimations shown in blue font), and the results are presented in Table 5.From the initially measured data, it is not trivial to determine the evaporator and condenser loads of GSHP, fulfilling the energy balance on both sides of the heat pump.That is why the methodology introduced in Section 2.4 was applied to manage the initial raw data (Table 5) using the introduction of daily uncertainty factors corresponding to cooling, BTES, heating, and dissipation.Consequently, daily loads for cooling, BTES field interaction, heating, and dissipation were adjusted by the corresponding uncertainty factors and shown in Figure 5. Additionally, all daily loads corresponding to GSHP, i.e., evaporator, condenser, and power demand, are presented in Figure 6.
It is interesting that, during the winter period (October to April), the daily heating demand achieves peaks of 20-21 MWh/day and the heat extracted from the ground is 11-12 MWh/day.On the other hand, space cooling is needed all year round (0.5-1 MWh/day), while the intensive use of AHU-cooling in summer (June, July, and August), mostly provided by active HP cooling, can increase the total cooling demand to 14-16 MWh/day.The intensity of summer active cooling (reflected in HP activity, condenser peaks up to 12-14 MWh/day) combined with low heating demand (1.5-2 MWh/day) can provoke daily peaks in dissipation loads up to 12-14 MWh/day.

Results and Analysis of the Evaporator Side
The results of the evaporator side between July 2019 and December 2020 are summarized in Table 6 and Figure 7.The total cooling loads from June to August 2020 account for almost 2/3 of overall annual cooling demand, while AHU-cooling demand over these three months is almost half of annual cooling loads.On the other hand, from October 2019 to April 2020, there is no AHU-cooling demand at all, although space cooling is needed during the whole year.Therefore, it can be concluded that summer operation lasts roughly for three months (June, July, and August) while winter mode covers seven months between October and April, with two transition months in between (May and September).
Similarly, heat from the BTES field is heavily extracted from October 2019 to April 2020 (again, May and September are transition months with lower extraction activity), whereas heat injection occurs during the summer months of June, July and August.Overall, the hourly calculated imbalance ratio during the last 12 months is 7 (1869 MWh heat extraction/268 MWh heat injection).Thus, the interactions with the BTES field are not balanced, which can potentially lead to overcooling in the long term.The uncertainty factors for cooling and BTES are still far from 1 from July to September 2019 and the standard deviation is high (Table 6).After that, they are more stable, i.e., close to 1 and with low standard deviation (due to testing of thermal energy meters, April and May 2020 are exceptions).The symmetry of Figure 7 also reflects the balance of all loads on the evaporator side.In summer, BTES loads injected into the ground (ground acting as heat sink) plus HP evaporator balance out all cooling loads.For the rest of the year (from September to May), evaporator loads are balanced out by total cooling and BTES loads extracted from the ground (being the ground used as heat source).GSHP power demand in addition to the evaporator loads is also shown in Figure 7 in order to represent the magnitude of HP condenser loads.

Results and Analysis of the Condenser Side
The monthly loads on the condenser side from July 2019 to December 2020 are presented in Table 7 and Figure 8.In Figure 8, condenser and DH loads are shown as negative in order to highlight the symmetry with all heating and dissipation loads due to the energy balance.As mentioned previously, winter mode can be clearly observed from October 2019 to May 2020, accounting for over 90% of annual heating demand and HP condenser loads.Moreover, space-and AHU-heating demands are more or less evenly distributed over the year, while snow melting is practically inexistent.As expected, dissipation loads are realized almost entirely during the summer operation in June, July, and August, when they impact hugely on the energy balance of the condenser side (Figure 8).Similar to Figure 7, the symmetry of Figure 8 highlights the energy balance on the condenser side of the heat pump: heat input is provided by the GSHP condenser and DH balanced out by total heating output and dissipation loads (in summer, acting as heat sink).Overall, the GSHP is utilized very intensively, covering some 95% of total annual heat input, with DH accounting only for 5%.However, the past winter (2019-2020) was exceptionally mild in Finland [22]; thus, probably in normal meteorological conditions, the share of DH of overall annual heat input is expected to be higher.
Overall, the uncertainty factors for heating and dissipation are close to 1.However, from July to September 2019, their standard deviation is high (Table 7).After that, the factors are quite stable, close to 1 and with low standard deviation.However, there was testing of all thermal energy meters in April and May 2020 and space-heating metering presented malfunctioning on 5 and 6 July 2020.

Results and Analysis of the Uncertainty Factors
The initial raw data and the estimations of relevant variables such as BTES field interaction and dissipation can present high uncertainty, which is reflected in the uncertainty factors from July to September 2019, separated from 1 (blue axis) and with wide standard deviation (Figure 9).However, since October 2019, all factors are closer to 1 and have small standard deviation.Again, April and May are exceptions due to continuous testing of all thermal energy meters as are July 5th and 6th, 2020 when there was malfunctioning of space-heat metering (recording negative values).In addition, the BTES circulation pumps were not working for over 20 h.Therefore, the initial annual period for performance evaluation of the system was selected from October 2019 to September 2020 (highlighted in yellow, Figure 9).During this period, the overall cooling demand, BTES, and HP evaporator loads were 642, 1726, and 2368 MWh, respectively, while HP condenser, total heating demand, DH share for heating, and dissipation loads accounted for 3253, 3242, 171, and 183 MWh, respectively.During this initial period, the imbalance ratio (BTES field) reached the extremely high value of 7.5 (1993 MWh heat extraction/267 MWh heat injection).

Performance of the Energy System
Monthly Seasonal Performance Factors (SPF) with SEPEMO boundary 2 were calculated using daily measured data (all thermal loads were adjusted by the uncertainty factors) for the entire period.All heating (SPF-H2), cooling (SPF-C2), and combined SPF-HC2 factors are shown in Figure 10, as is the heating COP (H1).
Cooling SPF factors increase greatly during the winter period (between 50 and 70) due to direct free-cooling from the ground supplied as space-cooling.However, spacecooling demand in winter is a tiny share of annual cooling demand (Figure 7).As expected, winter SPF-H2 is lower than COP(H1) since the former accounts for BTES circulation pumps.The SPF-H2 factor also takes into account the net heating demand, which in winter is a difference between HP condenser and DH input (there is no dissipation), while COP(H1) calculates the heat pump COP (ratio between condenser output and HP power input).During the summer months (June, July, and August) SPF-C2 is lower, between 6.2 and 7.3.Heat pump is intensively utilized to meet almost half of all summer cooling loads, providing active AHU-cooling and simultaneously enhancing the effect of "free heating" on the condenser side.That is why summer SPF-H2 factors are higher than HP COP (Figure 10).Overall, some two-thirds of annual cooling demand is generated by free-cooling (assuming all space-cooling as "free-cooling" provided by the ground source).The overall heat pump COP and SPF factors for heating, cooling, and combined (system boundary 2) are respectively 3.7, 3.5, 9.3, and 4.0 for the entire period (July 2019-December 2020), while the annualized values over the first evaluation period (October 2019-September 2020) are respectively 3.7, 3.5, 9.3, and 3.9.For comparison, Spitler and Gehlin reported 3.6 and 6.4 as median values for SPF-H2 and C2, respectively, in their exhaustive overview of long-term measured performance of 55 GSHP systems [8].They also provided values of 3.7 (H2) and 27 (C2) for their own case study: the new student center at Stockholm University in Sweden [8].The latter factor (C2) is so high because cooling was supplied entirely as free-cooling from the ground peaking in winter between 40 and 60, since source CPs were proportionally allocated to the supplied heating and cooling loads (similarly to the present case study).Another recent study [23] presented an 8/10-year performance of two GSHP-BTES office buildings in Germany, reporting 2.6 and 4 as average SPF H2 and C2, respectively.Similar research for a 4-year GSHP-BTES performance in a German office reported 3.8 and 8 as typical heating SPF H2 and typical cooling SPF C2, respectively [24].

GSHP Entering Fluid Temperature (EFT)
One of the most important parameters of GSHP's operation is the entering fluid temperature (EFT) in the evaporator since it reflects how well the BTES field performs and gives a clue about GSHP's efficiency (since, generally, for a given condenser supply temperature, the higher the EFT, the better the HP COP).The number of hours within the first evaluation period (October 2019-September 2020) with a binned 1 • C HP EFT is presented in Figure 11.Some 90% of the hours of HP EFT is within a narrow range between 3 and 9 • C, especially from 3 to 5 • C over 60% of the time.It may be normal during the winter period, but in summer operation, more hours with EFT between 10 and 16 • C are expected, since the HP evaporator works directly against the AHU-cooling exchanger.In winter operation (ODA up to 10-12 • C), the HP evaporator works against the BTES field, and as expected, the EFT is higher than the fluid temperature flowing to BTES (which is the fluid leaving HP evaporator).However, in summer, the HP evaporator operates against the AHU-cooling exchanger and the EFT is consistently lower than the fluid temperature leaving the space-cooling exchanger.The latter confirms our observation that, in summer, the HP evaporator generates AHU-cooling at lower temperatures than necessary, which is also reflected in the heat pump COP (Figure 12), breaking the clearly ascending trend in winter.
It is important to acknowledge also that COP depends on both temperatures on the evaporator side (heat source) and the condenser side (heat sink).The GSHP fluid temperatures on both sides are plotted in Figure 13, where also the measured vs. expected Lorentz/Carnot COP (calculated with equations ( 5) and ( 6), η Lor = 0.39/η Car = 0.43) are compared.Generally speaking, the higher the evaporator temperatures and the lower the condenser temperatures, the higher the COP and vice versa.During the winter period (October 2019-May 2020), the GSHP condenser must provide heating at higher temperatures (average supply/return 42/30 • C) while the intensive heat extraction from the BTES field provokes a temperature drop in the evaporator (average inlet/outlet 4/1 • C).The calculated winter COP (orange line) fits well with the measured COP (green line) showing low plateau with minimum values around 3.5 between November 2019 and March 2020.On the other hand, in summer operation (June, July, and August 2020), GSHP condenser operates at lower temperatures (average supply/return 36/29 • C) and the evaporator operates at higher temperatures (average inlet/outlet 10/7 • C).This enhances the calculated summer COP to around 5. However, this theoretical expectation is not fulfilled by the measured COP (only 3.7, very similar to the average winter COP), presenting a huge gap between June and September (Figure 13).

Effect of GSHP Partial Load Operation
The significant gap between the measured COP and expected (calculated) COP in summer operation-as discussed in the previous section-could probably be provoked by the HP partial load operation.Several studies [25][26][27] concluded that heat pump COP is impacted by a partial load operation and highlighted that more frequent "on/off" switching can reduce HP efficiency [25].As mentioned previously, the ANCC system's operation is based on cascade launching of HP modules (from 1 to 9 units) and each module can be additionally adjusted to deliver a fraction between 0% and 100% of its nominal capacity.This effect of partial load operation is studied in Figure 14, where 366 daily data points (representing the first evaluation period between October 2019 and September 2020) are plotted depending on the average daily percentage of HP operating modules.Data are additionally clustered in three groups: winter operation (258 days, ODA up to 12-14 • C), summer operation (78 days, ODA over 15 • C), and transition period (30 days, ODA 12-15 • C).
In Figure 14, the average daily COP is plotted on the left axis while the secondary (right) axis shows the average daily ODA temperature (with cruxes).During the winter period, HP COP steadily decreases as modules operate at higher rates (between 25 and 100%).However, it can be seen that COP is closely correlated with ODA temperature, which is expected, since COP is lower with colder weather (during the winter days with high heating demand, EFT is lower and HP supply temperature is higher; thus, COP goes down).In summer, when GSHP modules work at lower capacities, normally up to 30%, COP sharply decreases, especially for an HP capacity below 10-15%.Eventually, a maximum COP (plotted as the decile binned COP, brown line in Figure 14) is achieved between 25% and 50% of HP average capacity, which mostly corresponds to the start/end period of the heating operation (Figure 14, highlighted in yellow).

Conclusions
A new methodology for the measured data management of the ANCC energy system was introduced and discussed in the present research.Data acquisition for 18 months, from July 2019 to December 2020, was processed, while essential missing data before April/May 2020 was reconstructed using linear regression techniques.
The presented methodology is considered crucial for tackling the existing uncertainty related to all thermal energy meters difficult to evaluate since they have not been calibrated and has not even been properly commissioned yet.With the introduction of four uncertainty factors related to cooling demand, BTES loads, heating demand, and dissipation, the applied method adjusts all thermal measured data for consistency in order to achieve daily energy balance on both sides of the GSHP and, at the same time, determining the loads of the GSHP evaporator and condenser.It is acknowledged that the aforementioned uncertainty factors present huge spreading with wide standard deviation during the first months (from July to September 2019), implying high uncertainty in the measured/estimated data.After this first period, the uncertainty factors become closer to 1 and their standard deviation steadily decreases.Therefore, the initial annual period for system performance evaluation was selected between October 2019 and September 2020.
The results for all loads on both sides of GSHP (evaporator and condenser) were reported and analyzed for the entire 18-month period, as were HP fluid temperatures, COP, and SPF (for system boundary 2).The HP COP, heating SPF (H2), cooling SPF (C2), and combined SPF (HC2) for the entire period are 3.7, 3.5, 9.3, and 4.0, respectively.

Figure 2 .
Figure 2. General schematics of Aalto New Campus Complex (image courtesy of Tekla, Finland).

Figure 5 .
Figure 5. Daily loads after data management procedure.

Figure 7 .
Figure 7. Monthly loads on the evaporator side.

Figure 8 .
Figure 8. Monthly loads on the condenser side.

Figure 9 .
Figure 9. Uncertainty factors and standard deviation (proportional to the diameter of the bubbles).

3. 7 .
Heat Pump COP Daily values for heat pump fluid temperature vs. outdoor dry air temperature (ODA) are plotted in Figure 12.The HP entering fluid temperature (EFT) and the fluid temperature flowing to the BTES field are compared in the figure.

Figure 14 .
Figure 14.Heat pump's COP and ODA in partial load operation.

Table 1 .
Nomenclature of thermal energy meters in Aalto New Campus Complex (ANCC).

Table 3 .
HVAC pumps: parameters for linear regression models and accuracy.

Table 4 .
BTES field and dissipation: Parameters for regression models and accuracy.

Table 5 .
Raw measured data and estimations (blue font).

Table 6 .
Data management: results of the evaporator side.

Table 7 .
Data management: results of the condenser side.