Lowering greenhouse gas emissions in the built environment by combining ground source heat pumps, photovoltaics and battery storage

Ground source heat pumps (GSHPs) have been suggested to replace gas-based heating in urban environ- ments to reduce greenhouse gas emissions and help to comply with the Paris Agreement. The emission reduction from GSHP depends on the carbon intensity of the electricity generation mix. Moreover, grid capacity may be limiting the introduction of these high-electricity demand GSHP systems. Photovoltaics (PV) systems help to provide additional emission reductions for residential GSHP systems. Battery en- ergy storage systems can reduce the peak demand and allow for more GSHPs within the low voltage grid. We developed a techno-economic and environmental assessment model to quantify this impact of PV and batteries combined with residential GSHP systems. Measured demand data of 16 dwellings with GSHP and PV systems from the Netherlands were used. We show that PV can provide around 19% of the GSHP demand, while batteries enhance this by 53% and reduce the peak demand by 45%. Greenhouse gas emission of a GSHP with PV is reduced on average with 73 tCO 2 -eq, corresponding to a 80% reduction, over a 30-year lifetime. Dwellings with only a GSHP system have a net present values increase of around € 275 per tCO 2 -eq of avoided emission. This is reduced to € 230 per tCO 2 -eq when PV and storage is added to the system. Nevertheless, investment in GSHP systems today is not economically attractive for many dwellings. A sensitivity analysis showed that policies should focus on increasing natural gas tariffs, carbon taxation, investment subsidies or combinations of these routes to encourage sustainable heating.


Introduction
In the European Union (EU), around half of the buildings are provided with heat by fossil fuel boilers with an efficiency of 60% or lower and were installed before 1992 [1] . For example, in the Netherlands, 93.7% of the dwellings is heated using natural gas resulting in 11.1% of the total Dutch CO 2 emissions [2,3] . In line with the Paris Agreement, the Dutch government has set ambitious goals to phase out gas boilers by 2050 and replace them with other technologies to provide heat [4] .
Ground source heat pump systems (GSHPs), also known as ground coupled heat pumps, constitute a promising technology to replace fossil based heating systems [5] . GSHPs generate heat from electricity with high efficiencies and are seen as best available technology, especially in combination with renewable sources [1] . Emission reduction is largely dependent on the electricity generation mix in a specific country, GSHPs efficiencies and climate conditions [6][7][8][9] . Moreover, GSHPs can deliver flexibility such as demand response services to the electricity system [10] . These advantages led to policies that support investments in GSHPs. For example, in the USA a tax credit of 30% of GSHP investment costs presently exists [11] . In the Netherlands, a variable subsidy based on the installed GSHP capacity can be obtained [12] . In Europe, around 10 0,0 0 0 GSHP units are annually installed [13] . In the USA, over 560,0 0 0 units were installed by the end of 2014 [14] .
Currently, the Dutch electricity generation mix is relatively carbon intensive (490 gCO 2 -eq per kWh), which lowers the emission reduction potential of GSHP [15] . Rooftop photovoltaic (PV) systems are a worthwhile option to lower the electricity needs from the grid. Consequently, greenhouse gas (GHG) emissions can be lowered by avoiding electricity generation from fossil based power plants [16] . The direct use of PV electricity (referred to as self-consumption) is limited, due to a mismatch in time between PV production and electricity consumption. This can be in- creased by using a stationary battery energy storage system (BESS) that charges surplus PV electricity so that it can be used on later moments. The electrification of heating and the use of energy storage are highly recommended technologies to allow for more intermittent renewable energy in the electricity system worldwide [17,18] . The implementation of GSHPs could be restricted by the capacity of existing low voltage utility grid, due to the larger peak demand of dwellings with GSHP [19] . BESSs are suitable for peak shaving of the power demand and PV electricity production [20] . The latter application reduces potential energy losses due to PV curtailment requirements [21] . Consequently, more PV systems and GSHP systems can be installed on a local grid without expansion requirements. This reduces the societal grid costs and helps the profitability of BESSs [22] . Furthermore, less power generation capacity is required to meet peak electricity demand, especially for colder winter months. These associated benefits may result in a rapid deployment and cost decline of BESSs [23] . However, BESSs have charging and discharging losses that result in higher system emissions. These systems are only recommended if the share of renewable electricity generation is a large share of the total electricity generation [24] .

Literature review
Several studies assessed the technological, economical or environmental advantages of GSHP systems combined with PV and storage. PV self-consumption of dwellings with GSHP systems show a clear seasonal effect, with significantly higher selfconsumption in summer months than in winter months [25] . A GSHP control algorithm that used weather forecasts showed a limited increase of 7% in PV self-consumption for a Swedish building [26] . Another study using a residential dwelling from Switzerland found that controlling a heat pump on the availability of surplus PV power increases self-consumption with 1.5% [27] . Also, a few percent lower self-sufficiency by PV electricity was found when a heat pump was included for a residential dwelling in Germany [28] .
A German study developed a mixed linear programming model to optimize size the PV system, thermal and battery storage for residential dwellings with a heat pump [29] . It was found that the optimal battery storage size is mainly determined by the electricity demand and hardly influenced by the PV system size. Another study presented a multi objective model to reduce the GSHP consumption peak, using thermal storage and demand response for residential buildings in Belgium [30] . They found that the peak demand can be decreased with 2.5 kW per building, at an estimated capacity cost of 25 € /kW related to a lifetime of 25 years.
Combining a GSHP with a PV system was found to be a preferred option from an energy and economical point, then combining GSHP with a solar thermal collector. This is mainly because for each kWh of PV electricity produced, a multifold of heat can be produced with a GSHP system [31] . A study modelling residential heating options for buildings in Belgium found that GSHP reduce annual emissions with ≈ 1.5 tCO 2 . With the combination of a PV system, the annual emission reduction was 0.6-0.9 ton higher [32] . Also, higher CO 2 abatements cost were found for GSHP than air source heat pumps [9] . When surplus PV electricity is converted to heat using a GSHP and be exported to a district heating grid, then larger PV system capacities are economic suitable [33] . Also conversion of surplus PV electricity to heat was a better environmental option than exporting the electricity to the grid. This is related to the higher carbon intensity of district heating compared with the electricity generation mix [34] . A study conducted in the south of Spain compared 4 systems for provision of cooling heating and power demand to building. A system including an absorption chiller, an auxiliary heater, PV modules and solar thermal modules showed the best performance from an economic and environmental perspective [35] . Another study conducted in Italy found a 70-80% reduction in primary energy when a GSHP system is coupled with a PV system [36] .
The literature illustrates that GSHP systems combined with PV and storage show great potential for emissions reduction. Most studies focus on a single or a few research topics, and a broader integrated study combing technological economic and environmental impact was not found. These multidisciplinary studies are essential for a better understanding of the broader impact of PV and storage on GSHP systems. Furthermore, almost all studies use building models to assess the impact of these systems, since measurement data is difficult to obtain. If measured data was used, than only a single or a few buildings were included in the studies. Therefore the influence of real heating demand patterns is not well known.

Research aim
In our study, we present a broader integrated techno-economic and environmental impact for residential dwellings with a GSHP combined with PV and BESS. While most studies solely used a single or modelled consumption timeseries, our research used 16 measured residential electricity consumption and PV production timeseries with duration of 2-years at a 15 min time resolution.
We assessed the impact of 5 systems architectures by comparing the systems with a reference case over the 30-years lifetime. This reference case consists of dwellings with a condensing gas boiler (CGB) and no PV or storage installed. A battery control strategy was developed that aims to increase self-sufficiency and reduction of grid impact simultaneously. The economic impact was assessed and we provide policy options to increase the economic attractiveness of GSHP systems.
Avoided life cycle GHG emissions and payback periods for GHG emissions were determined for each system and compared with each other. The sensitivity of the electricity carbon intensity on the avoided emissions was assessed. Furthermore, we compared the results with surrounding western European countries. Consequently, our results could be used as an indication of emission reduction potential of GSHP combined with PV and batteries for these countries. The obtained knowledge is valuable for a broad range of users, from system owners, installers, distribution system operators and policy makers.

