Cost-Effectiveness of Carbon Emission Abatement Strategies for a Local Multi-Energy System—A Case Study of Chalmers University of Technology Campus

This paper investigates the cost-effectiveness of operation strategies which can be used to abate CO2 emissions in a local multi-energy system. A case study is carried out using data from a real energy system that integrates district heating, district cooling, and electricity networks at Chalmers University of Technology. Operation strategies are developed using a mixed integer linear programming multi-objective optimization model with a short foresight rolling horizon and a year of data. The cost-effectiveness of different strategies is evaluated across different carbon prices. The results provide insights into developing abatement strategies for local multi-energy systems that could be used by utilities, building owners, and authorities. The optimized abatement strategies include: increased usage of biomass boilers, substitution of district heating and absorption chillers with heat pumps, and higher utilization of storage units. The results show that, by utilizing all the strategies, a 20.8% emission reduction can be achieved with a 2.2% cost increase for the campus area. The emission abatement cost of all strategies is 36.6–100.2 (AC/tCO2), which is aligned with estimated carbon prices if the Paris agreement target is to be achieved. It is higher, however, than average European Emission Trading System prices and Sweden’s carbon tax in 2019.


Introduction
With growing amounts of distributed energy resources (e.g., distributed generation and energy storage) and introduction of new demands (e.g., electrical vehicles, heat pumps, etc.) in the distribution systems, local energy management systems have received attention because they provide, amongst others, the following benefits for energy systems: • Facilitating activation of small flexibilities in the system [1] • Reducing transmission and distribution losses [2] • Managing congestion on the distribution level [3] • Increasing awareness and engagement of end-users [4] Furthermore, multi-energy systems (MESs) have been suggested to increase the synergies in the energy system as a whole by integrating and managing different energy carriers such as electricity, district heating, district cooling, and natural gas simultaneously [5]. Subsequently, analyzing the combination of local energy management systems and MESs is an interesting emerging field of research.
Within the two research areas of local energy management systems and MESs, previous studies [3,5,6] have reviewed definitions, trends, challenges, and categorization of literature providing valuable insight into the topic. For example, Grosspeithsch et al. [6] categorized the literature into four categories: general overview, model and optimization, energy management and system analysis, and case study.
One feature of the model and optimization category is that energy systems have traditionally been modeled based solely on cost minimization objectives. However, multi-criterion optimization can help broaden decision making to consider cost, the environment, reliability, social impact, utilization of renewable energy, etc. [7]. As global concerns about greenhouse gas emissions increase, carbon emissions become an increasingly important criterion to be considered in optimizing the operation of local MESs. For instance, Majidi et al. [8] proposed a cost and emission framework to assess demand response programs, and Bracco et al. [9] developed a multi-objective model to evaluate the operation of a multi-energy system considering four different building types and three energy carriers (heat, gas, and electricity). Wang et al. [10] demonstrated that a multi-objective optimization will not give one single solution but rather a set of Pareto optimal solutions. Often, the objectives are conflicting and different approaches to solve the minimization problem exist, e.g., mixed integer linear programming with weighted sums [11], evolutionary algorithms [12], game theory [13], particle swarm optimization [14], genetic algorithms [15], etc.
Furthermore, an optimization model can have a short foresight (close to real-time) or a long foresight depending on the purpose and the characteristics of the energy technologies included in the system. Optimizations with long foresight result in a more optimal management of resources especially in energy systems with seasonal storage, conventionally dispatchable units, and perfect foresight. However, such long-term optimizations require long-term forecasts and can be computationally expensive as the size and complexity of the model increases [16].
On the other hand, optimization with a short foresight lowers the computational time (which would be of value when simulating complex systems) [17] and have lesser challenges with quality of forecasts. This is especially important for systems with a large share of renewable energy sources (RES) since as the share of intermittent RES increases in the system, their stochastic nature starts to affect forecasts, availability, and prices of energy carriers. Therefore, if the model represents a system including a large share of RES, or reacts in response to the energy prices from a system with a large share of RES, a close to real-time modeling approach with short foresight can represent agents' (energy technologies') behavior closer to reality [18].
As an example from the research area of modeling and optimization of MESs, Wu et al. [11] investigated the simultaneous optimization of annual cost andCO 2 emissions in the design and operation of a distributed energy network where distributed energy sources (DESs) can exchange heat with each other through pipelines. A mixed integer linear programming (MILP) model with a weighted sum approach is used in this multi-objective optimization. Di Somma et al. [19] also used a weighted sum approach in developing a multi-objective linear programming model considering both cost and emissions. The contribution various energy technologies have on achieving the objective function was evaluated by sensitivity analyses. A limitation of this study is that it was carried out only on one customer and not a community of customers. Falke et al. [20] developed a multi-objective model for design and operation of distributed energy systems using a heuristic optimization approach. The model decomposes into three sub-model of heating network planning, buildings' renovation planning, and operation simulation. However, cooling loads and district cooling is not considered in this study. Yan et al. [21] studied the operation optimization of multiple distributed energy systems where the emissions are considered in the form of monetary costs through a carbon tax. The DESs, in this model, can exchange electricity and thermal energy with each other and electricity can be sold back to the grid. Although emissions cost is considered through carbon tax in this study, the trade-off between the emissions and the monetary costs are not discussed. In [12], an evolutionary algorithm is used to solve a multi-objective isolated MES model with high share of renewables, including investment in RES as decision variables. The paper shows that different operational approaches may be beneficial for different seasons. In [13], a scheduling method for MES is proposed using game theory. The model was found to be robust and useful for solving problems with uncertainties in, for example, wind production. In [22], a MES model is developed which includes possible constraints in the energy flow within the MES. For the electricity network, this is accomplished using a DC load flow model, and a pipeline load flow model is used for the natural gas network.
In summary, the literature contains multi-objective optimization considering both cost and emissions with different energy carriers. However, to the best of the authors' knowledge, a study which specifically extracts emission abatement strategies from multi-objective optimization models and evaluates the abatement cost for these strategies is lacking. Such a study could provide insights regarding carbon pricing and investigate the possibilities of operating local MESs in a more environmentally responsible manner. In addition, no previous study has been found which multi-objectively optimizes the three energy carriers (electricity, district heating, and district cooling) using a short foresight rolling horizon over a whole year. The benefit of considering these three energy carriers is that we can capture the synergies for emission abatement through technologies such as heat pumps.
This study was carried out within an European Union project called Fossil-free Energy District [23] whose purpose was to build a testbed using Chalmers University's campus to implement and evaluate MESs alongside a local energy market for electricity, heating, and cooling. As local energy management systems, MESs and environmental aspects are gaining momentum in the research field of energy systems, this study aims mainly to investigate: • what operation strategies can be utilized to abate CO 2 emissions in a local multi-energy system; and • how cost-effective these strategies can be.
To the best of our knowledge, the combination of the following aspects is what distinguishes this study from other studies within the field:

