Emissions directly exported from GAMS and postprocessed

dispaset/GAMS/UCM_h.gms

@@ -117,7 +117,7 @@ Curtailment(n)                   [n.a]           Curtailment allowed or not {1 0
 Demand(mk,n,h)                   [MW]            Demand
 Efficiency(p2h,h)                [%]             Efficiency
 EmissionMaximum(n,p)             [tP]            Emission limit
-EmissionRate(u,p)                [tP\MWh]        P emission rate
+EmissionRate(au,p)               [tP\MWh]        P emission rate
 FlowMaximum(l,h)                 [MW]            Line limits
 FlowMinimum(l,h)                 [MW]            Minimum flow
 Fuel(u,f)                        [n.a.]          Fuel type {1 0}
@@ -532,6 +532,7 @@ EQ_SystemCost(i)..
          +sum(u,CostRampUpH(u,i) + CostRampDownH(u,i))
          +sum(u,CostVariable(u,i) * Power(u,i)*TimeStep)
          +sum(hu,CostVariable(hu,i) * Heat(hu,i)*TimeStep)
+         +sum(p2h,CostVariable(p2h,i) * Heat(p2h,i)*TimeStep)
          +sum(l,PriceTransmission(l,i)*Flow(l,i)*TimeStep)
          +sum(n,CostLoadShedding(n,i)*ShedLoad(n,i)*TimeStep)
          +sum(n_th, CostHeatSlack(n_th,i) * HeatSlack(n_th,i)*TimeStep)
@@ -551,6 +552,7 @@ EQ_SystemCost(i)..
          +sum(u,CostRampUpH(u,i) + CostRampDownH(u,i))
          +sum(u,CostVariable(u,i) * Power(u,i)*TimeStep)
          +sum(hu,CostVariable(hu,i) * Heat(hu,i)*TimeStep)
+         +sum(p2h,CostVariable(p2h,i) * Heat(p2h,i)*TimeStep)
          +sum(l,PriceTransmission(l,i)*Flow(l,i)*TimeStep)
          +sum(n,CostLoadShedding(n,i)*ShedLoad(n,i)*TimeStep)
          +sum(n_th, CostHeatSlack(n_th,i) * HeatSlack(n_th,i)*TimeStep)
@@ -1280,6 +1282,7 @@ ShadowPrice_RampDown_TC(au,h)
 OutputRampRate(au,h)
 OutputStartUp(au,h)
 OutputShutDown(au,h)
+OutputEmissions(n,p,z)
 ;
 
 OutputCommitted(au,z)=Committed.L(au,z);
@@ -1338,6 +1341,9 @@ OutputRampRate(hu,z) = - Heat.L(hu,z-1)$(ord(z) > 1) + Heat.L(hu,z);
 OutputStartUp(au,z) = StartUp.L(au,z);
 OutputShutDown(au,z) = ShutDown.L(au,z);
 
+OutputEmissions(n,p,z) = (sum(u,Power.L(u,z)*EmissionRate(u,p)*Location(u,n))
+                        + sum(hu,Heat.L(hu,z)*EmissionRate(hu,p)*Location(hu,n)))
+                        / (sum(u,Power.L(u,z)*Location(u,n)) + sum(hu,Heat.L(hu,z)*Location(hu,n)));
 
 EXECUTE_UNLOAD "Results.gdx"
 OutputCommitted,
@@ -1389,6 +1395,7 @@ ShadowPrice_RampDown_TC,
 OutputRampRate,
 OutputStartUp,
 OutputShutDown,
+OutputEmissions,
 status
 ;
 

dispaset/misc/gdx_handler.py

@@ -21,6 +21,7 @@
 import numpy as np
 import pandas as pd
 import logging
+import copy
 
 from .str_handler import shrink_to_64, force_str
 
@@ -306,8 +307,31 @@ def gdx_to_dataframe(data, fixindex=False, verbose=False):
                 logging.debug('Successfully loaded variable ' + symbol)
             elif dim == 1:
                 logging.warning('Variable ' + symbol + ' has dimension 0, which should not occur. Skipping')
-            elif dim > 3:
+            elif dim > 4:
                 logging.warning('Variable ' + symbol + ' has more than 2 dimensions, which is very tiring. Skipping')
+            elif dim == 4:
+                vars1 = set()
+                vars2 = set()
+                for element in data[symbol]:
+                    if not element[0] in vars1:
+                        vars1.add(element[0])
+                    if not element[1] in vars2:
+                        vars2.add(element[1])
+                vars1 = list(vars1)
+                vars2 = list(vars2)
+                pd_index = pd.MultiIndex.from_product([vars1, vars2], names=('Zones', 'Emissions'))
+                out[symbol] = pd.DataFrame(columns=pd_index, index=out['OutputPower'].index)
+                for element1 in data[symbol]:
+                    element = copy.deepcopy(element1)
+                    for var1 in vars1:
+                        for var2 in vars2:
+                            if (var2 == element[1]) and (var1 == element[0]):
+                                out[symbol].loc[element[2], (var1, var2)] = element[3]
+                out[symbol].index = out[symbol].index.astype(int)
+                out[symbol].sort_index(inplace=True)
+                out[symbol] = out[symbol].fillna(0)
+
+                logging.warning('Successfully loaded variable ' + symbol)
         else:
             logging.debug('Variable ' + symbol + ' is empty. Skipping')
     for symbol in out:

dispaset/postprocessing/data_handler.py

@@ -160,7 +160,7 @@ def get_sim_results(path='.', cache=None, temp_path=None, return_xarray=False, r
                    'OutputStorageSlack', 'OutputPtLDemand', 'OutputH2Output', 'OutputPowerMustRun',
                    'OutputReserve_2U', 'OutputReserve_2D', 'OutputReserve_3U', 'ShadowPrice_RampUp_TC',
                    'ShadowPrice_RampDown_TC', 'OutputRampRate', 'OutputStartUp', 'OutputShutDown','HeatShadowPrice',
-                   'H2ShadowPrice', 'OutputCurtailedHeat']
+                   'H2ShadowPrice', 'OutputCurtailedHeat', 'OutputEmissions']
 
     # Setting the proper index to the result dataframes:
     from itertools import chain