Method
This research assessed the impact photovoltaic and storage on residential dwellings using a techno-economic and environmental assessment model. First, electricity consumption and PV electricity production of 16 residential dwellings were measured and selected based on data availability and reliability. Also remaining model input parameters required were obtained from literature. Next, a battery storage model was developed which increased PV self-consumption and reduced the peak power on the local electricity grid. Then, technical, environmental and economic performance of each system design was assessed over their lifetime. We selected a 30-years lifetime based on the minimum expected lifetime or a renovation cycle of a residential dwelling [37] . Replacement cost of system components with a lifetime shorter than 30 years were also included. Finally, a sensitivity study on the input parameters was conducted. An overview of the used input data and model steps with corresponding sections and chapters is shown in Fig. 1 .

Energy production and consumption time series
An overview of the PV electricity production and energy consumption data selection is provided in this section. Also, key parameters containing statistics of PV electricity production and energy consumption from dwellings is provided. Further insights on the time series and patterns are given in Appendix A .
Rooftop PV system production and residential electricity consumption timeseries of 38 dwellings were measured in the Netherlands from 2013 until 2015. This was measured during the project called Your Energy Moment which assesses the influence of dynamic tariffs on residential demand. Data of this project has been validated and used in previous studies [38,39] . The electricity consumption of the dwelling, the GSHP system and PV electricity production were measured separately with a 15 min time step. The measured data from 1st of July 2013 until 30th of June 2015 was selected to obtain a period of 2 years. Only timeseries with a data availability of higher than 96% were selected. Furthermore, time series were manually analysed on incorrect measurements which reduced the dataset to 16 dwellings. Remaining missing data points were refilled using data from a similar time moment of the first available previous day or following day. As a result, each used timeseries contains 70,080 data points.
The dwellings are detached and semi-detached family houses, build in 2012, with highly insulated roofs, walls and windows, which are representative for typically newly build houses in the Netherlands. They contain PV systems that are oriented between 170 °and 210 °and have a module tilt of 10 °. Further technical and environmental performance statistics of the PV systems are provided in Fig A.1 .
Each dwelling holds a GSHP unit with a nominal output of 4 kW thermal. The GSHP provided the dwellings with three services: space heating, domestic hot water (DHW) provision and cooling. The GSHP system power requirements of the space heating are 1.1 kW and 1.8 kW for DHW production. The heat pumps are connected to a closed loop vertical soil heat exchanger that was designed to achieve a relative high brine temperature of 4 °C. Dwellings are heated using low temperature floor heating (28 °C). Combined with the high brine temperature, this results in a coefficient of performances (COP) of 5.7 for space heating. DHW production requires a higher temperature (60 °C), which results in a lower COP of 1.9 [40] . DHW is stored in 190 L storage tanks, which are refilled if reaching a lower water level. The GSHP can also cool the inside temperature with a few degrees using the floor heating system. Therefore, cool water from the ground source is circulated through the floor heating system. If cooling is applied, the GSHP uses ≈ 0.15 kW for the circulation pumps. In addition, the GSHP system uses electricity for standby use.
The minimum and maximum GSHP power consumption for each time step were measured. Consequently, we could separate the GSHP consumption time series into four profiles for each function, based on the maximum energy consumption measured for each time step. Electricity consumption below 0.08 kW was allocated to standby usages. Consumption between 0.08 and 0.2 kW was allocated to cooling use. The consumption between 0.2 and 1.4 kW was allocated to heating and above 1.2 kW allocated to DHW production. The heat demand was derived from the electricity demand required for space heating and DHW production. The natural gas demand was determined from the heat demand using the efficiencies of the condensing gas boiler, which can be found in Table 1 . Table 1 System input assumptions used in this study. The PV system, battery energy storage system, ground source heat pump system, condensing gas boiler, and remaining system input assumptions are separated by dashed lines. Replacement is denoted ad repl. and emissions for manufacturing and installation of system components is denoted as mfg emis.

Parameter
Value Unit Sources Statistics of the residential electricity consumption and the share of the GSHP functions are presented in Table A.1 . Also the estimated natural gas consumption statistics are shown in Table A.2 . The two years measurement period was scaled to a 30-years period by duplicating the 2 years period. We assumed no future change in electricity consumption from the dwellings. More details on the hourly, daily and seasonal fluctuation of the electricity consumption are presented in Figs. A.2 and A.3 . The winter months in the used time period (December, January and February) had higher temperatures than the long term average of 3.4 °C. The winter of 2013-2014 had an average winter temperature of 6 °C and the winter of 2014-2015 had an average temperature of 4.1 °C. The overall temperatures for 2014 and 2015 were respectively 11.7 °C and 10.9 °C compared to a long term-average of 10.1 °C [41] . Solar radiation during the used time periods was also a few percent higher [42] .

Model input assumptions and explanation
The model assumed input values to determine the technical, environmental and economic performance of the assessed systems. An overview of the model input assumptions for the various sys-tem components is shown in Table 1 . A detailed explanation of theses model assumptions is given in Appendix B .
The dwellings had no BESS installed, therefore we modelled the battery energy storage behaviour. An AC (alternating current) coupled lithium-based PV-battery storage system was assumed. AC coupled systems are very suitable to retrofit existing PV systems with storage, and are widely used in literature [46] . The behaviour of the BESS was modelled using a novel algorithm aimed for two applications simultaneously: increase in self-consumption and reduction of grid impact. Self-consumption was enhanced by storing surplus PV energy to use this on later moments. Grid impact was reduced by peak shaving of the imported and exported peak. The novelty of this algorithm is to predict the annual peak shaving potential by including the expected self-consumption enhancement. This potential depends on the used battery capacity and battery inverter rating.
This algorithm was written in Python (v3.5) and uses electricity production and demand profiles of the dwellings. The battery charge and discharge flows were simulated for 30 years using a time step of 5 min. This time step was used to obtain a good accuracy of battery storage capacity degradation. Battery degradation consists of cycle degradation and calendric degradation and was determined annually. The diminished storage capacity was subtracted from the original storage capacity for each year. After 15 years the battery storage is replaced and the capacity is set back to the full storage capacity. The used battery degradation model is explained in detail in a previous study [64] . Perfect forecasts of the PV energy production and the electricity consumption patterns were assumed. These were used to predict the moments to charge or discharge the battery storage system and consequently reduce the power impact on the local grid. An overview of the energy model algorithm with the description and explanation of the model steps is presented in Appendix C . Furthermore, the model is clarified using two days as example, graphical supported by

Calculation of self-sufficiency and grid impact
The self-sufficiency ratio was used to assess the impact of the PV self-consumption. The self-sufficiency ratio indicates the share of electricity demand ( P demand ) that is fulfilled by direct or indirect PV self-consumption. The direct self-consumption consists of PV produced power which is directly used in the dwelling ( P direct SC ). The indirect self-consumption contains power that is discharged from the battery and delivered to the dwelling ( P B discharge ). The self-consumed power is aggregated over the total lifetime of the system from the first time step of the first year (t = 1) until the last time step of the final year ( t end ), with a 5 min time step see Eq. (1) .
The specific self-sufficiency ratio ( SSR specific ) indicates the share of electricity consumption of each of the four GSHP functions that could be fulfilled by the direct or indirect self-consumption. The GSHP self-consumed power was divided into four temporal subsets (T F ), based on moments with power demand for each GSHP function. This self-consumed power was aggregated and divided by the electricity consumption of each function ( P demand, GSHP, t ) see Eq. (2) .
The impact of the grid was assessed using four parameters. The maximum import power peak and export power peak to the grid ( P G ) over the lifetime of the system was determined. This peak is a valuable indicator for distribution system operators to assess the impact of the residential system on their network. The other two parameters are the import peak to average ratio (IPAR) and export peak to average ratio (EPAR). The IPAR is defined as the ratio between the maximum import power peak and the average imported power from the grid. The export peak to average ratio EPAR is defined as the maximum exported power to the grid and the averaged exported power. These parameters are an indicator for the variability and magnitude of the power flows on the network, see Eq. (3) .