•
Defining emission abatement strategies in MESs' operation and evaluating cost-effectiveness of these strategies across different carbon prices • Optimizing a local multi-energy system over a year with a short foresight rolling time horizon • Multi-objective optimization of three energy carriers: district heating (DH), district cooling and electricity.
In this study, the emission abatement strategies are introduced by a short foresight rolling horizon, multi-objective optimization model which is formulated with the weighted sum approach, considering costs and emissions as the objectives. The problem is solved with a rolling time horizon of ten hours to assess the real-time behavior of different technologies. Afterwards, the cost of these strategies is evaluated against the carbon prices to provide insight on the cost effectiveness of these strategies. The study includes a coupled system of district heating, district cooling, and electricity networks. Moreover, a wide range of distributed energy technologies such as biomass combined heat and power (CHP), biomass boilers, electrical heat pumps (HP), absorption chillers, battery energy storage (BES), building inertia thermal energy storage (BITES), cold water basins (CWB), and solar photovoltaic panels (PV) are considered.

The Paper's Structure
The rest of this paper is organized as follows. Section 2 presents the Chalmers University campus energy system model, its general formulation, and its input data. Section 3 presents the emission abatement strategies and their corresponding costs. In Section 4, the incentives behind each strategy and their emission abatement costs are discussed. Finally, the study is concluded and summarized in Section 5.