dispaset/preprocessing/build.py

@@ -717,11 +717,18 @@ def build_single_run(config, profiles=None, PtLDemand=None, MTS=0):
         found = False
         for FuelEntry in FuelEntries:
             if Plants_merged['Fuel'][unit] == FuelEntry:
-                parameters['CostVariable']['val'][unit, :] = FuelPrices[FuelEntries[FuelEntry]][c] / \
-                                                             Plants_merged['Efficiency'][unit] + \
-                                                             Plants_merged['EmissionRate'][unit] * \
-                                                             FuelPrices['PriceOfCO2'][c]
-                found = True
+                if Plants_merged['Technology'][unit] == 'ABHP':
+                    parameters['CostVariable']['val'][unit, :] = FuelPrices[FuelEntries[FuelEntry]][c] / \
+                                                                 1.55 + \
+                                                                 Plants_merged['EmissionRate'][unit] * \
+                                                                 FuelPrices['PriceOfCO2'][c]
+                    found = True
+                else:
+                    parameters['CostVariable']['val'][unit, :] = FuelPrices[FuelEntries[FuelEntry]][c] / \
+                                                                 Plants_merged['Efficiency'][unit] + \
+                                                                 Plants_merged['EmissionRate'][unit] * \
+                                                                 FuelPrices['PriceOfCO2'][c]
+                    found = True
         # Special case for biomass plants, which are not included in EU ETS:
         if Plants_merged['Fuel'][unit] == 'BIO':
             parameters['CostVariable']['val'][unit, :] = FuelPrices['PriceOfBiomass'][c] / \

dispaset/preprocessing/reserves.py

@@ -91,15 +91,13 @@ def probabilistic_reserve(config, allunits, load, AvailabilityFactors, zone,
         wton = pd.Series(0, index=load.index)
 
     # 5% of the windoff forecast
-
     if units['Technology'].str.contains('WTOF').any():
         tmp = units.index[units['Technology'] == 'WTOF']
         wtof = AvailabilityFactors[tmp].agg(func).sum(axis=1)
     else:
         wtof = pd.Series(0, index=load.index)
 
     # 5% of solar forecast
-
     if units['Technology'].str.contains('PHOT').any():
         tmp = units.index[units['Technology'] == 'PHOT']
         phot = AvailabilityFactors[tmp].agg(func).sum(axis=1)