Calculation of investment attractiveness
The investment attractiveness of each system depends on the initial investment cost of a system ( I total ) and the annual cash flows (CFs) over the lifetime of the systems. The investment cost depends on the expenses for the heating system ( I HS ), the cost of the PV system( I PV ) and investment cost of the storage system ( I BESS ), see Eq. (4) . The cost of the PV and storage systems depends on the installed system capacities.
The total cash flow depends on the annual ( y ) energy cost from natural gas and electricity. The cash flow from natural gas ( CF NG ) depends on the imported natural gas ( NG import ) and the natural gas tariff ( π NG ). The cash flow of electricity ( CF Elec ) depends on the imported electricity ( Elec import ) and exported electricity ( Elec export ) and the consumption tariff ( π cons ) and the feed-in tariff ( π feed-in ).
The attractiveness of the investment in the system architectures was evaluated using the net present value (NPV), cost reduction, internal rate of return (IRR) and discounted payback period (DPBP), defined according Eq. (6) . The NPV gives a perspective of the current value by including the future risk and returns of the investment. The NPV provides an absolute number that can help in the decision process for the system selection of dwelling owners. A positive NPV represents a feasible investment, while a negative NPV shows an unwise investment. The future risk and return is included with the discount rate ( r ) and used to calculate the diminishing value of future returns over the economic lifetime ( L econ ) of the system. The NPV was calculated by subtracting the reference NPV ( NPV reference ) minus the actual NPV ( NPV actual ). The reference NPV included cash flows and investment cost from the reference case, and the actual NPV included the cash flow and investment cost of the analysed system architecture. The cost reduction ( C reduction ) is the relative reduction in NPV of the actual system compared cash flow to the reference system. The IRR gives an indication of the discount rate required to obtain a feasible investment. This was obtained by solving the discount rate with an NPV identical to zero. Thus the IRR shows the influence of the discount rate on the investment. The DPBP presents the time period until an economic investment is recovered. This was found by solving the year in which the NPV is identical to zero. System payback periods are frequently used in the residential PV and GSHP market, therefore DPBP was selected as an economic indicator. A maximum lifetime of 50 years was included in the mathematical solver. We included lifetimes longer than the economic lifetime of the systems to give a clear depiction of the DPBP diversity.

Calculation of avoided GHG life cycle emissions
The life cycle emissions of each system mainly depends on two emissions categories. The first category are emissions released during manufacturing and installation of the specific system components ( GHG mfg ). These emissions are from the heating system ( GHG mfg HS ), PV system ( GHG mfg PV ), and battery storage system ( GHG mfg BESS ). Emissions from PV systems manufactured in China were selected since more than half of PV module globally are currently produced in China [65] . The emissions from the heating system are per unit, emissions from PV per kWp capacity and BESS per kWh storage capacity, see Eq. (7) .
The second category contains emissions associated with the imported energy ( GHG E import ) and potential avoided emissions due to the exported energy ( GHG E export ). The emissions associated with imported energy consist of the emissions from electricity imported and emissions from natural gas combustion. Emissions from imported and exported electricity are calculated by multiplying the annual electricity import and export with the average annual emission factor of electricity (EFE). An annual linear reduction from the current emission factor towards zero for a certain future year was assumed. The emissions associated with exported electricity depend on the exported electricity multiplied by the annual emission factor, see Eq. (8) . Emissions from the natural gas consumption were included when a system contained a CGB to provide space heating and DHW. A emission factor of 50.93 kg CO 2 /GJ HHV was used to obtain the emissions from natural gas consumption.
A larger share of renewable PV production capacity could result in PV feed-in limitations. Moreover, exported PV produced by the dwellings could replace other renewable sources and these future marginal GHG emissions of the electricity systems are difficult to predict. Therefore, we assessed the avoided life cycle GHG emissions using two system boundaries: from an electricity system perspective ( GHG system ) and from a dwelling perspective ( GHG dwelling ). In the first perspective, all exported PV electricity is allocated to replace emissions from electricity generated by other sources. From the dwelling perspective, exported PV electricity does not account for additional avoided emissions, see Eq. (9) .
The avoided emissions ( GHG avoided ) are defined as the difference in emissions of the reference system ( GHG reference ) and the actual system ( GHG actual ) architecture. The reference system contains a CGB for heating demand that has similar standby consumption as the GSHP. Also the reference system included the emissions from import electricty and natural gas. The reduction is defined as the ratio between the avoided emissions and the reference emissions, see Eq. (10) .
The GHG payback period (GHGPBP) is an indication of the time period to offset the additional emissions from manufacturing and installation of the actual system ( GHG mfg actual ) compared to the reference system ( GHG mfg reference ). The avoided emissions excluding the emissions from manufacturing ( GHG avoided excl mfg ) from the reference system and the actual system were determined for each year. These were calculating according to Eq. (8) and Eq. (9) , excluding the emissions from manufacturing. The GHG payback period was found by selecting the time step for which the cumulative avoided emissions excluding manufacturing (CAE) are identical to the extra emissions due to manufacturing, according Eq. (11) CAE , y = Distribution of the self-sufficiency ratios for the four system configurations (a) and distributions of specific self-sufficiency ratios for each of the heat pump functions (b). The impact of power flows from the electricity grid is shown by the import peak distributions (c) and the import peak to average ration (d). The impact of power flows to the grid is presented by the export peak (e) and the export peak to average ratio (f). The left part of the violin plot shows the distribution of a PV system only and the right part of the distribution with a storage system. Distributions are shown for a 30-year lifetime. Mean values of distributions of the 16 dwellings are indicated by solid lines, and the 25% and 75% percentiles by dotted lines.

Results
This section provided the technical, economic and environmental impact of the GSHP systems combined with PV and storage. The results of each individual system was visualized using violin plots [66] . This type of graphical illustration gives a quick indication of the distribution of results obtained from each dwelling and system. The results are presented for normalized PV system capacities of 1 kWp for each MWh of annual electricity consumption and battery storage capacities of 1 kWh per MWh annual consumption. Technical, economic and environmental impact of smaller and larger PV system sizes and battery storage capacities are presented in Appendix D .