Multi-Energy System Model of Chalmers University of Technology Campus
The model developed and used in this study is modeled after the energy system of Chalmers University of Technology's Johannesburg campus. The model is open-source and available on Github [24]. The model is implemented using General Algebraic Modeling System (GAMS) [25] and Matlab [26]. Only the objective function and the structure of the model are presented in detail in the main body of this article; however, the formulations for the energy technologies are presented in Appendix A. Figure 1 shows the flow chart of the dispatch model structure. The flow chart consists of three main modules, which are the measurement module, the MatLab module, and the GAMS module.

Demand profiles Local generation
Weather data BITS characterisitcs Prices Other constraints

Load and radiation forcasting
Pre-processing of data Exporting data to GAMS Post-processing of GAMS outputs The first module or the measurement module reads input data to the dispatch model, including demand data, generation data, BITES data, price data, solar irradiance, outdoor temperature data, and energy storage data.
The second module of the model consists of two parts. In the first part, forecasts of the input data to the dispatch model over a given time horizon are made and the input data are then pre-processed to be exported to the GAMS optimization model. The second part is post-processing of the GAMS outputs.
The third module is the GAMS optimization model. This model takes forecast data, measurement data, and network data from the Matlab module and simulates the optimal dispatch of the local generating units, the storage units, and demand side flexibility for a given objective. The optimization and dispatch of units is performed for the desired forecast horizon, e.g., 10 h. The model is designed to perform the optimizations with a rolling time horizon, i.e., after the first optimization, the forecast is updated with new data and a new optimization is performed in order to redispatch the demand and generation units. Between each optimization, the results are stored to be used for evaluation and comparison with the actual dispatch of the system.
To better understand the formulations and results, a simplified overview of the structure of the campus's energy system is provided in Figure 2. In this figure, HP R , HP S , and HP W1,2 are the four heat pumps in the system. HP R is a refrigeration heat pump, HP S is a heat pump used in the summer time, and HP W1,2 are the heat pumps used during the heating season. Note that the absorption chiller can only receive heating from the external district heating network. Further note that the CHP unit has the flexibility to be used purely for heating or partly heating with electricity production. B+FGC represents campus's biomass boiler (B) with a flue gas condenser (FGC). There are three storage types in the system. Thermal energy is stored in buildings' thermal inertia (BITES), a cold water basin (CWB) stores chilled water, and battery energy storage (BES) is used to store electricity. In Chalmers' campus energy system model, district heating, district cooling, and electricity networks are coupled by means of heat pumps, absorption chillers, and a combined heat and power (CHP) unit:

•
Heat pumps: This technology can be placed between the three energy carriers (district heating, district cooling, and electricity) by assigning the hot source to district heating network and cold sink to district cooling network (see Figure 3).

•
Absorption chiller: This technology can provide flexibility in the system by using heat from district heating network to produce cooling for the district cooling.

•
Combined Heat and Power: The biomass boiler and turbine can provide electricity and heating simultaneously to the electricity and district heating networks. Flexibility is introduced by the power-to-heat ratio (alpha) that varies the amount of fuel to be converted to electricity or heating. Although coupling of the heat pumps to the district heating and cooling systems does provide flexible and efficient production of either cooling or heating, it also imposes limitations on the load factor of heat pumps. This limitation is necessary to maintain balance between energy production and consumption in the district heating and district cooling networks. For example, if the heating demand is high while cooling demand is low, full capacity usage of heat pumps would not be possible to prevent the formation of frost in the district cooling network. This trade-off between energy efficiency and load factor is discussed more in Section 4.
In this study, a close to real-time simulation (one hour ahead) with a rolling time horizon of 10 h was chosen. This forecast horizon and optimization window were chosen so as to be technically and realistically achievable in actual energy systems. Figure 4 shows how the optimizations are carried out under the rolling-time horizon.

