Passed
Pull Request — dev (#924)
by
unknown
01:39
created

data.datasets.sanity_checks.sanitycheck_eGon2035_electricity()   F

Complexity

Conditions 16

Size

Total Lines 230
Code Lines 97

Duplication

Lines 64
Ratio 27.83 %

Importance

Changes 0
Metric Value
eloc 97
dl 64
loc 230
rs 1.7781
c 0
b 0
f 0
cc 16
nop 0

How to fix   Long Method    Complexity   

Long Method

Small methods make your code easier to understand, in particular if combined with a good name. Besides, if your method is small, finding a good name is usually much easier.

For example, if you find yourself adding comments to a method's body, this is usually a good sign to extract the commented part to a new method, and use the comment as a starting point when coming up with a good name for this new method.

Commonly applied refactorings include:

Complexity

Complex classes like data.datasets.sanity_checks.sanitycheck_eGon2035_electricity() often do a lot of different things. To break such a class down, we need to identify a cohesive component within that class. A common approach to find such a component is to look for fields/methods that share the same prefixes, or suffixes.

Once you have determined the fields that belong together, you can apply the Extract Class refactoring. If the component makes sense as a sub-class, Extract Subclass is also a candidate, and is often faster.

1
"""
2
This module does sanity checks for both the eGon2035 and the eGon100RE scenario
3
separately where a percentage error is given to showcase difference in output
4
and input values. Please note that there are missing input technologies in the
5
supply tables.
6
Authors: @ALonso, @dana, @nailend, @nesnoj
7
"""
8
9
from sqlalchemy import Numeric
10
from sqlalchemy.sql import and_, cast, func, or_
11
import numpy as np
12
import pandas as pd
13
14
from egon.data import config, db, logger
15
from egon.data.datasets import Dataset
16
from egon.data.datasets.electricity_demand_timeseries.cts_buildings import (
17
    EgonCtsElectricityDemandBuildingShare,
18
    EgonCtsHeatDemandBuildingShare,
19
from egon.data import config, db
20
from egon.data.datasets import Dataset
21
22
from egon.data.datasets.emobility.motorized_individual_travel.db_classes import (
23
    EgonEvCountMunicipality,
24
    EgonEvCountMvGridDistrict,
25
    EgonEvCountRegistrationDistrict,
26
    EgonEvMvGridDistrict,
27
    EgonEvPool,
28
    EgonEvTrip,
29
)
30
from egon.data.datasets.emobility.motorized_individual_travel.helpers import (
31
    DATASET_CFG,
32
    read_simbev_metadata_file,
33
)
34
from egon.data.datasets.etrago_setup import (
35
    EgonPfHvLink,
36
    EgonPfHvLinkTimeseries,
37
    EgonPfHvLoad,
38
    EgonPfHvLoadTimeseries,
39
    EgonPfHvStore,
40
    EgonPfHvStoreTimeseries,
41
)
42
from egon.data.datasets.scenario_parameters import get_sector_parameters
43
44
TESTMODE_OFF = (
45
    config.settings()["egon-data"]["--dataset-boundary"] == "Everything"
46
)
47
48
49
class SanityChecks(Dataset):
50
    def __init__(self, dependencies):
51
        super().__init__(
52
            name="SanityChecks",
53
            version="0.0.5",
54
            dependencies=dependencies,
55
            tasks={
56
                etrago_eGon2035_electricity,
57
                etrago_eGon2035_heat,
58
                residential_electricity_annual_sum,
59
                residential_electricity_hh_refinement,
60
                cts_electricity_demand_share,
61
                cts_heat_demand_share,
62
                sanitycheck_emobility_mit,
63
            },
64
        )
65
66
67
def etrago_eGon2035_electricity():
68
    """Execute basic sanity checks.
69
70
    Returns print statements as sanity checks for the electricity sector in
71
    the eGon2035 scenario.
72
73
    Parameters
74
    ----------
75
    None
76
77
    Returns
78
    -------
79
    None
80
    """
81
82
    scn = "eGon2035"
83
84
    # Section to check generator capacities
85
    logger.info(f"Sanity checks for scenario {scn}")
86
    logger.info(
87
        "For German electricity generators the following deviations between "
88
        "the inputs and outputs can be observed:"
89
    )
90
91
    carriers_electricity = [
92
        "other_non_renewable",
93
        "other_renewable",
94
        "reservoir",
95
        "run_of_river",
96
        "oil",
97
        "wind_onshore",
98
        "wind_offshore",
99
        "solar",
100
        "solar_rooftop",
101
        "biomass",
102
    ]
103
104
    for carrier in carriers_electricity:
105
106
        if carrier == "biomass":
107
            sum_output = db.select_dataframe(
108
                """SELECT scn_name, SUM(p_nom::numeric) as output_capacity_mw
109
                    FROM grid.egon_etrago_generator
110
                    WHERE bus IN (
111
                        SELECT bus_id FROM grid.egon_etrago_bus
112
                        WHERE scn_name = 'eGon2035'
113
                        AND country = 'DE')
114
                    AND carrier IN ('biomass', 'industrial_biomass_CHP',
115
                    'central_biomass_CHP')
116
                    GROUP BY (scn_name);
117
                """,
118
                warning=False,
119
            )
120
121
        else:
122
            sum_output = db.select_dataframe(
123
                f"""SELECT scn_name,
124
                 SUM(p_nom::numeric) as output_capacity_mw
125
                         FROM grid.egon_etrago_generator
126
                         WHERE scn_name = '{scn}'
127
                         AND carrier IN ('{carrier}')
128
                         AND bus IN
129
                             (SELECT bus_id
130
                               FROM grid.egon_etrago_bus
131
                               WHERE scn_name = 'eGon2035'
132
                               AND country = 'DE')
133
                         GROUP BY (scn_name);
134
                    """,
135
                warning=False,
136
            )
137
138
        sum_input = db.select_dataframe(
139
            f"""SELECT carrier, SUM(capacity::numeric) as input_capacity_mw
140
                     FROM supply.egon_scenario_capacities
141
                     WHERE carrier= '{carrier}'
142
                     AND scenario_name ='{scn}'
143
                     GROUP BY (carrier);
144
                """,
145
            warning=False,
146
        )
147
148
        if (
149
            sum_output.output_capacity_mw.sum() == 0
150
            and sum_input.input_capacity_mw.sum() == 0
151
        ):
152
            logger.info(
153
                f"No capacity for carrier '{carrier}' needed to be"
154
                f" distributed. Everything is fine"
155
            )
156
157
        elif (
158
            sum_input.input_capacity_mw.sum() > 0
159
            and sum_output.output_capacity_mw.sum() == 0
160
        ):
161
            logger.info(
162
                f"Error: Capacity for carrier '{carrier}' was not distributed "
163
                f"at all!"
164
            )
165
166
        elif (
167
            sum_output.output_capacity_mw.sum() > 0
168
            and sum_input.input_capacity_mw.sum() == 0
169
        ):
170
            logger.info(
171
                f"Error: Eventhough no input capacity was provided for carrier"
172
                f"'{carrier}' a capacity got distributed!"
173
            )
174
175
        else:
176
            sum_input["error"] = (
177
                (sum_output.output_capacity_mw - sum_input.input_capacity_mw)
178
                / sum_input.input_capacity_mw
179
            ) * 100
180
            g = sum_input["error"].values[0]
181
182
            logger.info(f"{carrier}: " + str(round(g, 2)) + " %")
183
184
    # Section to check storage units
185
186
    logger.info(f"Sanity checks for scenario {scn}")
187
    logger.info(
188
        "For German electrical storage units the following deviations between"
189
        "the inputs and outputs can be observed:"
190
    )
191
192
    carriers_electricity = ["pumped_hydro"]
193
194
    for carrier in carriers_electricity:
195
196
        sum_output = db.select_dataframe(
197
            f"""SELECT scn_name, SUM(p_nom::numeric) as output_capacity_mw
198
                         FROM grid.egon_etrago_storage
199
                         WHERE scn_name = '{scn}'
200
                         AND carrier IN ('{carrier}')
201
                         AND bus IN
202
                             (SELECT bus_id
203
                               FROM grid.egon_etrago_bus
204
                               WHERE scn_name = 'eGon2035'
205
                               AND country = 'DE')
206
                         GROUP BY (scn_name);
207
                    """,
208
            warning=False,
209
        )
210
211
        sum_input = db.select_dataframe(
212
            f"""SELECT carrier, SUM(capacity::numeric) as input_capacity_mw
213
                     FROM supply.egon_scenario_capacities
214
                     WHERE carrier= '{carrier}'
215
                     AND scenario_name ='{scn}'
216
                     GROUP BY (carrier);
217
                """,
218
            warning=False,
219
        )
220
221
        if (
222
            sum_output.output_capacity_mw.sum() == 0
223
            and sum_input.input_capacity_mw.sum() == 0
224
        ):
225
            logger.info(
226
                f"No capacity for carrier '{carrier}' needed to be "
227
                f"distributed. Everything is fine"
228
            )
229
230
        elif (
231
            sum_input.input_capacity_mw.sum() > 0
232
            and sum_output.output_capacity_mw.sum() == 0
233
        ):
234
            logger.info(
235
                f"Error: Capacity for carrier '{carrier}' was not distributed"
236
                f" at all!"
237
            )
238
239
        elif (
240
            sum_output.output_capacity_mw.sum() > 0
241
            and sum_input.input_capacity_mw.sum() == 0
242
        ):
243
            logger.info(
244
                f"Error: Eventhough no input capacity was provided for carrier"
245
                f" '{carrier}' a capacity got distributed!"
246
            )
247
248
        else:
249
            sum_input["error"] = (
250
                (sum_output.output_capacity_mw - sum_input.input_capacity_mw)
251
                / sum_input.input_capacity_mw
252
            ) * 100
253
            g = sum_input["error"].values[0]
254
255
            logger.info(f"{carrier}: " + str(round(g, 2)) + " %")
256
257
    # Section to check loads
258
259
    logger.info(
260
        "For German electricity loads the following deviations between the"
261
        " input and output can be observed:"
262
    )
263
264
    output_demand = db.select_dataframe(
265
        """SELECT a.scn_name, a.carrier,  SUM((SELECT SUM(p)
266
        FROM UNNEST(b.p_set) p))/1000000::numeric as load_twh
267
            FROM grid.egon_etrago_load a
268
            JOIN grid.egon_etrago_load_timeseries b
269
            ON (a.load_id = b.load_id)
270
            JOIN grid.egon_etrago_bus c
271
            ON (a.bus=c.bus_id)
272
            AND b.scn_name = 'eGon2035'
273
            AND a.scn_name = 'eGon2035'
274
            AND a.carrier = 'AC'
275
            AND c.scn_name= 'eGon2035'
276
            AND c.country='DE'
277
            GROUP BY (a.scn_name, a.carrier);
278
279
    """,
280
        warning=False,
281
    )["load_twh"].values[0]
282
283
    input_cts_ind = db.select_dataframe(
284
        """SELECT scenario,
285
         SUM(demand::numeric/1000000) as demand_mw_regio_cts_ind
286
            FROM demand.egon_demandregio_cts_ind
287
            WHERE scenario= 'eGon2035'
288
            AND year IN ('2035')
289
            GROUP BY (scenario);
290
291
        """,
292
        warning=False,
293
    )["demand_mw_regio_cts_ind"].values[0]
294
295
    input_hh = db.select_dataframe(
296
        """SELECT scenario, SUM(demand::numeric/1000000) as demand_mw_regio_hh
297
            FROM demand.egon_demandregio_hh
298
            WHERE scenario= 'eGon2035'
299
            AND year IN ('2035')
300
            GROUP BY (scenario);
301
        """,
302
        warning=False,
303
    )["demand_mw_regio_hh"].values[0]
304
305
    input_demand = input_hh + input_cts_ind
306
307
    e = round((output_demand - input_demand) / input_demand, 2) * 100
308
309
    logger.info(f"electricity demand: {e} %")
310
311
312
def etrago_eGon2035_heat():
313
    """Execute basic sanity checks.
314
315
    Returns print statements as sanity checks for the heat sector in
316
    the eGon2035 scenario.
317
318
    Parameters
319
    ----------
320
    None
321
322
    Returns
323
    -------
324
    None
325
    """
326
327
    # Check input and output values for the carriers "other_non_renewable",
328
    # "other_renewable", "reservoir", "run_of_river" and "oil"
329
330
    scn = "eGon2035"
331
332
    # Section to check generator capacities
333
    logger.info(f"Sanity checks for scenario {scn}")
334
    logger.info(
335
        "For German heat demands the following deviations between the inputs"
336
        " and outputs can be observed:"
337
    )
338
339
    # Sanity checks for heat demand
340
341
    output_heat_demand = db.select_dataframe(
342
        """SELECT a.scn_name,
343
          (SUM(
344
          (SELECT SUM(p) FROM UNNEST(b.p_set) p))/1000000)::numeric as load_twh
345
            FROM grid.egon_etrago_load a
346
            JOIN grid.egon_etrago_load_timeseries b
347
            ON (a.load_id = b.load_id)
348
            JOIN grid.egon_etrago_bus c
349
            ON (a.bus=c.bus_id)
350
            AND b.scn_name = 'eGon2035'
351
            AND a.scn_name = 'eGon2035'
352
            AND c.scn_name= 'eGon2035'
353
            AND c.country='DE'
354
            AND a.carrier IN ('rural_heat', 'central_heat')
355
            GROUP BY (a.scn_name);
356
        """,
357
        warning=False,
358
    )["load_twh"].values[0]
359
360
    input_heat_demand = db.select_dataframe(
361
        """SELECT scenario, SUM(demand::numeric/1000000) as demand_mw_peta_heat
362
            FROM demand.egon_peta_heat
363
            WHERE scenario= 'eGon2035'
364
            GROUP BY (scenario);
365
        """,
366
        warning=False,
367
    )["demand_mw_peta_heat"].values[0]
368
369
    e_demand = (
370
        round((output_heat_demand - input_heat_demand) / input_heat_demand, 2)
371
        * 100
372
    )
373
374
    logger.info(f"heat demand: {e_demand} %")
375
376
    # Sanity checks for heat supply
377
378
    logger.info(
379
        "For German heat supplies the following deviations between the inputs "
380
        "and outputs can be observed:"
381
    )
382
383
    # Comparison for central heat pumps
384
    heat_pump_input = db.select_dataframe(
385
        """SELECT carrier, SUM(capacity::numeric) as Urban_central_heat_pump_mw
386
            FROM supply.egon_scenario_capacities
387
            WHERE carrier= 'urban_central_heat_pump'
388
            AND scenario_name IN ('eGon2035')
389
            GROUP BY (carrier);
390
        """,
391
        warning=False,
392
    )["urban_central_heat_pump_mw"].values[0]
393
394
    heat_pump_output = db.select_dataframe(
395
        """SELECT carrier, SUM(p_nom::numeric) as Central_heat_pump_mw
396
            FROM grid.egon_etrago_link
397
            WHERE carrier= 'central_heat_pump'
398
            AND scn_name IN ('eGon2035')
399
            GROUP BY (carrier);
400
    """,
401
        warning=False,
402
    )["central_heat_pump_mw"].values[0]
403
404
    e_heat_pump = (
405
        round((heat_pump_output - heat_pump_input) / heat_pump_output, 2) * 100
406
    )
407
408
    logger.info(f"'central_heat_pump': {e_heat_pump} % ")
409
410
    # Comparison for residential heat pumps
411
412
    input_residential_heat_pump = db.select_dataframe(
413
        """SELECT carrier, SUM(capacity::numeric) as residential_heat_pump_mw
414
            FROM supply.egon_scenario_capacities
415
            WHERE carrier= 'residential_rural_heat_pump'
416
            AND scenario_name IN ('eGon2035')
417
            GROUP BY (carrier);
418
        """,
419
        warning=False,
420
    )["residential_heat_pump_mw"].values[0]
421
422
    output_residential_heat_pump = db.select_dataframe(
423
        """SELECT carrier, SUM(p_nom::numeric) as rural_heat_pump_mw
424
            FROM grid.egon_etrago_link
425
            WHERE carrier= 'rural_heat_pump'
426
            AND scn_name IN ('eGon2035')
427
            GROUP BY (carrier);
428
    """,
429
        warning=False,
430
    )["rural_heat_pump_mw"].values[0]
431
432
    e_residential_heat_pump = (
433
        round(
434
            (output_residential_heat_pump - input_residential_heat_pump)
435
            / input_residential_heat_pump,
436
            2,
437
        )
438
        * 100
439
    )
440
    logger.info(f"'residential heat pumps': {e_residential_heat_pump} %")
441
442
    # Comparison for resistive heater
443
    resistive_heater_input = db.select_dataframe(
444
        """SELECT carrier,
445
         SUM(capacity::numeric) as Urban_central_resistive_heater_MW
446
            FROM supply.egon_scenario_capacities
447
            WHERE carrier= 'urban_central_resistive_heater'
448
            AND scenario_name IN ('eGon2035')
449
            GROUP BY (carrier);
450
        """,
451
        warning=False,
452
    )["urban_central_resistive_heater_mw"].values[0]
453
454
    resistive_heater_output = db.select_dataframe(
455
        """SELECT carrier, SUM(p_nom::numeric) as central_resistive_heater_MW
456
            FROM grid.egon_etrago_link
457
            WHERE carrier= 'central_resistive_heater'
458
            AND scn_name IN ('eGon2035')
459
            GROUP BY (carrier);
460
        """,
461
        warning=False,
462
    )["central_resistive_heater_mw"].values[0]
463
464
    e_resistive_heater = (
465
        round(
466
            (resistive_heater_output - resistive_heater_input)
467
            / resistive_heater_input,
468
            2,
469
        )
470
        * 100
471
    )
472
473
    logger.info(f"'resistive heater': {e_resistive_heater} %")
474
475
    # Comparison for solar thermal collectors
476
477
    input_solar_thermal = db.select_dataframe(
478
        """SELECT carrier, SUM(capacity::numeric) as solar_thermal_collector_mw
479
            FROM supply.egon_scenario_capacities
480
            WHERE carrier= 'urban_central_solar_thermal_collector'
481
            AND scenario_name IN ('eGon2035')
482
            GROUP BY (carrier);
483
        """,
484
        warning=False,
485
    )["solar_thermal_collector_mw"].values[0]
486
487
    output_solar_thermal = db.select_dataframe(
488
        """SELECT carrier, SUM(p_nom::numeric) as solar_thermal_collector_mw
489
            FROM grid.egon_etrago_generator
490
            WHERE carrier= 'solar_thermal_collector'
491
            AND scn_name IN ('eGon2035')
492
            GROUP BY (carrier);
493
        """,
494
        warning=False,
495
    )["solar_thermal_collector_mw"].values[0]
496
497
    e_solar_thermal = (
498
        round(
499
            (output_solar_thermal - input_solar_thermal) / input_solar_thermal,
500
            2,
501
        )
502
        * 100
503
    )
504
    logger.info(f"'solar thermal collector': {e_solar_thermal} %")
505
506
    # Comparison for geothermal
507
508
    input_geo_thermal = db.select_dataframe(
509
        """SELECT carrier,
510
         SUM(capacity::numeric) as Urban_central_geo_thermal_MW
511
            FROM supply.egon_scenario_capacities
512
            WHERE carrier= 'urban_central_geo_thermal'
513
            AND scenario_name IN ('eGon2035')
514
            GROUP BY (carrier);
515
        """,
516
        warning=False,
517
    )["urban_central_geo_thermal_mw"].values[0]
518
519
    output_geo_thermal = db.select_dataframe(
520
        """SELECT carrier, SUM(p_nom::numeric) as geo_thermal_MW
521
            FROM grid.egon_etrago_generator
522
            WHERE carrier= 'geo_thermal'
523
            AND scn_name IN ('eGon2035')
524
            GROUP BY (carrier);
525
    """,
526
        warning=False,
527
    )["geo_thermal_mw"].values[0]
528
529
    e_geo_thermal = (
530
        round((output_geo_thermal - input_geo_thermal) / input_geo_thermal, 2)
531
        * 100
532
    )
533
    logger.info(f"'geothermal': {e_geo_thermal} %")
534
535
536
def residential_electricity_annual_sum(rtol=1e-5):
537
    """Sanity check for dataset electricity_demand_timeseries :
538
    Demand_Building_Assignment
539
540
    Aggregate the annual demand of all census cells at NUTS3 to compare
541
    with initial scaling parameters from DemandRegio.
542
    """
543
544
    df_nuts3_annual_sum = db.select_dataframe(
545
        sql="""
546
        SELECT dr.nuts3, dr.scenario, dr.demand_regio_sum, profiles.profile_sum
547
        FROM (
548
            SELECT scenario, SUM(demand) AS profile_sum, vg250_nuts3
549
            FROM demand.egon_demandregio_zensus_electricity AS egon,
550
             boundaries.egon_map_zensus_vg250 AS boundaries
551
            Where egon.zensus_population_id = boundaries.zensus_population_id
552
            AND sector = 'residential'
553
            GROUP BY vg250_nuts3, scenario
554
            ) AS profiles
555
        JOIN (
556
            SELECT nuts3, scenario, sum(demand) AS demand_regio_sum
557
            FROM demand.egon_demandregio_hh
558
            GROUP BY year, scenario, nuts3
559
              ) AS dr
560
        ON profiles.vg250_nuts3 = dr.nuts3 and profiles.scenario  = dr.scenario
561
        """
562
    )
563
564
    np.testing.assert_allclose(
565
        actual=df_nuts3_annual_sum["profile_sum"],
566
        desired=df_nuts3_annual_sum["demand_regio_sum"],
567
        rtol=rtol,
568
        verbose=False,
569
    )
570
571
    logger.info(
572
        "Aggregated annual residential electricity demand"
573
        " matches with DemandRegio at NUTS-3."
574
    )
575
576
577
def residential_electricity_hh_refinement(rtol=1e-5):
578
    """Sanity check for dataset electricity_demand_timeseries :
579
    Household Demands
580
581
    Check sum of aggregated household types after refinement method
582
    was applied and compare it to the original census values."""
583
584
    df_refinement = db.select_dataframe(
585
        sql="""
586
        SELECT refined.nuts3, refined.characteristics_code,
587
                refined.sum_refined::int, census.sum_census::int
588
        FROM(
589
            SELECT nuts3, characteristics_code, SUM(hh_10types) as sum_refined
590
            FROM society.egon_destatis_zensus_household_per_ha_refined
591
            GROUP BY nuts3, characteristics_code)
592
            AS refined
593
        JOIN(
594
            SELECT t.nuts3, t.characteristics_code, sum(orig) as sum_census
595
            FROM(
596
                SELECT nuts3, cell_id, characteristics_code,
597
                        sum(DISTINCT(hh_5types))as orig
598
                FROM society.egon_destatis_zensus_household_per_ha_refined
599
                GROUP BY cell_id, characteristics_code, nuts3) AS t
600
            GROUP BY t.nuts3, t.characteristics_code    ) AS census
601
        ON refined.nuts3 = census.nuts3
602
        AND refined.characteristics_code = census.characteristics_code
603
    """
604
    )
605
606
    np.testing.assert_allclose(
607
        actual=df_refinement["sum_refined"],
608
        desired=df_refinement["sum_census"],
609
        rtol=rtol,
610
        verbose=False,
611
    )
612
613
    logger.info("All Aggregated household types match at NUTS-3.")
614
615
616
def cts_electricity_demand_share(rtol=1e-5):
617
    """Sanity check for dataset electricity_demand_timeseries :
618
    CtsBuildings
619
620
    Check sum of aggregated cts electricity demand share which equals to one
621
    for every substation as the substation profile is linearly disaggregated
622
    to all buildings."""
623
624
    with db.session_scope() as session:
625
        cells_query = session.query(EgonCtsElectricityDemandBuildingShare)
626
627
    df_demand_share = pd.read_sql(
628
        cells_query.statement, cells_query.session.bind, index_col=None
629
    )
630
631
    np.testing.assert_allclose(
632
        actual=df_demand_share.groupby(["bus_id", "scenario"])[
633
            "profile_share"
634
        ].sum(),
635
        desired=1,
636
        rtol=rtol,
637
        verbose=False,
638
    )
639
640
    logger.info("The aggregated demand shares equal to one!.")
641
642
643
def cts_heat_demand_share(rtol=1e-5):
644
    """Sanity check for dataset electricity_demand_timeseries
645
    : CtsBuildings
646
647
    Check sum of aggregated cts heat demand share which equals to one
648
    for every substation as the substation profile is linearly disaggregated
649
    to all buildings."""
650
651
    with db.session_scope() as session:
652
        cells_query = session.query(EgonCtsHeatDemandBuildingShare)
653
654
    df_demand_share = pd.read_sql(
655
        cells_query.statement, cells_query.session.bind, index_col=None
656
    )
657
658
    np.testing.assert_allclose(
659
        actual=df_demand_share.groupby(["bus_id", "scenario"])[
660
            "profile_share"
661
        ].sum(),
662
        desired=1,
663
        rtol=rtol,
664
        verbose=False,
665
    )
666
667
    logger.info("The aggregated demand shares equal to one!.")
668
669
670
def sanitycheck_emobility_mit():
671
    """Execute sanity checks for eMobility: motorized individual travel
672
673
    Checks data integrity for eGon2035, eGon2035_lowflex and eGon100RE scenario
674
    using assertions:
675
      1. Allocated EV numbers and EVs allocated to grid districts
676
      2. Trip data (original inout data from simBEV)
677
      3. Model data in eTraGo PF tables (grid.egon_etrago_*)
678
679
    Parameters
680
    ----------
681
    None
682
683
    Returns
684
    -------
685
    None
686
    """
687
688
    def check_ev_allocation():
689
        # Get target number for scenario
690
        ev_count_target = scenario_variation_parameters["ev_count"]
691
        print(f"  Target count: {str(ev_count_target)}")
692
693
        # Get allocated numbers
694
        ev_counts_dict = {}
695
        with db.session_scope() as session:
696
            for table, level in zip(
697
                [
698
                    EgonEvCountMvGridDistrict,
699
                    EgonEvCountMunicipality,
700
                    EgonEvCountRegistrationDistrict,
701
                ],
702
                ["Grid District", "Municipality", "Registration District"],
703
            ):
704
                query = session.query(
705
                    func.sum(
706
                        table.bev_mini
707
                        + table.bev_medium
708
                        + table.bev_luxury
709
                        + table.phev_mini
710
                        + table.phev_medium
711
                        + table.phev_luxury
712
                    ).label("ev_count")
713
                ).filter(
714
                    table.scenario == scenario_name,
715
                    table.scenario_variation == scenario_var_name,
716
                )
717
718
                ev_counts = pd.read_sql(
719
                    query.statement, query.session.bind, index_col=None
720
                )
721
                ev_counts_dict[level] = ev_counts.iloc[0].ev_count
722
                print(
723
                    f"    Count table: Total count for level {level} "
724
                    f"(table: {table.__table__}): "
725
                    f"{str(ev_counts_dict[level])}"
726
                )
727
728
        # Compare with scenario target (only if not in testmode)
729
        if TESTMODE_OFF:
730
            for level, count in ev_counts_dict.items():
731
                np.testing.assert_allclose(
732
                    count,
733
                    ev_count_target,
734
                    rtol=0.0001,
735
                    err_msg=f"EV numbers in {level} seems to be flawed.",
736
                )
737
        else:
738
            print("    Testmode is on, skipping sanity check...")
739
740
        # Get allocated EVs in grid districts
741
        with db.session_scope() as session:
742
            query = session.query(
743
                func.count(EgonEvMvGridDistrict.egon_ev_pool_ev_id).label(
744
                    "ev_count"
745
                ),
746
            ).filter(
747
                EgonEvMvGridDistrict.scenario == scenario_name,
748
                EgonEvMvGridDistrict.scenario_variation == scenario_var_name,
749
            )
750
        ev_count_alloc = (
751
            pd.read_sql(query.statement, query.session.bind, index_col=None)
752
            .iloc[0]
753
            .ev_count
754
        )
755
        print(
756
            f"    EVs allocated to Grid Districts "
757
            f"(table: {EgonEvMvGridDistrict.__table__}) total count: "
758
            f"{str(ev_count_alloc)}"
759
        )
760
761
        # Compare with scenario target (only if not in testmode)
762
        if TESTMODE_OFF:
763
            np.testing.assert_allclose(
764
                ev_count_alloc,
765
                ev_count_target,
766
                rtol=0.0001,
767
                err_msg=(
768
                    "EV numbers allocated to Grid Districts seems to be flawed."
769
                ),
770
            )
771
        else:
772
            print("    Testmode is on, skipping sanity check...")
773
774
        return ev_count_alloc
775
776
    def check_trip_data():
777
        # Check if trips start at timestep 0 and have a max. of 35040 steps
778
        # (8760h in 15min steps)
779
        print("  Checking timeranges...")
780
        with db.session_scope() as session:
781
            query = session.query(
782
                func.count(EgonEvTrip.event_id).label("cnt")
783
            ).filter(
784
                or_(
785
                    and_(
786
                        EgonEvTrip.park_start > 0,
787
                        EgonEvTrip.simbev_event_id == 0,
788
                    ),
789
                    EgonEvTrip.park_end
790
                    > (60 / int(meta_run_config.stepsize)) * 8760,
791
                ),
792
                EgonEvTrip.scenario == scenario_name,
793
            )
794
        invalid_trips = pd.read_sql(
795
            query.statement, query.session.bind, index_col=None
796
        )
797
        np.testing.assert_equal(
798
            invalid_trips.iloc[0].cnt,
799
            0,
800
            err_msg=(
801
                f"{str(invalid_trips.iloc[0].cnt)} trips in table "
802
                f"{EgonEvTrip.__table__} have invalid timesteps."
803
            ),
804
        )
805
806
        # Check if charging demand can be covered by available charging energy
807
        # while parking
808
        print("  Compare charging demand with available power...")
809
        with db.session_scope() as session:
810
            query = session.query(
811
                func.count(EgonEvTrip.event_id).label("cnt")
812
            ).filter(
813
                func.round(
814
                    cast(
815
                        (EgonEvTrip.park_end - EgonEvTrip.park_start + 1)
816
                        * EgonEvTrip.charging_capacity_nominal
817
                        * (int(meta_run_config.stepsize) / 60),
818
                        Numeric,
819
                    ),
820
                    3,
821
                )
822
                < cast(EgonEvTrip.charging_demand, Numeric),
823
                EgonEvTrip.scenario == scenario_name,
824
            )
825
        invalid_trips = pd.read_sql(
826
            query.statement, query.session.bind, index_col=None
827
        )
828
        np.testing.assert_equal(
829
            invalid_trips.iloc[0].cnt,
830
            0,
831
            err_msg=(
832
                f"In {str(invalid_trips.iloc[0].cnt)} trips (table: "
833
                f"{EgonEvTrip.__table__}) the charging demand cannot be "
834
                f"covered by available charging power."
835
            ),
836
        )
837
838
    def check_model_data():
839
        # Check if model components were fully created
840
        print("  Check if all model components were created...")
841
        # Get MVGDs which got EV allocated
842
        with db.session_scope() as session:
843
            query = (
844
                session.query(
845
                    EgonEvMvGridDistrict.bus_id,
846
                )
847
                .filter(
848
                    EgonEvMvGridDistrict.scenario == scenario_name,
849
                    EgonEvMvGridDistrict.scenario_variation
850
                    == scenario_var_name,
851
                )
852
                .group_by(EgonEvMvGridDistrict.bus_id)
853
            )
854
        mvgds_with_ev = (
855
            pd.read_sql(query.statement, query.session.bind, index_col=None)
856
            .bus_id.sort_values()
857
            .to_list()
858
        )
859
860
        # Load model components
861
        with db.session_scope() as session:
862
            query = (
863
                session.query(
864
                    EgonPfHvLink.bus0.label("mvgd_bus_id"),
865
                    EgonPfHvLoad.bus.label("emob_bus_id"),
866
                    EgonPfHvLoad.load_id.label("load_id"),
867
                    EgonPfHvStore.store_id.label("store_id"),
868
                )
869
                .select_from(EgonPfHvLoad, EgonPfHvStore)
870
                .join(
871
                    EgonPfHvLoadTimeseries,
872
                    EgonPfHvLoadTimeseries.load_id == EgonPfHvLoad.load_id,
873
                )
874
                .join(
875
                    EgonPfHvStoreTimeseries,
876
                    EgonPfHvStoreTimeseries.store_id == EgonPfHvStore.store_id,
877
                )
878
                .filter(
879
                    EgonPfHvLoad.carrier == "land transport EV",
880
                    EgonPfHvLoad.scn_name == scenario_name,
881
                    EgonPfHvLoadTimeseries.scn_name == scenario_name,
882
                    EgonPfHvStore.carrier == "battery storage",
883
                    EgonPfHvStore.scn_name == scenario_name,
884
                    EgonPfHvStoreTimeseries.scn_name == scenario_name,
885
                    EgonPfHvLink.scn_name == scenario_name,
886
                    EgonPfHvLink.bus1 == EgonPfHvLoad.bus,
887
                    EgonPfHvLink.bus1 == EgonPfHvStore.bus,
888
                )
889
            )
890
        model_components = pd.read_sql(
891
            query.statement, query.session.bind, index_col=None
892
        )
893
894
        # Check number of buses with model components connected
895
        mvgd_buses_with_ev = model_components.loc[
896
            model_components.mvgd_bus_id.isin(mvgds_with_ev)
897
        ]
898
        np.testing.assert_equal(
899
            len(mvgds_with_ev),
900
            len(mvgd_buses_with_ev),
901
            err_msg=(
902
                f"Number of Grid Districts with connected model components "
903
                f"({str(len(mvgd_buses_with_ev))} in tables egon_etrago_*) "
904
                f"differ from number of Grid Districts that got EVs "
905
                f"allocated ({len(mvgds_with_ev)} in table "
906
                f"{EgonEvMvGridDistrict.__table__})."
907
            ),
908
        )
909
910
        # Check if all required components exist (if no id is NaN)
911
        np.testing.assert_equal(
912
            model_components.drop_duplicates().isna().any().any(),
913
            False,
914
            err_msg=(
915
                f"Some components are missing (see True values): "
916
                f"{model_components.drop_duplicates().isna().any()}"
917
            ),
918
        )
919
920
        # Get all model timeseries
921
        print("  Loading model timeseries...")
922
        # Get all model timeseries
923
        model_ts_dict = {
924
            "Load": {
925
                "carrier": "land transport EV",
926
                "table": EgonPfHvLoad,
927
                "table_ts": EgonPfHvLoadTimeseries,
928
                "column_id": "load_id",
929
                "columns_ts": ["p_set"],
930
                "ts": None,
931
            },
932
            "Link": {
933
                "carrier": "BEV charger",
934
                "table": EgonPfHvLink,
935
                "table_ts": EgonPfHvLinkTimeseries,
936
                "column_id": "link_id",
937
                "columns_ts": ["p_max_pu"],
938
                "ts": None,
939
            },
940
            "Store": {
941
                "carrier": "battery storage",
942
                "table": EgonPfHvStore,
943
                "table_ts": EgonPfHvStoreTimeseries,
944
                "column_id": "store_id",
945
                "columns_ts": ["e_min_pu", "e_max_pu"],
946
                "ts": None,
947
            },
948
        }
949
950
        with db.session_scope() as session:
951
            for node, attrs in model_ts_dict.items():
952
                print(f"    Loading {node} timeseries...")
953
                subquery = (
954
                    session.query(
955
                        getattr(attrs["table"], attrs["column_id"])
956
                    )
957
                    .filter(attrs["table"].carrier == attrs["carrier"])
958
                    .filter(attrs["table"].scn_name == scenario_name)
959
                    .subquery()
960
                )
961
962
                cols = [
963
                    getattr(attrs["table_ts"], c) for c in attrs["columns_ts"]
964
                ]
965
                query = session.query(
966
                    getattr(attrs["table_ts"], attrs["column_id"]), *cols
967
                ).filter(
968
                    getattr(attrs["table_ts"], attrs["column_id"]).in_(
969
                        subquery
970
                    ),
971
                    attrs["table_ts"].scn_name == scenario_name,
972
                )
973
                attrs["ts"] = pd.read_sql(
974
                    query.statement,
975
                    query.session.bind,
976
                    index_col=attrs["column_id"],
977
                )
978
979
        # Check if all timeseries have 8760 steps
980
        print("    Checking timeranges...")
981
        for node, attrs in model_ts_dict.items():
982
            for col in attrs["columns_ts"]:
983
                ts = attrs["ts"]
984
                invalid_ts = ts.loc[ts[col].apply(lambda _: len(_)) != 8760][
985
                    col
986
                ].apply(len)
987
                np.testing.assert_equal(
988
                    len(invalid_ts),
989
                    0,
990
                    err_msg=(
991
                        f"{str(len(invalid_ts))} rows in timeseries do not "
992
                        f"have 8760 timesteps. Table: "
993
                        f"{attrs['table_ts'].__table__}, Column: {col}, IDs: "
994
                        f"{str(list(invalid_ts.index))}"
995
                    ),
996
                )
997
998
        # Compare total energy demand in model with some approximate values
999
        # (per EV: 14,000 km/a, 0.17 kWh/km)
1000
        print("  Checking energy demand in model...")
1001
        total_energy_model = (
1002
            model_ts_dict["Load"]["ts"].p_set.apply(lambda _: sum(_)).sum()
1003
            / 1e6
1004
        )
1005
        print(f"    Total energy amount in model: {total_energy_model} TWh")
1006
        total_energy_scenario_approx = ev_count_alloc * 14000 * 0.17 / 1e9
1007
        print(
1008
            f"    Total approximated energy amount in scenario: "
1009
            f"{total_energy_scenario_approx} TWh"
1010
        )
1011
        np.testing.assert_allclose(
1012
            total_energy_model,
1013
            total_energy_scenario_approx,
1014
            rtol=0.1,
1015
            err_msg=(
1016
                "The total energy amount in the model deviates heavily "
1017
                "from the approximated value for current scenario."
1018
            ),
1019
        )
1020
1021
        # Compare total storage capacity
1022
        print("  Checking storage capacity...")
1023
        # Load storage capacities from model
1024
        with db.session_scope() as session:
1025
            query = session.query(
1026
                func.sum(EgonPfHvStore.e_nom).label("e_nom")
1027
            ).filter(
1028
                EgonPfHvStore.scn_name == scenario_name,
1029
                EgonPfHvStore.carrier == "battery storage",
1030
            )
1031
        storage_capacity_model = (
1032
            pd.read_sql(
1033
                query.statement, query.session.bind, index_col=None
1034
            ).e_nom.sum()
1035
            / 1e3
1036
        )
1037
        print(
1038
            f"    Total storage capacity ({EgonPfHvStore.__table__}): "
1039
            f"{round(storage_capacity_model, 1)} GWh"
1040
        )
1041
1042
        # Load occurences of each EV
1043
        with db.session_scope() as session:
1044
            query = (
1045
                session.query(
1046
                    EgonEvMvGridDistrict.bus_id,
1047
                    EgonEvPool.type,
1048
                    func.count(EgonEvMvGridDistrict.egon_ev_pool_ev_id).label(
1049
                        "count"
1050
                    ),
1051
                )
1052
                .join(
1053
                    EgonEvPool,
1054
                    EgonEvPool.ev_id
1055
                    == EgonEvMvGridDistrict.egon_ev_pool_ev_id,
1056
                )
1057
                .filter(
1058
                    EgonEvMvGridDistrict.scenario == scenario_name,
1059
                    EgonEvMvGridDistrict.scenario_variation
1060
                    == scenario_var_name,
1061
                    EgonEvPool.scenario == scenario_name,
1062
                )
1063
                .group_by(EgonEvMvGridDistrict.bus_id, EgonEvPool.type)
1064
            )
1065
        count_per_ev_all = pd.read_sql(
1066
            query.statement, query.session.bind, index_col="bus_id"
1067
        )
1068
        count_per_ev_all["bat_cap"] = count_per_ev_all.type.map(
1069
            meta_tech_data.battery_capacity
1070
        )
1071
        count_per_ev_all["bat_cap_total_MWh"] = (
1072
            count_per_ev_all["count"] * count_per_ev_all.bat_cap / 1e3
1073
        )
1074
        storage_capacity_simbev = count_per_ev_all.bat_cap_total_MWh.div(
1075
            1e3
1076
        ).sum()
1077
        print(
1078
            f"    Total storage capacity (simBEV): "
1079
            f"{round(storage_capacity_simbev, 1)} GWh"
1080
        )
1081
1082
        np.testing.assert_allclose(
1083
            storage_capacity_model,
1084
            storage_capacity_simbev,
1085
            rtol=0.01,
1086
            err_msg=(
1087
                "The total storage capacity in the model deviates heavily "
1088
                "from the input data provided by simBEV for current scenario."
1089
            ),
1090
        )
1091
1092
        # Check SoC storage constraint: e_min_pu < e_max_pu for all timesteps
1093
        print("  Validating SoC constraints...")
1094
        stores_with_invalid_soc = []
1095
        for idx, row in model_ts_dict["Store"]["ts"].iterrows():
1096
            ts = row[["e_min_pu", "e_max_pu"]]
1097
            x = np.array(ts.e_min_pu) > np.array(ts.e_max_pu)
1098
            if x.any():
1099
                stores_with_invalid_soc.append(idx)
1100
1101
        np.testing.assert_equal(
1102
            len(stores_with_invalid_soc),
1103
            0,
1104
            err_msg=(
1105
                f"The store constraint e_min_pu < e_max_pu does not apply "
1106
                f"for some storages in {EgonPfHvStoreTimeseries.__table__}. "
1107
                f"Invalid store_ids: {stores_with_invalid_soc}"
1108
            ),
1109
        )
1110
1111
    def check_model_data_lowflex_eGon2035():
1112
        # TODO: Add eGon100RE_lowflex
1113
        print("")
1114
        print("SCENARIO: eGon2035_lowflex")
1115
1116
        # Compare driving load and charging load
1117
        print("  Loading eGon2035 model timeseries: driving load...")
1118
        with db.session_scope() as session:
1119
            query = (
1120
                session.query(
1121
                    EgonPfHvLoad.load_id,
1122
                    EgonPfHvLoadTimeseries.p_set,
1123
                )
1124
                .join(
1125
                    EgonPfHvLoadTimeseries,
1126
                    EgonPfHvLoadTimeseries.load_id == EgonPfHvLoad.load_id,
1127
                )
1128
                .filter(
1129
                    EgonPfHvLoad.carrier == "land transport EV",
1130
                    EgonPfHvLoad.scn_name == "eGon2035",
1131
                    EgonPfHvLoadTimeseries.scn_name == "eGon2035",
1132
                )
1133
            )
1134
        model_driving_load = pd.read_sql(
1135
            query.statement, query.session.bind, index_col=None
1136
        )
1137
        driving_load = np.array(model_driving_load.p_set.to_list()).sum(axis=0)
1138
1139
        print(
1140
            "  Loading eGon2035_lowflex model timeseries: dumb charging "
1141
            "load..."
1142
        )
1143
        with db.session_scope() as session:
1144
            query = (
1145
                session.query(
1146
                    EgonPfHvLoad.load_id,
1147
                    EgonPfHvLoadTimeseries.p_set,
1148
                )
1149
                .join(
1150
                    EgonPfHvLoadTimeseries,
1151
                    EgonPfHvLoadTimeseries.load_id == EgonPfHvLoad.load_id,
1152
                )
1153
                .filter(
1154
                    EgonPfHvLoad.carrier == "land transport EV",
1155
                    EgonPfHvLoad.scn_name == "eGon2035_lowflex",
1156
                    EgonPfHvLoadTimeseries.scn_name == "eGon2035_lowflex",
1157
                )
1158
            )
1159
        model_charging_load_lowflex = pd.read_sql(
1160
            query.statement, query.session.bind, index_col=None
1161
        )
1162
        charging_load = np.array(
1163
            model_charging_load_lowflex.p_set.to_list()
1164
        ).sum(axis=0)
1165
1166
        # Ratio of driving and charging load should be 0.9 due to charging
1167
        # efficiency
1168
        print("  Compare cumulative loads...")
1169
        print(f"    Driving load (eGon2035): {driving_load.sum() / 1e6} TWh")
1170
        print(
1171
            f"    Dumb charging load (eGon2035_lowflex): "
1172
            f"{charging_load.sum() / 1e6} TWh"
1173
        )
1174
        driving_load_theoretical = (
1175
            float(meta_run_config.eta_cp) * charging_load.sum()
1176
        )
1177
        np.testing.assert_allclose(
1178
            driving_load.sum(),
1179
            driving_load_theoretical,
1180
            rtol=0.01,
1181
            err_msg=(
1182
                f"The driving load (eGon2035) deviates by more than 1% "
1183
                f"from the theoretical driving load calculated from charging "
1184
                f"load (eGon2035_lowflex) with an efficiency of "
1185
                f"{float(meta_run_config.eta_cp)}."
1186
            ),
1187
        )
1188
1189
    print("=====================================================")
1190
    print("=== SANITY CHECKS FOR MOTORIZED INDIVIDUAL TRAVEL ===")
1191
    print("=====================================================")
1192
1193
    for scenario_name in ["eGon2035", "eGon100RE"]:
1194
        scenario_var_name = DATASET_CFG["scenario"]["variation"][scenario_name]
1195
1196
        print("")
1197
        print(f"SCENARIO: {scenario_name}, VARIATION: {scenario_var_name}")
1198
1199
        # Load scenario params for scenario and scenario variation
1200
        scenario_variation_parameters = get_sector_parameters(
1201
            "mobility", scenario=scenario_name
1202
        )["motorized_individual_travel"][scenario_var_name]
1203
1204
        # Load simBEV run config and tech data
1205
        meta_run_config = read_simbev_metadata_file(
1206
            scenario_name, "config"
1207
        ).loc["basic"]
1208
        meta_tech_data = read_simbev_metadata_file(scenario_name, "tech_data")
1209
1210
        print("")
1211
        print("Checking EV counts...")
1212
        ev_count_alloc = check_ev_allocation()
1213
1214
        print("")
1215
        print("Checking trip data...")
1216
        check_trip_data()
1217
1218
        print("")
1219
        print("Checking model data...")
1220
        check_model_data()
1221
1222
    print("")
1223
    check_model_data_lowflex_eGon2035()
1224
1225
    print("=====================================================")
1226