Self-sufficiency ratio and grid impact
SSR distributions of dwellings and their specific GSHP functions are visualized using violin plots in Fig. 2 (a) and (b). PV systems with GSHP show an average SSR of 25%, with a range between 21% and 31%. Adding BESSs increases this SSR on average with 18% point, ranging from 39% until 49%. Dwellings without GSHP show a higher SSR of 31% on average for PV only and 49% with BESS. Clearly, adding a GSHP reduces SSR. This is simply due to the fact that PV production occurs mainly in the summer months, whereas most GSHP demand consumption occurs during colder winter months (see Figs. A. 2

and A.3 ).
Each of the specific heat pump functions shows a different contribution to the SSR ( Fig. 2 (b)). Cooling demand shows the largest specific SSR of which 53% directly is provided by the PV system. This number can be increased to 94% when a storage system is added. The SSR of space heating is limited to 9% in case of PV only and 14% when a storage system is used, since most space heating occurs during colder winter months. DHW is temporarily stored in a tank and therefore DHW production has a demand side management potential. The average SSR of DHW increased by 12% with storage and by 23% with twice the storage capacity, see Based on these numbers, we estimated that a quarter of DHW production has the potential to be shifted to moments with excess PV electricity. This will increase the direct self-consumption of DHW production and reduce the storage requirements for batteries. In total, 19% of the GSHP consumption can be provided by direct PV electricity and an additional 10% by BESS.
The influences of the systems on the electricity grid are presented in Fig. 2 (c)-(f). Dwellings with GSHP have an average 0.8 kW higher import peak than dwellings without GSHP, however the maximum peak demand of the GSHP system is 1.8 kW. Hence, power peaks from GSHP are not occurring at the same moment as power peaks of remaining dwelling appliances. Deploying BESS systems leads to an average reduction in import peak of 45% for dwellings with GSHP and of 30% for dwellings without GSHP. Dwellings with GSHP systems show an average reduction of import peak to average ratio from 15 to 11 when BESS are used. The reduction of exported peak PV power is around 38% for both dwellings with and without GSHP. Export peak to average ratio are reduced from 3.3 to 2.1 when storage was included in the system. Larger storage capacities lead to larger reduction of the import and export peaks, see Fig D.2 . Also shown here is that on average 0.09% of the total charged energy was pre-charged to reduce peak demand and only 0.04% of discharged energy was pre-discharged to store future PV peaks. Furthermore, SSR are 0.07% lower due to the peak shaving application. Moreover Fig D.2 shows the impact of three potential PV feed-in limitations. These results show that BESS have a large potential in enabling more GSHP and PV systems on existing low voltage grids.

Investment attractiveness
The investment attractiveness for 16 dwellings and 5 systems architectures is presented in Fig. 3 . A large distribution range in net present value is shown for dwellings with a GSHP, highly influenced by the dwellings heat consumption. The GSHP system shows a net present value range from € −1 . 6 · 10 3 to 6 · 10 3 , with a mean of ≈ € 1100. 7 of the 16 dwellings with solely a GSHP system have a negative NPV. A PV system combined with GSHP results in a positive average NPV of € 5.7 · 10 3 . A BESS reduces the average NPV around € 1100, showing that investing in storage is not profitable.
The average cost reduction of a PV system with GSHP is larger than for only a PV system. However, highest internal rate of returns are shown for PV systems only, with an average of 7.5%. The investment in a GSHP system is two to three times larger than the investments in a PV system. Consequently, the internal rate of return is higher for solely a PV system than for a PV system with GSHP. The GSHP systems show an average IRR of 2.6% and GSHP systems with PV an average of and 4.2%. Discounted payback period of GSHP systems are between 18 to 37 years, whereas solely PV systems show an average of 12 years. Distribution range for DPBP by only GSHP is substantially larger than the distributions seen from only a PV system. The revenue of GSHP is caused by the avoidance of natural gas consumption and is highly dependent on the occupancy and behaviour of the residents as well as the dwelling properties. The DPBP decreases with larger PV system or BESS capacities (see Fig. D.3 .)

Avoided life cycle GHG emissions
Distributions of avoided life cycle GHG emissions, GHG emission reduction and GHG payback periods are presented in Fig. 4 . Avoided emissions from an electricity system perspective show a relative large distribution spread, especially for systems with GSHP. Avoided emissions of GSHP systems are between 30.6 and 59.2, with an average of 42.4 tCO 2 -eq for the 16 dwellings. Systems with PV have an average of 72.7 tCO 2 eq, demonstrating that PV systems contributes with 30.3 tCO 2 eq, on average of avoided emissions. Avoided emissions are larger for dwellings which include a GSHP due to the larger PV capacities. BESS systems lower the avoided emissions from an electricity system perspective, mainly caused by charge and discharge losses. The higher avoided emissions by storage are caused by the higher SSR of 10%. Larger PV system capacities have more avoided emissions from an electricity system perspective, see Fig. D.4 . Emissions are reduced by 47% using a GSHP and 80% with the use of a PV system. GHGPBP are lowest for a GSHP system only with 2.5 years of time required to offset the invested emissions from manufacturing and installation. All system configurations have an average GHGPBP below 6 years from a system perspective.
Avoided emissions were also determined from a dwelling perspective, which only includes emissions avoided within a dwelling. This perspective is currently not realistic since surplus PV electricity can be exported to the grid without limit, and is shown purely as a theoretical indicator. Avoided emissions show an average of 44.3 tCO 2 -eq for GSHP with PV, and an additional 7.3 t when storage is included. The higher avoided emissions by storage are caused by the higher SSR of 10%. Larger PV system capacities have more avoided emissions from an electricity system perspective, but less avoided emissions from a dwelling perspective, see

Techno-economic correlations with avoided emissions
Observable correlations between technological or economic parameters with the avoided life cycle GHG emissions parameters are given in Fig. 5 . The contribution of the PV system to the avoided life cycle GHG emissions from a dwelling perspective is shown using the self-sufficiency ratio. A clear linear trend is observed between the SSR and the avoided GHG emissions by a PV system. A minimal SSR of 22% is required to obtain net positive emissions. For SSR > 22%, avoided emissions from a dwelling perspective increase with ≈ 80 gCO 2 -eq per Wp for each percentage point SSR.
Another clear correlation is shown between the net present value and the total avoided emissions from a system perspective. The NPV increases with ≈ € 275 for each additional ton of avoided tCO 2 -eq. Around 40 tCO 2 -eq should be avoided to obtain benefits from emissions reduction. Dwellings with a GSHP system plus PV and storage show an NPV increase of around € 230 for each additional ton of avoided tCO 2 -eq. Complementing these dwellings with storage shows a comparable slope. Dwellings with PV show a NPV increase of ≈ € 180 for each ton of avoided emissions. This is decreased to € 170 when storage is added. These numbers show that avoiding emissions with PV systems obtain lower benefits as avoiding emissions using GSHP systems per tCO 2 -eq. However, PV systems show only positive NPV values, whereas GSHP system needs a minimum of avoided emissions to obtain a positive NPV number.
Discounted payback periods are decreasing with avoided emissions. Dwellings with GSHP systems only show a decrease of 0.7 year for each tonne. The dwellings with the largest avoided emissions of 58 tCO 2 -eq show the lowest DPBPs of 18 years. GSHP systems with PV and storage have a DPBP reduction of ≈ 0.25 year per tCO 2 -eq. Dwellings with PV only show a decrease in DPBP of 0.03 year per ton tCO 2 -eq, or 12 days.

Discounted payback period sensitivity
A sensitivity analysis of five input parameters on the DPBPs of GSHP systems and GSHP with PV and with PV + storage is presented in Fig. 6 . Dwellings with only a GSHP system need an increase of 1.2% per year in natural gas tariff to obtain DPBP < 30 years for all dwellings. Higher electricity tariffs lead to higher DPBP and are even sharply increase when these changes are above 2%. Yet, the impact of higher electricity tariffs is greatly reduced when a PV system is added to the dwelling. The value of selfconsumed energy increases with higher electricity tariffs, therefore PV systems obtain more revenue. This is even more reduced when a battery storage system is added due to the increased selfconsumption. Consequently, the DPBP of GSHP with PV and BESS are only slightly increasing with higher electricity tariffs.
A CO 2 taxation on the used electricity and gas consumption of around € 70 per tCO 2 reduces the DPBP < 30 year for all dwellings. An additional 11% cost reduction results in similar DPBP. This is a realistic cost reduction for the future GSHP mass market scenario [56] . Yet, an increase in the discount rate, from 2% to 3%, results in an average DPBP to 33 years. In general, the PV and PV battery systems decrease the DPBPs of the overall system. Two extended sensitivity analyse on the DPBPs are given in Appendix E . First, the combined impact of natural gas tariffs and electricity tariffs is presented. For example, it is shown that for dwellings with a GSHP, a 1% increase in both electricity tariff and natural gas tariff will keep the average DPBP between 25 and 30 years. Furthermore, lower electricity tariff are positive for GSHP payback periods but negative for PV and battery storage payback periods. Consequently, GSHP combined with PV and storage level out the influences of higher or lower electricity tariffs. Second, the combined impact of a change in GSHP investment subsidies with a change of natural gas tariff, electricity tariff, CO 2 taxation or discount rate is shown. This shows that a combination of GSHP investment subsidy with a decreasing natural gas tariff, or increasing CO 2 taxation shows a steep reduction in DPBP. For example, an additional 10% additional subsidy and a € 65/ tonne CO 2 will reduce the average DPBP for GSHP systems to 20 years.