. . .
Step i  In each step, the dispatch of units is decided based on the assumed forecast horizon ( f h). In this way, the results are more closely based on real-time decision making and the solver will not be able to see the hours with less accurate forecasts. The steps are continued until the whole simulation period (one full year in this study) is covered.

Model formulation of Chalmers Campus Energy System
In this section, the main formulations of the model such as objective function, energy balances, and general constraints of the energy technologies are presented. More detailed formulations of the energy technologies can be found in Appendix A.

Objective Function
To consider both economic and environmental objectives, the weighted sum method has been used. The objective function is calculated based on Equation (1).
In Equation (1), C t and E t represent the cost and the emissions of the system at hour t. The total cost and emissions are summed and minimized for each time interval of i with the span of f h hours.
α and β are the weighting factors for each objective. The cost weighting factor α is fixed at 1 while the environmental weighting factor β is changed from zero to 1. Thus, in other words, β can be seen as an additional carbon tax on the emissions. This approach for changing α and β is chosen to make the carbon emission reductions comparable with carbon prices. The units for C t and E t are SEK and gCO 2 .
In Equation (2), C t includes fixed costs of the assets (C f ix ), variable costs (C var ), revenues from exporting electricity and heat to the external networks (R exp ), and costs for importing from external networks. Subscript j represents the energy technologies shown in Figure 2. The cost for electricity or district heating input for technologies is counted in C imp and not in C var nor C f uel .
E t , which is calculated according to Equation (3), accounts for the emission for the imports from the outer grids and the assets inside the MES.
P and Q represent the electric power and the heating power. EF exg is the external grid's emission factor. E j accounts for emissions from each energy technology.

Energy Balances
Demand-supply balances are presented in Equations (4)-(6) for electricity, heating, and cooling. In the heating system (Equation (5)), over-production is possible with help of a cooling tower which is available in the campus' energy system. The dem subscript represents the demand of the campus.

Energy Technologies
As shown in Figure 2, there are different generation and storage technologies in the system. The formulation of energy technologies are provided in Appendix A and their characteristics are presented in Appendix B. In this section, only the general explanation for the technologies is provided.
All the generation units are limited by their maximum capacity. For the case study at the Chalmers campus, some site-specific limitations were set to better represent the actual operating condition at the campus. For example, during weekday working hours, the biomass boiler and CHP units are not dispatchable as they are assumed to be used for research purposes at those times. Moreover, exporting district heating to the external grid is limited to the capacity of biomass boiler because of the network structure. Electricity production from the CHP unit is flexible, so that the heat production of its boiler can feed the turbine or go directly to the internal district heating network. The absorption chiller is only operational during the cooling season and can only receive heating from the external district heating network due to the network's design. The production of heat pumps is calculated based on Equations (7) and (8), where the cooling and heating is produced at the same time. COP k and COP q are the coefficient of performance for cooling and heating, respectively. HP S is used only in cooling season while HP W1−2 are available in heating seasons. HP R is a refrigeration machine and only 20% of its heat production can be used in the system due to the network configuration.
For both storage units, the cold water basin and battery energy storage, charge/discharge efficiency has been considered. Moreover, these technologies are constrained by their charging/ discharging power and energy storage capacity. BITES has been modeled based on a methodology using two interlinked storages, a shallow and a deep storage [27].