Avoided life cycle GHG emissions sensitivity
Avoided system emission depends on the carbon intensity of electricity from the grid. This intensity is expected to decrease in the next decades due to a larger share of renewable electricity generation capacity. Yet, the gradient of this decrease depends on policy and technological development. We assessed the impact of this reduction by assuming a linear reduction of 2016 GHG emissions to zero emissions for a given year, presented in Fig. 7 . A faster decrease in emission intensity results in lower avoided emissions by PV and BESS, but more emissions are avoided using a GSHP. From a dwelling perspective, the emissions due to manufacturing of solely a PV system are not recouped before 2040. A sixth scenario was added which presents the impact of a PV system made in the EU. In this case, the emissions due to PV system manufacturing are already recovered by 2029, due to the lower emissions during production of PV in the EU compared to China [16] .

Discussion
This research assessed the techno-economic and environmental impact of GSHP systems combined with PV and storage. We showed that PV systems contribute to one fifth of the total GSHP electrify demand. Battery energy storage reduces the GSHP impact on the local grid. Also GSHP systems can greatly reduce the life cycle GHG emissions of dwellings by replacing natural gas-fired boilers.

Comparison with previous studies
PV self-sufficiency ratios are comparable as found in previous studies [23,24,69] . Thus, the models used in these studies are representative for real applications. We found that investment attractiveness of systems is highly dependent on the total avoided life cycle GHG emissions. An NPV increase of € 275 per tCO 2 -eq was observed for GSHP systems, which are in a similar range as found previously [9] . At least 40 tCO 2 -eq of life cycle GHG emission should be avoided to obtain a positive NPV for these systems. With a PV and storage system, the NPV decreases to € 230 per tCO 2 -eq. No comparable studies were found which assessed this combined influence.
Results on the environmental impact are comparable with previous studies that investigated avoided emissions using GSHP systems from German, Belgium or the UK. These countries have similar climate and dwelling conditions as the Netherlands. We found that GSHP systems reduce life cycle GHG emissions between 18  and 56% compared to a natural gas boiler. This corresponds to annual avoided emissions ranging between 1 and 2 tCO 2 -eq for GSHP systems. These results are comparable with previous studies that investigated savings for Germany, Belgium or UK [6,7,9] . If PV is  [67] . Also, we assumed a 10% larger PV production for France and kept remaining model parameters similar [68] . added to the systems, than an additional 1 ton of avoided CO 2eq emissions can be achieved, similar as found before for Belgium [32] .
We did not find studies that assessed emission reductions from GSHP systems combined with PV and storage. A large share of heat demand in dwellings is provided by fossil fuels in France, Germany and United Kingdom. This heat demand can be fulfilled by GSHP systems [1] . Therefore we provide a first indication of this potential for these countries, shown in Fig. 8 . This shows that United Kingdom and France would have respectively 22% and 45% more avoided emissions for GSHP systems respectively, compared to the Netherlands.
Yet the impact of solely PV is lower, with even negative emissions for France due to the large share of nuclear electricity generation. However, if we would assume that PV systems were produced in Europe, then the PV system footprint per kWh would be almost similar as the emission factor of France, see Fig. A.1 . Consequently, negative emissions would be reduced to −1.1 t of CO 2 -eq. Moreover, Germany has a higher electricity tariff of ≈ 0.30 € for each kWh [70] . This results in higher electricity costs for GSHP but also larger revenues from PV and BESS systems. We analysed the impact of this tariff and found that a solely GSHP system would be far from profitable. However, GSHP with PV and storage would have a lower DPBP than presented in our research.

Implementation challenges
Dwellings used in this study use low temperature floor heating, resulting in highly efficient GSHP operation and thus avoid more emissions compared to using a conventional dwelling heating system. Existing dwellings must be renovated and insulated to a certain level to effectively use a GSHP system [71] . Sufficient ground area and proper subsurface characteristics should be available for the ground heat source [72] . As a result, for densely populated urban areas, other heating systems could be more beneficial from an environmental and economic perspective. For example, heating using solar water heating, bio gas combustion or district heating systems are excellent alternatives [73] . Moreover, sufficient roof area should be available to install the PV system. The assumed PV systems sizes for the dwellings with GSHP are between 3.5 and 7.8 kWp. Assuming a 0.2 Wp/m 2 capacity factor, the required roof area are between 17.5 m 2 and 39 m 2 . This could be a limitation in densely populated areas. Nevertheless, 2 million dwellings in the Netherlands are detached or semi-detached, that have a high prospective to install GSHP with sufficient roof space for a PV system. The annually avoided emissions for all these dwellings would be 4.8 Mt, reducing the total annual CO 2 emissions with 2.9% in the Netherlands [3] .
Other cost reduction options excluded in this research are also available. The decrease of power flows by BESS could potentially reduce grid connection cost, yet a flexible or dynamic capacity tariff structure should be available to obtain monetary benefits [74] . In addition, BESS revenues can be increased by provision of energy arbitrage or frequency control restoration [75] . Also, collaborative planning and investing in PV, storage and GSHP systems is highly recommended to decrease cost.

Limitations and further research
Our research has several limitations that could affect the findings. We used data of 2 years thereby including the relative colder winter period of the beginning of 2015. However, future winters are expected to be warmer due to global warming [76] . This would lower heating demand, and therefore slightly reduce the investment attractiveness of GSHP systems. Moreover, annual electricity prices and emission factors were assumed, but future residential energy tariffs could change every 15 min due to an increased share of variable renewable energy generation.
It is expected that emission factors from power generation will have a higher variability due to the larger share of renewables. Avoided emission of PV systems could be lower when power of these systems will not replace fossil fuels. Consequently, storage is required to keep avoiding emissions, and could even increase avoided emissions. Especially, when battery storage systems could discharge energy to the dwelling on moments when a large share of fossil fuel fired power plants are in the power generation mix. Therefore new battery control strategies should be developed that include the marginal emissions factors. Future research should focus on the role of marginal emissions factors and dynamic tariffs.

Conclusion
This study assessed the technical, economical and avoided life cycle GHG emissions of GSHP systems with PV and battery energy storage. We used measured data of 16 dwellings and assessed the performance over a 30-year lifetime.
Dwellings with a PV system show a SSR of around 31%, while the use of a GSHP system reduces this to 25%. PV can supply 19% of the total GSHP demand, while this can be increased to 29% using battery energy storage. Moreover, storage reduces the peak demand with 45% of dwellings with a GSHP, which enables more GSHP systems on the low voltage grid. The investment attractiveness for dwellings with a GSHP shows DPBPs between 19 and 34 years. PV and storage have significant lower DPBPs helping to reduce the overall systems DPBP.
Avoided life cycle GHG emissions from GSHP systems are between 31 and 59 tCO 2 -eq for the 16 dwellings. The NPV increases with € 275 per tonne of avoided emissions. Adding a PV system to the dwelling increases the average avoided life cycle GHG emissions with 30 tCO 2 -eq. Also, this lowers the NPV increase per tonne of avoided emissions to € 230. Battery energy storage only avoids emissions from a dwelling perspective. Therefore, storage is currently not recommended to reduce emissions since all PV electricity can be exported.
The sensitivity study provides recommendations to improve the investment attractiveness of GSHP systems. Policies should focus on higher natural gas tariffs or include CO 2 taxation to encourage less natural gas use and switch to GSHP systems. Also, additional GSHP investment subsidies combined with these policies show promising results to obtain economic feasibility investments in GSHP system.

Acknowledgement
This work is part of the research programme Transitioning to a More Sustainable Energy System (Grant No. 022.004.023), which is financed by the Netherlands Organisation for Scientific Research (NWO). We are grateful to Daphne Geelen from Enexis for providing data.

Appendix A. Input time series statistics
Technical and environmental performance statistics of the PV systems over the assessed lifetime are shown in Fig A.1 . Used PV system capacities for dwellings without a GSHP have an average capacity of 3.3 kWp, varying between 1.7 and 4.9 kWp. Dwellings with GSHP have an average capacity of 5.9 kWp. The annual specific yield of the systems ranges between 791 and 886 kWh per kWp, with an average of 853 kWh per kWp. Note that this performance includes a 0.5% annual degradation. The PV systems show a good performance which is comparable with other Dutch PV systems [68] . GHG emissions are on average 67 gCO 2 -eq per kWh for PV modules manufactured in China and 37 gCO 2 -eq per kWh for PV modules manufactured in Europe.
Statistics of residential electricity consumption and contribution of the GSHP functions are given in Table A.1 . Annual average electricity consumption without GSHP is 3139 kWh, which corresponds to the average consumption of Dutch households [77] . The electricity consumption shows a relative large range, between 1557 kWh and 4727 kWh. This is mainly influenced by the household appliances and composition. Dwellings with GSHP show an average electricity consumption of 5920 kWh. The contribution of a GSHP shows a more narrow distribution, from 1933 kWh to 3760 kWh.  GSHPs contribute to an average of 48% of the total electricity demand, ranging between 35% and 60%. The largest share of the GSHP consumption is used for space heating with an average of 27%. Cooling and standby consumption shows smallest shares of respectively 1.7 and 2.6%. The GSHP consumption is mainly influenced by space heating consumption, which has the largest range between a minimum share of 17.1% and a maximum share of 41.1%.
The estimated natural gas consumption for DHW production and space heating is shown in Table A.2 . 17% of the natural gas consumption is used for DHW production and 83% for space heat-

Table A1
Annual electricity consumption statistics of the 16 dwellings included in this study. The absolute values are given in the left columns and the shares of the total is given in the right columns.

Average
Max ing. The total annual natural gas consumption is between 30 and 58 GJ, with an average of 42 GJ. This corresponds to an annual natural gas demand of 1490 m 3 , which agrees with the average natural gas consumption of Dutch households [77] . The average stacked daily energy consumption of the 16 residential dwellings is shown in Fig A.2 . The conventional electricity consumption is relative constant with an average of ≈ 9 kWh. The energy consumption for space heating shows a high volatility, related to the outdoor temperature. Cold winter days have the most significant impact on the variation in daily electricity consumption of the dwellings. Besides, the winter months of 2015 required more space heating than the winter months of 2014, due to lower temperatures. This shows the relevance of using multiple year consumption data for this type of research.
Distributions of hour of the day and monthly electricity demand from the 16 dwellings and from the specific GHSP functions are shown using a violin plot in Fig. A.3 . Conventional electricity consumption shows a clear trend with higher consumption during daytime and peaks between 5 until 6 pm. The consumption of the GSHP shows a clear increase of consumption over all hours. The average standby consumption shows a small increase during daytime, indicating that the GSHP is less used for other functions during these moments. The cooling demand shows a decrease during daytime, and increases during night-time when residents are sleeping. Also, the dwellings are highly insulated and therefore the heating of the dwellings by the environments encounters a time lag. DHW and space heating have an order of magnitude higher demand than the standby and cooling consumptions. DHW electricity demand shows a low consumption during the night, and peaks in the morning and evening hours, mainly caused by showering. Space heating shows high demand during the night and decreases and low during daytime. The distribution range of electricity consumption for space heating and DHW is highly dependent on the dwelling and its residents. The conventional consumption is quite constant over the months, whereas the heat pump consumption shows a seasonal effect. Higher consum ption is seen by the cooling demand in the summer months, whereas space heating shows peaks in winter months. Standby consumption shows a small increase in the summer month, caused by the lower utilization of the heat pump system for other functions. The DHW demand shows a small increase of demand during winter months. The majority of electricity demand is caused by space heating which increases the overall demand of GSHP system in winter months.

Appendix B. Explanation model input assumptions
This section provided an additional explanation of the input parameters used, also shown in Table 1 .

B.1. PV system assumptions
Currently, net-metering policy is in place for residential PV systems. This policy enables system owners to offset the cost of electricity from the utility with their own generated PV electricity. Consequently, most system owners installed a PV system capacity that is able to cover the annual electricity consumption. Therefore, the installed PV system size was set to 1 kWp for each MWh of annual electricity consumption [21] . The electricity production of the PV system was reduced with 0.5% per year by PV system degradation [43] . PV system costs of 1500 € /kWp were assumed, including all components and installations [44] . In the Netherlands, value added tax (21%) of the PV system investment can be reclaimed, which results in a net investment cost of 1200 € /kWp. Annual maintenance cost (O&M) were set to 1% of the PV system investment costs [45] . The PV inverter is replaced after 15 years [46] . Re-placement cost of the PV inverter are 100 € /kWp [44] . Emissions from manufacturing PV systems with multi crystalline silicon modules were expected to be 1590 gCO 2 -eq for each Wp when made in China, and 824 g CO 2 -eq for each Wp when made in Europe [47] . Emissions from PV produced in Europe are significantly lower due to lower emission intensity of the European (ENTSO-E) electricity mix.

B.2. Battery energy storage system assumptions
The battery energy storage capacity was set to 1 kWh for each MWh of annual electricity consumption, based on previous research on optimal storage system designs for residential dwellings [46] . A commonly installed battery inverter rating of 0.5 kW per kWh of storage capacity was used [48] . Thus, two hours are required to completely charge the battery capacity from 0% until 100%. The battery state of charge (SOC) range was set between 0% and 100% of the battery storage capacity, so the full capacity potential could be assessed. A constant battery (direct current) DC-DC efficiency of 96% was assumed for battery charging and discharging, almost similar as the commercial available Tesla Powerwall [48] . The battery inverter efficiency curve from SMA Sunny Boy Storage was used to model the AC-DC and DC-AC conversion [78] . A BESS standby consumption of 0.1% of the rated inverter power was assumed. A battery cycle lifetime of 50 0 0 full equivalent cycles and a calendric lifetime of 15 years were used based on a previous study [49] . We expected battery storage pack cost of 300 € /kWh and storage balance of system (BOS) and EPC (Engineering, Procurement, and Construction) costs of 200 € /kW [50] . Annual maintenance cost (O&M) were set to 1% of the investment cost of battery energy storage system [45] . The battery inverter and battery cells are replaced after 15 years [46] . Replacement cost of the battery inverter is assumed 100 € /kW [44] . Replacement cost for the battery energy storage pack are 100 € /kWh, based on learning curves [23] . For the production of Li-Ion battery energy storage systems 110 gCO 2 -eq for each Wh of Li-Ion storage was selected [53] . 124 gCO 2 -eq was used for the production for each W of battery inverter, based on production requirements for a PV inverter [47] .