Input Data
Input data for energy prices, marginal emissions, demand profiles, weather data ,and energy technologies are presented in this section.
Both the district heating and the electricity prices considered in this study are hourly dynamic prices. Electricity prices are obtained from Nordpool day ahead market for region SE3 year 2016-2017. To present a more realistic price on the consumer side, the sum of the spot price, network tariffs, average profits, and (when applicable) electricity certificates are used. Currently, the pricing system in the local district heating network is based on seasonal heating prices with a peak demand charge and an efficiency rebate [28]. However, to create more dynamic interactions between energy carriers, hourly district heating prices are modeled based on a typical generation mix in a Swedish district heating system (see Appendix B).
The marginal emissions of the external electricity grid are obtained from Electricity Map [29]. This methodology uses a statistical method to predict the marginal emissions for every hour based on historical data. For district heating, the marginal emissions are calculated based on the typical marginal unit in the generation mix which is in turn determined by the hourly district heating prices (see Appendix B).
Demand profiles for electricity, heating, and cooling are provided by the primary real estate owners on the campus, Akademiska Hus and Chalmers Fastigheter. The demand profiles are presented in Appendix B.
Weather data for outdoor temperature are from a local weather station and solar radiations are based on [30].
The capacity, performance, and similar energy technology parameters are obtained from the MES operator at Chalmers campus which are presented in Appendix B.

Key Results
In this section, the results are presented showing the trade-off between cost and emissions and the corresponding abatement strategies. Discussion of the results can be found in Section 4.
The trade-off between cost and emissions are presented in Figure 5. It can be seen that, at lower β values, higher emissions reductions can be achieved with smaller changes in the cost. The relative changes in cost and emissions in different beta values are presented in Table 1.

Emission Abatement Strategies
To show the emission abatement strategies, the aggregated annual change in usage of each technology is presented in Figure 6 for heating, cooling, electricity, and storage technologies. The values represent the difference in usage (kWh/year) at each β value, compared to the cost optimum case (where β = 0). The lines for coupling technologies (e.g., heat pumps and the absorption chiller) can represent consumption or production in different energy carriers. For example, the heat pump lines in the heating and cooling graphs are the production of heat pumps while in the electricity graph they show the electricity consumed by the heat pumps. Note the very different scales used in the graph for each energy carrier.
Moreover, at different β values, there is a mix of different strategies for abatement. The β spectrum is divided into five phases (Phases (I)-(V)) as in Figure 6 based on when these strategies start and end. The range of β in each phase can be seen in Table 1. The dominant strategies at each phase are summarized in Figure 7. The observed strategies in each phase are: • Phase (I): In this phase, although the COP of HP w2 is higher than HP w1 , the MES optimizer decides to replace HP w2 with HP w1 . Moreover, the electricity and heat production from the CHP unit is increased which leads to a decrease in electricity and district heating imports. Another strategy started in this phase is the increase in usage of storage. Furthermore, it's observed that the operation of CHP starts to move towards relatively greater heat production than electricity and consequently more electricity imports.    Figure 6. The change in the total energy per year (consumed, produced or stored) of each technology compared to the cost optimal case versus beta for heating, cooling, electricity, and storage.  Figure 7. CO 2 abatement strategies at different phases. Figure 8 is an example of the simulated dispatch of the campus energy system at three example beta values of 0, 0.001, and 1. A few different abatement strategies can be observed in this figure such as increase in heat production from the biomass boiler, increase and then decrease of electricity production from the CHP unit, increase in heating exports, replacement of HP W 2 with HP W 1, replacement of HP W 2 with HP R , and increase in utilization of the cold water basin and building inertia thermal energy storage. The strategy related to increase in using battery energy storage cannot be noticed in this figure due to its small scale compared to the electricity demand. Moreover, replacement of the absorption chiller with heat pumps cannot be seen because the absorption chiller is not in operation during the winter.

Cost of Strategies
The cost abatement in each phase is presented in the form of total abatement cost and marginal abatement cost in Table 2. Total abatement cost at each emission weighting factor (β) is calculated based on Equation (9). It includes the total abated CO 2 over the year divided by added costs compared to the cost minimization case. Marginal abatement cost of each β is calculated by Equation (10), which is the same as total abatement cost but, instead, compared to the previous β.  In Table 2, it can be seen that, as β increases, the cost effectiveness of strategies decreases as it gets more costly to avoid further emissions. Moreover, the difference between total abatement cost and marginal abatement cost shows that most of the emission reductions have been carried out in the cheaper Phases (II) and (III). This can be seen in Figure 5b.

Discussion
In this section, first the incentives behind each strategy are explained and then the cost of CO 2 emission abatement is discussed.

Incentives Behind The Strategies
It has been observed that one of the strategies was switching from HP w2 to HP w1 and later to HP R . This strategy is interesting because the COP values for HP w2 are higher than the unit it is substituted with. The reason for this switching is connected to the strict coupling of district heating and district cooling system. As shown in Table A1, COP q COP c for HP w2 is lower than the other two which means less heating can be obtained for each unit of cooling produced. From analyzing the input data by Equation (11), we know that there is 39% higher potential for emission savings in the heating network than in the electricity network.
Therefore, the optimizer decides to use the heat pumps with higher COP q COP c and produce as much heat as possible with cooling production as low as possible to prevent violating energy balance of internal cooling network and thus forming frost in the district cooling network. Moreover, if the district cooling network is cooled extra, the temperature would fall, and this would cause the COPs to be reduced and also cause non-linearity in the model. This behavior therefore shows the system to be limited in the low grade heat that is available from the district cooling system, which artificially constrains the dispatch of the available heat pumps.
This energy system would therefore benefit from bore holes or other low grade heat sources which would lead to more dispatch of the higher efficiency heat pumps.
Note, furthermore, that, although the district heating emissions in summer are lower than electricity in the summer, the absorption chiller is replaced by HP S . This is because the COP of the heat pump is higher than that of the absorption chiller, so much so that it compensates for the lower emissions in the district heating system. In this study, the COP c of the absorption chiller is assumed to be 0.5 compared to the assumed COP value of 1.8 for HP S .
We also see that the emissions in the external district heating system are higher than the biomass boiler emissions. This causes increased production and even export from the boiler to the district heating network to reduce overall emissions. This strategy plays a large role in emission abatement since 93.3% of the total emission in the emission minimization case (β = 1) are achieved during the phases which this strategy is utilized (i.e., Phases (II) and (III)).
As mentioned above, the heat production from the CHP unit for all β values is higher than the cost minimization case because there is a high potential to reduce emissions in the heating system. Regarding electricity production from the CHP, however, the increase in usage is observed only up to the end of Phase (II). From Phase (III) onwards, the electricity production is replaced by heat production because β is large enough to outweigh the loss of revenues from electricity production and more emissions can be saved on heating than electricity.
The use of storage technologies, mainly begins to increase from Phase (III). The incentive behind this trend is that storage is used to gain upon the variations in emission levels between consecutive hours. This happens later in the phases because the round-trip storage losses are considerable compared to the abatement gained from mitigating the volatility in emission levels.

Cost-Effectiveness of Abatement Strategies
The cost of abatement for MES can be evaluated at different magnitudes, locations and from different viewpoints by comparing three different carbon pricing schemes. The first is the prices from carbon markets (e.g., European Emission Trading System (EU ETS)), the second is the carbon taxes and the third is the emission abatement cost of similar pilot projects.
The average carbon price on EU ETS in 2019 was 24.9 (A C/tCO 2 ) [31]. According to the current EU ETS regulations [32], "combustion of fuels in installations with a total rated thermal input exceeding 20 MW" are considered under the cap. This means the installations on this campus would not be regulated under the current system. However, if the regulations were to change in the future to include smaller installations in the cap and trade system, then the prices of carbon allowances would have to increase to the total abatement cost levels in Table 2 for abatement strategies to be cost efficient. Otherwise, it would be more cost efficient to buy emission allowances. According to current EU ETS prices, the emission abatement strategies would not currently be cost-effective. However, according to Stiglitz et al. [33], the carbon prices to achieve Paris agreement's goal need to be 36-73 (A C/tCO 2 ) by 2020 and 45-91 (A C/tCO 2 ) by 2030. This suggests that the discussed abatement strategies up to Phase (III) can be economically competitive measures if the Paris agreement target is to be achieved.
In this study, the carbon taxes are included in the pricing of energy carriers. Therefore, the introduced abatement strategies can be economically beneficial if the tax levels are increased to or above the total abatement cost levels of each phase plus the current considered tax level. The carbon tax is applicable to the installations which are not under the EU ETS system, such as small MESs.
Another comparison can be with similar pilot projects aiming to reduce emissions. Klimatklivet is a program in Sweden, which supports such projects [34]. The focus of this program is the non-ETS sector and the funding is appointed to projects which have the highest abatement per investment. Iseberg et al. [34] evaluated Klimatklivet projects and found that the abatement cost is within a range of 10-80 (A C/tCO 2 ) over the other policies such as a carbon tax. These values show that the discussed abatement strategies until Phase (III) would be cost competitive with other pilot projects.