B.3. Ground source heat pump system assumptions
A constant coefficient of performances was assumed based on the relative constant temperatures of the ground sources. A COP for space heating of 5.7 and a COP of 1.9 were assumed based on the technical specification of the heat pump [40] . The largest cost component of a 4 kW GSHP is the drilling and installation of the vertical ground loop, which is estimated at € 60 0 0 [54] . The EPC cost is expected to be € 40 0 0. The cost for the heat pump are estimated € 10 0 0 for each kW of nominal thermal output, thus euro 40 0 0 for a 4 kW thermal unit [9] . This results in a total of a GSHP of € 14,0 0 0 for the 4 kW unit. Heat pump systems in The Netherlands are subsidized depending on the thermal power and system type [12] . The investment subsidy for this system is € 2800, which is similar to an investment subsidy of 20% of the total GSHP investment. Annual maintenance costs for GSHP are estimated 50% lower than the maintenance cost of CGB systems, specifically 50 € [55] . The lifetime of the heat pump is expected 20 years [9,55] . We expect that future heat pump costs will be reduced with 10% and that no subsidies can be used for the replacement of a heat pump. This results in heat pump replacement costs of 3600 € per unit [56] . Emissions for manufacturing, transportation and installation of the 4 kW GSHP are estimated to be 1760 kg of CO 2 -eq per unit [7] .

B.4. Condensing gas boiler assumptions
The efficiency of the condensing gas boiler depends on the required temperature of the heated water. These requirements are set similar as the temperatures delivered by the ground source heat pump system, specifically 28 °C for space heating and 60 °C for domestic hot water. Corresponding efficiency are 95% for space heating and 85% for domestic hot water [32] . The efficiencies are based on a higher heating value of Dutch natural gas of 35.17 MJ/Nm 3 [79] . The cost of the reference residential CGB are € 1500 [4] . CGB have a lifetime of 15 years and a similar replacement cost of € 1500 was selected. The CGB requires annual services cost of € 100 which is estimated based on the average service contract in the Netherlands for 2015 [57] .

B.5. Remaining systems assumptions
Costs of replacing components that are not mentioned above were included in the maintenance costs, as well as the EPC cost of replacing components. We assumed no salvage value for the system components. Also, we expect that emissions from manufacturing are 25% lower for the replaced components, except for the GSHP system. For this component we assumed 440 gCO 2 -eq per unit since only the heat pump is replaced.
A discount rate of 2% was selected based on the currently low interest rates [59] . The consumption tariff was set to 0.178 € /kWh, based on the average household retail prices for electricity from 2014 until 2016 [60] . The net metering policy in the Netherlands is currently in place and probably will be replaced with a feed-in premium by 2020. This feed-in premium will be gradually reduced from 0.12 € /kWh in 2020 to 0 € /kWh by 2036 [61] . Next to the feed-in premium, residential consumers will receive a fee based on the wholesale electricity price. We assumed that this is similar as the expected wholesale electricity price. A linear increase from 0.032 € /kWh in 2020 until 0.044 € /kWh in 2030 was assumed, based on the wholesale market prices projections for 2030 in the Netherlands [80] . From 2030 onwards, we assumed a 0.5% increase for this wholesale electricity price. The natural gas consumption price was set to 21.84 € /GJ based on the average household retail prices for natural gas from 2014 until 2016 [60] . An 0.5% increase in gas and electricity tariffs was selected due to additional energy taxes for supporting investments of renewable energy generation [62] . An average CO 2 emissions intensity factor of 490 gCO 2 -eq per kWh for 2016 was used [15] . A linear reduction of CO 2 emissions to zero in the year 2050 was assumed based on the Dutch energy agreement for sustainable growth [63] . The emission factor for natural gas is 56.6 kg CO 2 /GJ using a low heating value of 31.65 MJ/Nm 3 [81] . In this research we used the high heating value, thus the emission factor was converted to correspond with the HHV. This results in an emission factor of 50.93 kg CO 2 /GJ.

Appendix C. Detailed energy model algorithm explanation
The used energy model algorithm is explained with a schematic overview, shown by Fig. C.1 . Before the model algorithm is executed, the power export limit ( P EL ) and power import limit ( P IL ) were defined. These are pre-defined limitations for exporting and  importing power respectively to and from the grid. These limits depend on the maximum peak shaving capacity of the battery storage systems. This depends on the used battery inverter rating and the battery storage capacity. The inverter rating determines the maximum height of the peak that can be shaved. The storage capacity determines the maximum duration of the peak that can be shaved. The PV peak power and peak demand that could be stored by BESS was calculated for each day of a year. This calculation included pre-charging of energy to the battery to shave the demand power peak on a later moment. Also pre-discharging of energy to shave the PV production power peak was included. The power export and import limits were determined from the maximum PV peak and peak demand that could be stored.
Each model time step starts with the assessment or surplus of PV power is available or power demand from the building is requested, see grey rectangle with dashed lines. The difference between PV power and electricity demand is given by ( P pot ). If more PV electricity is produced than consumed in the dwelling, then the power difference will be compared with the power export limit. When the power difference is larger than the export limit, then the potential SOC for the next time step ( E B pot ) is determined. This depends on the current and expected charge potential for the next time step. The potential battery SOC is compared with the battery SOC reserved to charge the PV peak power that is larger than the export limit ( E B REL ). The reserved amount of charged or discharged energy is limited by the maximum battery SOC ( E B max ) and minimum SOC ( E B min ), as well by the battery inverter rating.
If the potential SOC is higher than the reserved SOC for PV peak charging, then the battery will only charge the PV peak dif-ference. Else, the battery inverter will charge all surplus PV electricity. If P pot is smaller than ( P EL ), then the battery SOC reserved for PV peak charging is compared with the current SOC ( E B cur ). When the current SOC is larger, then the battery is discharged until the current SOC is similar as the reserved SOC for PV peak charging.
When more electricity is consumed than produced by the PV system, then P pot is compared to the power import limit ( P IL ). When this limit is smaller, then the potential SOC for the next time step is determined, which in this case depends on the current and expected discharge potential. The potential SOC is compared with the battery SOC reserved to discharge the peak demand ( E B RIL ). If the potential SOC is lower than the revered SOC, then the battery will only discharge the peak demand. Else, the battery inverter aims to discharge all required electricity. If P pot is smaller than P IL , then the battery SOC reserved for load peak discharging is compared with the current SOC E B cur . When the current SOC is below the reserved SOC level, then the battery is pre-charged until the reserved SOC is reached. If the current SOC is larger than the reserved SOC level, then the model assesses or pre-discharging is required to charge PV peaks on later moments. This is conducted when the solar elevation ( G ELV ) is > 0 and if the current SOC is larger than the revered SOC for PV peak charging E B REL . The model time step ends with an action to charge, discharge or do nothing. Afterwards the algorithm continues to the next time-step, back to the assessment of P pot .
The model behaviour for a summer and a winter day is graphically explained in Fig. C.2 . The 13th of June has a larger forecasted PV production than a forecasted electricity consumption. Consequently, battery storage capacity is reserved to charge the excess PV power peak, visualized by the reserved capacity to reduce the exported peak power ( E B REL ). To obtain this reserved capacity, a small amount of storage is discharged to the grid in the early morning around 6.00. Consequently PV peak production is used for charging and stored in the battery. The battery is fully charged at 15.00, and discharged later in the evening and night. The 14th of June has less PV production forecasted, therefore a lower amount of capacity is reserved for PV peak shaving. No electricity from the grid was used during these two summer days.
The winter days show a different consumption pattern with significantly higher electricity consumption by the GSHP, mainly used for space heating. A demand power peak is forecasted around 12.00 on 14th of December. Consequently, battery energy storage is required to supply this demand peak. The battery storage is precharged on 13th of December to provide the peak demand of the next day. The small amounts of PV produced electricity are directly used or stored.