Conclusions
This study set out to find operation strategies which can be utilized for carbon emission abatement in a local multi-energy system and investigate the cost-effectiveness against carbon prices in carbon markets, carbon tax scheme, and cost-effectiveness of similar pilot projects.
The results of the case study show that, by utilizing all the abatement strategies, a 20.8% emission reduction can be achieved with a 2.2% increase in cost. The abatement strategies which were extracted from the optimization results include: more usage of biomass boilers in heat production, substitution of district heating and absorption chillers with heat pumps, and higher utilization of storage units. It should be noted that the system was shown to be limited in the low grade heat that was available from the district cooling system, which artificially constrained the dispatch of the available heat pumps. This system would therefore benefit from bore holes or other low grade heat sources which would lead to more dispatch of the higher efficiency heat pumps. Furthermore, the utilization of the CHP unit has shown to be sensitive to the relative weighting of emissions vs. cost in the objective function. The relative share of electricity production from the CHP unit is also shown to decrease at higher emissions weighting factors (above 10 −3 ) due to the relatively higher emissions in the district heating system compared to the electricity system. This analysis demonstrates that across all abatement strategies the total carbon dioxide abatement cost is 36.6-100.2 (A C/tCO 2 ), which is higher than both the average carbon price in EU ETS and carbon tax prices in Sweden in 2019, but at the same level as similar pilot projects in Sweden. Furthermore, all strategies up to Phase (IV) are cost competitive if the carbon prices are adjusted to reach the Paris agreement target.
The results from this study can provide insight that carbon prices might not yet be high enough to make the multi-objective operation of similar local MESs cost-effective. However, the suggested abatement strategies for similar local MESs can be cost-competitive in achieving Paris agreement's goal if carbon prices were aligned to reach the agreement's goal. Of course general conclusions can hardly be made from only one case study, and the presented conclusions are only based on this case study.
Future studies can investigate more generalized abatement strategies and cost-effectiveness evaluation for local multi-energy systems. A similar methodology can be utilized for case studies in different countries and with different characteristics, and the results can be compared to provide a broader insight into the research questions. Moreover, studying abatement strategies considering future scenarios for emissions and prices in electricity and district heating systems can lead to further insights in the operation planning of local multi-energy systems.
Supplementary Materials: The following are available online at http://www.mdpi.com/1996-1073/13/7/1626/ s1, the simulation model's code, input data files for the model, and the detailed calculations for emission factors. Further details are provided in "Readme.txt" at the address above.

Abbreviations
The following abbreviations are used in this manuscript:

MES
As mentioned in Section 2.1.3, the biomass boiler and the CHP unit are not always dispatchable. Therefore, they are regularly fixed to specific levels (Equation (A2)). Moreover, HP S , HP W1 , HP W2 , and the absorption chiller are only in operation during specific periods of the year, as also represented by Equation (A2).

. Biomass Boiler (B) and Flue Gas Condenser (FGC)
The heat produced in the boiler (Equation (A3)) is related to the fuel input (Q f uel,B,t ) and the efficiency of the boiler (η B ) . Besides, the heat recovered from the flue gas condenser (Equation (A4)) is related to the heat not absorbed by boiler and the efficiency of the condenser (η f gc ).
Since the flue gas condenser and the biomass boiler are in the same unit, their heat is aggregated and presented as one unit in this study (Equation (A5)).
Moreover, the changes in boiler's output is limited to the ramp-up (RU B ) and ramp-down (RD B ) limits (Equation (A6)).
As shown in Figure 2, the export to external district heating can only come from the biomass boiler and therefore is limited to its production (Equation (A7)).
discharging (η dis ). Charging and discharging efficiencies are considered to be equal in each unit. Moreover, for the cold water basin, instead of P, K is used in Equation (A16).
Appendix A.6.1. Cold Water Basin (CWB) Due to local system's limitations, discharging of the cold water basin is only limited to the cooling load of the building it is situated in (K dem,y,t ). This constraint is presented in Equation (A17).
Moreover, it has been observed that only considering the charging and discharging efficiencies will not stop simultaneous charge and discharge of the cold water basin. The reason is that the cooling demand is lower than heating demand in the district as a whole and there is not enough cooling storage to handle the excess produced cooling. Therefore, the optimizer would theoretically curtail cooling by charging and discharging at the same time, causing additional losses in the cooling system.
To prevent simultaneous charging and discharging, binary variables are used (Equation (A18)) where b is the binary variable and M is a big number (10 6 kW in this case study).
For the battery units, an additional constraint is considered on the minimum state of charge (SOC min ). This constraint (Equation (A19)) is to prevent high degradation costs.

Appendix B. Input Data
In this section, the input data for energy technologies, marginal emissions, marginal prices, and energy demand of the campus are provided. The parameters for energy technologies are presented in Table A1.
The heating and cooling seasons, which are considered for the heat pumps' and absorption chiller's operation is presented in Table A2.
Results from simulation models are very dependent on their assumptions and therefore other emission factors or marginal prices can affect the simulation results. With the developed model, similar studies can be performed to assess efficient abatement strategies in other locations or with different assumptions. The marginal emissions and marginal prices for external electricity and district heating networks that are used in this study are presented in Figures A1 and A2. The emission factors used for calculating the marginal emissions are presented in Tables A3 and A4. All electricity efficiencies and emission factors are based on [35,36], except for oil and waste incineration where the work of Bertoldi et al. [37] is the source. Pumped hydro and 'Unknown' are taken from [29]. Pumped hydro is based on average Swedish emissions during times when hydro is being charged. 'Unknown' is based on the average generation plant in the Swedish electrical grid. The emission factors for CHP units are calculated based on the Benefit Allocation Method [38].  In the district heating system, the emission factor (Table A3) for the marginal unit is considered to be the marginal emission factor for every hour. The emission factors are primarily based on the energy allocation method which may differ from current practice among district heating operators in Sweden. More details on the calculation of the emissions factors can be found in Supplementary Materials based on the following references [35][36][37][38][39]. The hourly marginal unit data are taken from a district heating operator in Sweden. Based on the information from the district heating provider, for this dataset, the most expensive unit at each hour has been assumed to be the marginal unit. However, the common practice for marginal unit allocation in the external district heating system is not necessarily the most expensive unit because of other constraints such as the ramp-up and ramp-down limitations of the units. However, due to limited data access on the external district heating network, this assumption was considered reasonable. For calculating the marginal emissions for the electricity grid, a statistical method is used to predict the marginal emissions at every hour based on historical data [40]. Firstly, the probability that each technology is on the margin based on a consumed unit of electricity in Sweden is taken from Electricity Map [29]. Secondly, the probabilities are multiplied by the relevant emission factor of each technology and summed. The emission factors for each technology can be seen in Table A4. The differing approaches for calculating marginal emissions in the electricity and district heating networks are due to the respective size of these systems and therefore the ability to know precisely which unit is on the margin. The district heating network of a city is usually smaller than the electricity network and limited to the boundaries of the city. Identification of the exact marginal unit can therefore be easily estimated. However, the electricity network is much larger and tangled with production units in other countries through import and export across regional and country borders. Since the exact marginal unit cannot be identified, the electricity system's marginal emissions are determined by statistical methods, as described above.
Energy demand of the campus for electricity, heating, and cooling for the simulation period is presented in Figure A3. The simulation period is selected based on availability and quality of real data from the campus energy system. The simulations are run for the period from 1 March 2016 to 28 February 2017.