D.1. Self-sufficiency impact
Self-sufficiency ratios and specific self-sufficiency ratios for smaller and larger PV systems and storage capacities are shown in Fig. D.1 . Self-sufficiency ratios are increasing with larger PV system capacities. A higher SSR is observed for dwellings without a GSHP than dwellings with a GSHP. The influence of the individual demand pattern on the SSR decreases with larger PV system capacities. Furthermore, the distribution range increases due to the larger influence of the demand patterns on the SSR. Battery energy storage results a larger increase in SSR for bigger PV systems than for smaller PV systems. A 1 kWh storage system increases the SSR with 14% for a 0.5 kWp PV system with GSHP. This increase is 20% in SSR for a 1.5 kWp PV system with GSHP. Larger storage sizes, of 2 kWh per MWh of annual demand, increase only the SSR of PV systems with a capacity of 1 or 1.5 kWp.
The specific self-sufficiency ratio of space heating is relatively low compared to the other GSHP consumption components. PV production occurs mainly in the summer months, whereas most space heating demand request occurs in the winter months ( Fig. A.3 ). A 1 kWh storage system does not show a slightly improvement of the SSR. Specific SSRs for domestic hot water are significant larger than space heating and go up to 58% with a 1.5 kWp PV system and a 2 kWh storage system. DHW production occurs largely during evening hours, which makes storage ideal to improve SSR. Cooling demand occurs in the summer months, and thus shows the highest specific SSR. Also, cooling reaches the highest specific SSR when storage is deployed.

D.2. Grid impact
The impact of smaller and larger storage capacities on the maximum import peak, maximum export peak, import peak to average ratio and export peak to average ratio are shown Fig. D.2 . Also three other indicators are presented. The pre-charged percentage shows the ratio of energy that is charged in the battery used for demand peak provision to the total charged energy. The predischarged percentage shows the ratio of energy used to discharge the battery for PV peak storing to the total discharged energy. The curtailment loss ratio (CLR) is an indicator of the share of PV en-  Distribution of the influence of battery storage capacities on lowering demand power peaks from the grid and reducing PV peaks to the grid. The reduction of the demand peak on the grid is given by the import peak (a), import peak to average ratio (b), pre-charged ratio (c). The reduction of the PV system peaks on the grid is shown using the export peak (d), export peak to average ratio (e), pre-discharged ratio (f). The impact of storage on the curtailment loss ratios for potential feed-in limitations are shown for a 0.1 kW/kWp feed-in limit (g), a 0.2 kW/kWp feed-in limit (h) and a 0.5 kW/kWp feed-in limit (i). Distributions are shown for dwellings without and with a heat pump. A normalized PV system size of 1 kWp pv for each MWh demand was selected. Note that the battery storage capacities are not equally dispersed. ergy which cannot be exported due to a potential future feed-in limit. This is the PV energy above a certain feed-in limit divided by the total produced PV energy. The PV feed-in limit is given in kW per kWp of installed PV capacity.
With a 2 kWh storage capacity, the import peak for dwellings without GSHP is reduced to an average of 2.8 kW and with a GSHP to 3.1 kW. The reduction using a 1 kWh storage capacity is larger for dwellings with GSHP than dwellings without, yet the opposite is shown for smaller and larger storage sizes. For smaller storage sizes, the absolute storage size for systems with a GSHP are larger, hence more peak demand can be shaved. For larger storage sizes, the peak demand reduction of dwellings with a GSHP is already reaching a limit. This limit is due to the width of the peak, which is increasing significantly when the peak is further reduced. IPAR for dwellings without a GSHP is greatly reduced with storage, from 15.0 with no storage to 8.1 with 2 kWh of storage. Dwellings with a GSHP show a smaller reduction, from 8.1 to 3.8. The pre-charged distribution shows an average value of around 0.3% for a 2 kWh storage capacity.
The export peaks are reduces in a similar way for both dwellings without and with a GSHP. A remarkable increase in EPAR from a 1 kWh battery to a 2 kWh battery is observed for dwellings without GSHP. Larger storage capacities reduce the peak PV power, but also reduce the average exported power, since more electricity is used for self-consumption. If the average power is reduced more than the peak power, then a higher EPAR occurs. A significantly higher percentage of a 2 kWh battery storage capacity is pre-discharge for dwellings with GSHP than dwellings without GSHP. Dwellings with GSHP have a relative higher storage capacity since their electricity consumption is larger. Consequently, more energy is still stored in the battery before the next morning, es-pecially in the summer months with relative low night consumption from the GSHP. As a result, more electricity is pre-discharged which is required to obtain empty storage capacity for PV peak charging.
The curtailment loss ratios are shown for three potential PV feed-in limitations. CLR are higher for dwellings with a GSHP then systems without. The most restricted feed-in limit (0.1 kW/kWp) shows an average PV energy loss of 40% without GSHP and 45% with a GSHP. With a 2 kWh storage capacity, this can be reduced to 17% and 22% for respectively dwellings without and with a GSHP. Moreover, with a feed-in limit of 0.5 kW per kWp, a 2 kWh storage capacity can almost completely abolish the potential curtailment losses.

D.3. Impact on discounted payback periods
The impact of smaller and larger PV system and battery storage capacities on the discounted payback periods is presented in Fig. D.3 . DPBPs for dwellings without GSHP are in increasing with PV system size when no storage is used. Dwellings with GSHP show a slight reduction in average DPBP with an increase of PV capacity. The results also show that storage systems are not profitable under all scenarios. If a storage system is added, then it should be designed based on the installed PV capacity. If the PV system size is relative small, then an oversized storage capacity will greatly increase the DPBP.

D.4. Impact on avoided life cycle GHG emissions
The avoided life cycle GHG emissions from an electricity system perspective and dwelling perspective for smaller and larger PV  system and storage capacities are presented in Fig. D.4 . Avoided emissions from an electricity system perspective increase linearly with larger PV system capacities, but also show a small decrease with bigger storage capacities. Emissions from a dwelling perspective show a different behaviour. When no storage is included, then larger PV systems lead to a reduction of avoided GHG emissions. With storage, the largest emission reductions are obtained with a 1 kWp PV system. This is especially visible for a 2 kWh storage capacity per MWh demand. Smaller PV systems have a higher direct self-consumption, therefore the avoided emissions with self-consumption by storage are lower. Larger PV systems have much more electricity production resulting in more storage and avoided emissions. However, these larger PV systems (2 kWp) also have higher emissions from manufacturing, which results in lower net avoided emissions compared to a 1 kWp PV system.

Appendix E. Extended sensitivity analysis on the discounted payback period
The combined impact of a change in natural gas tariff and a change in electricity tariff is presented in Fig. E.1 . A combined increase in electricity tariff and decrease in natural gas tariff results in a rapid increase of DPBP. An annual increase of 1.4% in electricity tariff shows an average DPBP of the GSHP of 35 years. Yet, a 2.2% annual increase in natural gas tariff decreases the average DPBP of the GSHP of 20 years. The change in electricity tariff has a lower influence on dwellings with only a PV system than for dwellings with PV and storage. The average DPBP of the latter system will be lower than 15 years with an electricity tariff increase of 2.9% per year. The value of self-consumption increases with higher electricity tariffs, thus the revenues of storage are higher. Lower electricity tariff are positive for GSHP payback periods, but negative for PV and battery storage payback periods. Consequently, GSHP combined with PV and storage level out the influences of higher or lower electricity tariffs, as can be seen in subplot (c).
The influences of a higher GSHP investment subsidy combined with four other parameters are shown in Fig. E.2 . The combination of a higher GSHP investment subsidy and a decrease in natural gas tariff show a strong reduction in DPBP. For example, a 10% increase in subsidy combined with a natural gas tariff increase of 1.3% results in average DPBP of 20 years. Also, the combination with CO 2 taxation is promising to obtain lower DPBP. With a significantly increase in investment subsidies, higher discount rates and electricity tariffs become feasible to obtain DPBP below 25 years.  Extended analyses on the impact of the GSHP investment subsidy combined with the impact of the natural gas tariff (a), the electricity gas tariff (b), the CO 2 taxation (c) and the discount rate (d) on the average discounted payback period of the 16 residential dwellings with a GSHP system only. The red dot indicates the reference scenario value.