{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Advanced features\n",
"\n",
"* File name: advanced_features.ipynb\n",
"* Last edited: 2020-06-25\n",
"* Created by: Stefan Bruche (TU Berlin)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"```python\n",
"import copy\n",
"import pandas as pd\n",
"import aristopy as ar\n",
"\n",
"# Get data from csv-file \n",
"data = pd.read_csv('testdata.csv', sep=';', decimal='.', index_col=[0])\n",
"\n",
"# USER-INPUT: \"hpts\" (hours per time step)\n",
"hpts = 4\n",
"\n",
"# Aggregate data with \"hours_per_time_step\" for simplification of this notebook\n",
"c_elec = ar.utils.aggregate_time_series(data['c_elec [EUR/MWh]'], hpts, 'mean')\n",
"q_demand = ar.utils.aggregate_time_series(data['Q_demand [MWh]'], hpts, 'sum')\n",
"\n",
"# Create energy system instance\n",
"es = ar.EnergySystem(\n",
" number_of_time_steps=8760//hpts, hours_per_time_step=hpts,\n",
" interest_rate=0.05, economic_lifetime=20, logging=ar.Logger(\n",
" default_log_handler='file', default_log_level='INFO', \n",
" write_screen_output_to_logfile=True))\n",
"\n",
"# Add a gas source, a heat and electricity sink and a heat bus for collecting\n",
"gas_source = ar.Source(\n",
" ensys=es, name='gas_source', commodity_cost=20, outlet=ar.Flow('Fuel'))\n",
"\n",
"heat_sink = ar.Sink(\n",
" ensys=es, name='heat_sink', inlet=ar.Flow('Heat', link='heat_bus'),\n",
" commodity_rate_fix=ar.Series('q_demand', q_demand))\n",
"\n",
"elec_sink = ar.Sink(\n",
" ensys=es, name='elec_sink', inlet=ar.Flow('Elec'),\n",
" commodity_revenues=ar.Series('c_elec', c_elec))\n",
"\n",
"heat_bus = ar.Bus(\n",
" ensys=es, name='heat_bus',\n",
" inlet=ar.Flow('Heat', var_name='Heat_inlet'),\n",
" outlet=ar.Flow('Heat', var_name='Heat_outlet'))\n",
"\n",
"gas_boiler = ar.Conversion(\n",
" ensys=es, name='gas_boiler', basic_variable='Heat',\n",
" inlet=ar.Flow('Fuel', 'gas_source'), outlet=ar.Flow('Heat', 'heat_bus'),\n",
" capacity_max=250, capex_per_capacity=60e3,\n",
" user_expressions='Heat == 0.9 * Fuel')\n",
"\n",
"# Add a group of 5 gas engine units with fixed size\n",
"gas_engine = ar.Conversion(\n",
" ensys=es, name='gas_engine', basic_variable='Elec',\n",
" inlet=ar.Flow('Fuel', 'gas_source'),\n",
" outlet=[ar.Flow('Heat', 'heat_bus'), ar.Flow('Elec', 'elec_sink')],\n",
" has_existence_binary_var=True, has_operation_binary_var=True,\n",
" additional_vars=ar.Var('Load', ub=1), scalar_params={'dt': hpts},\n",
" capacity=20, min_load_rel=0.2,\n",
" user_expressions=['Load <= BI_OP',\n",
" 'Elec == (20 * Load) * dt',\n",
" 'Heat <= (17 * Load + 2.5 * BI_OP) * dt',\n",
" 'Fuel == (37 * Load + 7 * BI_OP) * dt'],\n",
" capex_per_capacity=1000e3, opex_per_capacity=12.5e3, opex_operation=6,\n",
" instances_in_group=5)\n",
"\n",
"# Add a heat storage component\n",
"heat_storage = ar.Storage(\n",
" ensys=es, name='heat_storage',\n",
" inlet=ar.Flow('Heat', link='heat_bus', var_name='Heat_charge'),\n",
" outlet=ar.Flow('Heat', link='heat_bus', var_name='Heat_discharge'),\n",
" has_existence_binary_var=True,\n",
" maximal_module_number=5, capacity_per_module=50, # max: 5 * 50 MWh\n",
" capex_per_capacity=12e3, opex_per_capacity=0.015*12e3, # 12,000 €/MWh CAPEX\n",
" self_discharge=0.01, charge_rate=1/6, discharge_rate=1/6,\n",
" opex_operation=1e-9, use_inter_period_formulation=False)\n",
"\n",
"# Use a backdoor to the pyomo model instance and add a user-defined constraint\n",
"# specifying the installed nominal thermal capacity needs to be at least 220 MW.\n",
"def minimal_installed_capacity_heat(ensys, m):\n",
" expr = 0\n",
" for comp in ensys.components.values():\n",
" if comp.name == 'gas_boiler':\n",
" expr += comp.block.CAP\n",
" elif comp.group_name == 'gas_engine':\n",
" expr += 19.5 * comp.block.BI_EX\n",
" return expr >= 220\n",
"es.add_constraint(name='min_cap_heat', has_time_set=False,\n",
" rule=minimal_installed_capacity_heat)\n",
"\n",
"\n",
"# Create a deepcopy of the original (full scale) model instance for later use\n",
"original_model = copy.deepcopy(es)\n",
"\n",
"\n",
"# Conventional optimization:\n",
"# ~~~~~~~~~~~~~~~~~~~~~~~~~~\n",
"# Run the optimization of the original model instance\n",
"es.optimize(solver='gurobi', time_limit=300, optimization_specs='mipgap=5e-3')\n",
"\n",
"\n",
"# Time series clustering:\n",
"# ~~~~~~~~~~~~~~~~~~~~~~~\n",
"# Cluster time series data to 20 typical days.\n",
"es.cluster(number_of_typical_periods=20,\n",
" number_of_time_steps_per_period=24//hpts)\n",
"\n",
"# Use aggregated data to run the optimization\n",
"es.optimize(use_clustered_data=True)\n",
"\n",
"# Plot results for a single period and scaled to full year with hourly resolution\n",
"plotter = ar.Plotter('results.json')\n",
"plotter.plot_operation('heat_storage', 'Heat', ylabel='Thermal energy [MWh]',\n",
" file_name='heat_storage', plot_single_period_with_index=7)\n",
"plotter.plot_operation('heat_bus', 'Heat', ylabel='Thermal energy [MWh]', \n",
" file_name='heat_bus', scale_to_hourly_resolution=True)\n",
"\n",
"\n",
"# Two-stage model run:\n",
"# ~~~~~~~~~~~~~~~~~~~~\n",
"# Construct (declare) a persistent model of the original (full scale) instance\n",
"original_model.declare_model(use_clustered_data=False, declare_persistent=True)\n",
"\n",
"# Cluster time series data, declare the model and relax the integrality of the\n",
"# binary operation variables\n",
"es.cluster(number_of_typical_periods=20,\n",
" number_of_time_steps_per_period=24 // hpts)\n",
"es.declare_model(use_clustered_data=True)\n",
"es.relax_integrality(include_time_dependent=True, \n",
" include_existence=False, include_modules=False)\n",
"\n",
"# Iteratively run the model (1st stage and 2nd stage)\n",
"for icc in range(3):\n",
"\n",
" # Run optimization with clustered data\n",
" es.optimize(use_clustered_data=True, declare_model=False,\n",
" tee=False, results_file=f'1st_stage_results_{icc}.json')\n",
"\n",
" # Export the configuration to a DataFrame and store it in an Excel-File\n",
" config = es.export_component_configuration(\n",
" write_to_excel=True, file_name=f'1st_stage_config_{icc}.xlsx')\n",
"\n",
" # Import the optimal configuration and fix design variables\n",
" original_model.import_component_configuration(config)\n",
"\n",
" # Run the already declared persistent, full scale model\n",
" original_model.optimize(declare_model=False, tee=False,\n",
" results_file=f'2nd_stage_results_{icc}.json')\n",
"\n",
" # Reset variables to their original (unfixed) values\n",
" original_model.reset_component_variables()\n",
"\n",
" # Add an integer-cut constraint for the gas engines to the aggregated model\n",
" # => Exclude previously found optimal design from the solution space\n",
" es.add_design_integer_cut_constraint(which_instances=['gas_engine'])\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Data and model\n",
"\n",
"Read the input data from a csv-file for one year in hourly resolution (8760 time steps) with pandas."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"# Import the required packages (jupyter magic only required for jupyter notebooks)\n",
"%reload_ext autoreload\n",
"%autoreload 2\n",
"%matplotlib inline\n",
"\n",
"import copy\n",
"import pandas as pd\n",
"import aristopy as ar\n",
"\n",
"# Get data from csv-file \n",
"data = pd.read_csv('testdata.csv', sep=';', decimal='.', index_col=[0])"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
c_elec [EUR/MWh]
\n",
"
Q_demand [MWh]
\n",
"
\n",
" \n",
" \n",
"
\n",
"
01.01.2018 00:00
\n",
"
-5.27
\n",
"
51.828932
\n",
"
\n",
"
\n",
"
01.01.2018 01:00
\n",
"
-29.99
\n",
"
51.123716
\n",
"
\n",
"
\n",
"
01.01.2018 02:00
\n",
"
-56.65
\n",
"
49.202127
\n",
"
\n",
"
\n",
"
01.01.2018 03:00
\n",
"
-63.14
\n",
"
49.060032
\n",
"
\n",
"
\n",
"
01.01.2018 04:00
\n",
"
-64.62
\n",
"
52.669893
\n",
"
\n",
"
\n",
"
01.01.2018 05:00
\n",
"
-67.00
\n",
"
59.055207
\n",
"
\n",
"
\n",
"
01.01.2018 06:00
\n",
"
-72.54
\n",
"
64.771925
\n",
"
\n",
"
\n",
"
01.01.2018 07:00
\n",
"
-76.01
\n",
"
66.424269
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" c_elec [EUR/MWh] Q_demand [MWh]\n",
"01.01.2018 00:00 -5.27 51.828932\n",
"01.01.2018 01:00 -29.99 51.123716\n",
"01.01.2018 02:00 -56.65 49.202127\n",
"01.01.2018 03:00 -63.14 49.060032\n",
"01.01.2018 04:00 -64.62 52.669893\n",
"01.01.2018 05:00 -67.00 59.055207\n",
"01.01.2018 06:00 -72.54 64.771925\n",
"01.01.2018 07:00 -76.01 66.424269"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# print the first 8 rows\n",
"data.head(8)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Aggregate the data for simplification purposes in this notebook by summing/averaging over 4 hours. There is a built-in function in the file utils.py to assist in completing this task."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(array([-38.7625, -70.0425, -64.1025, ..., 61.3 , 64.1375, 43.5775]),\n",
" array([201.21480811, 242.92129394, 306.69530448, ..., 432.8322987 ,\n",
" 463.4687946 , 381.27112154]))"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# USER-INPUT: \"hpts\" (hours per time step)\n",
"hpts = 4\n",
"\n",
"# Aggregate data with \"hours_per_time_step\" for simplification of this notebook\n",
"c_elec = ar.utils.aggregate_time_series(data['c_elec [EUR/MWh]'], hpts, 'mean')\n",
"q_demand = ar.utils.aggregate_time_series(data['Q_demand [MWh]'], hpts, 'sum')\n",
"\n",
"# Show aggregated data on screen\n",
"c_elec, q_demand"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"As in the first example, first, an energy system instance must be created as the main model container. This time an additional logger is added to store information about the modeling and solving process in a log-file. Moreover, a gas source and electricity and heat sinks are added."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Create energy system instance\n",
"es = ar.EnergySystem(\n",
" number_of_time_steps=8760//hpts, hours_per_time_step=hpts,\n",
" interest_rate=0.05, economic_lifetime=20, logging=ar.Logger(\n",
" default_log_handler='file', default_log_level='INFO', \n",
" write_screen_output_to_logfile=True))\n",
"\n",
"# Add a gas source, a heat and electricity sink and a heat bus for collecting\n",
"gas_source = ar.Source(\n",
" ensys=es, name='gas_source', commodity_cost=20, outlet=ar.Flow('Fuel'))\n",
"\n",
"heat_sink = ar.Sink(\n",
" ensys=es, name='heat_sink', inlet=ar.Flow('Heat', link='heat_bus'),\n",
" commodity_rate_fix=ar.Series('q_demand', q_demand))\n",
"\n",
"elec_sink = ar.Sink(\n",
" ensys=es, name='elec_sink', inlet=ar.Flow('Elec'),\n",
" commodity_revenues=ar.Series('c_elec', c_elec))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The bus component collects all heat flows at the inlet and transports them (with or without losses) to the outlet. \n",
"\n",
"
\n",
"\n",
"**Note:** \n",
"\n",
"The same commodity is present at the inlet and the outlet. Therefore, different names for the created variables must be specified in the corresponding flows.\n",
"
"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"heat_bus = ar.Bus(\n",
" ensys=es, name='heat_bus',\n",
" inlet=ar.Flow('Heat', var_name='Heat_inlet'),\n",
" outlet=ar.Flow('Heat', var_name='Heat_outlet'))\n",
"\n",
"gas_boiler = ar.Conversion(\n",
" ensys=es, name='gas_boiler', basic_variable='Heat',\n",
" inlet=ar.Flow('Fuel', 'gas_source'), outlet=ar.Flow('Heat', 'heat_bus'),\n",
" capacity_max=250, capex_per_capacity=60e3,\n",
" user_expressions='Heat == 0.9 * Fuel')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, we add five gas engine instances with identical specifications that are simultaneously created and arranged in a group (see \"instances_in_group\"). The five components can work independently, but have an order for their binary existence and operation variables (see API of the [Conversion class](conversion.rst))."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"# Add maximal 5 gas engine units with fixed size\n",
"gas_engine = ar.Conversion(\n",
" ensys=es, name='gas_engine', basic_variable='Elec',\n",
" inlet=ar.Flow('Fuel', 'gas_source'),\n",
" outlet=[ar.Flow('Heat', 'heat_bus'), ar.Flow('Elec', 'elec_sink')],\n",
" has_existence_binary_var=True, has_operation_binary_var=True,\n",
" additional_vars=ar.Var('Load', ub=1), scalar_params={'dt': hpts},\n",
" capacity=20, min_load_rel=0.2,\n",
" user_expressions=['Load <= BI_OP',\n",
" 'Elec == (20 * Load) * dt',\n",
" 'Heat <= (17 * Load + 2.5 * BI_OP) * dt',\n",
" 'Fuel == (37 * Load + 7 * BI_OP) * dt'],\n",
" capex_per_capacity=1000e3, opex_per_capacity=12.5e3, opex_operation=6,\n",
" instances_in_group=5)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The performance modeling in the user expressions of the gas engines involves using binary operating variables (default name \"BI_OP\" is defined in file utils.py). Due to the y-axis intercept, the resulting efficiency curves show non-linear characteristics.\n",
"\n",
"
\n",
"\n",
"**Note:** \n",
"\n",
"If required, we could also introduce time-varying parameters in the user-expressions to describe the dependency of the conversion rates on, e.g., ambient conditions.\n",
"
"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import matplotlib.pyplot as plt\n",
"import numpy as np\n",
"\n",
"idx = np.linspace(0.2, 1, 20)\n",
"plt.plot(idx, (20*idx)/(37*idx+7)*100, label='el. efficiency')\n",
"plt.plot(idx, (17*idx+2.5)/(37*idx+7)*100, label='th. efficiency')\n",
"plt.xlabel('Load [-]'), plt.ylabel('Efficiency [%]'), plt.legend(loc=4)\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"A heat storage component is also included in the design. The storage tank has a modular design and can consist of up to 5 modules of 50 MWh. It should be noted that, similar to the bus component, different variable names must be specified for the inlet and outlet flows. The small value for the argument \"opex_operation\" is set to prevent simultaneous loading and unloading of the storage. For all other keyword arguments, please consult the documentation of the [storage class API](storage.rst)."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"# Add a heat storage component\n",
"heat_storage = ar.Storage(\n",
" ensys=es, name='heat_storage',\n",
" inlet=ar.Flow('Heat', link='heat_bus', var_name='Heat_charge'),\n",
" outlet=ar.Flow('Heat', link='heat_bus', var_name='Heat_discharge'),\n",
" has_existence_binary_var=True,\n",
" maximal_module_number=5, capacity_per_module=50, # max: 5 * 50 MWh\n",
" capex_per_capacity=12e3, opex_per_capacity=0.015*12e3, # 12,000 €/MWh CAPEX\n",
" self_discharge=0.01, charge_rate=1/6, discharge_rate=1/6,\n",
" opex_operation=1e-9, use_inter_period_formulation=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Finally, we introduce a constraint that the total installed thermal capacity of the gas boiler and all gas engines must be at least 220 MW. To add this condition to the model, we use a back door to the main pyomo model instance, accessible via the method ```add_constraint``` of the energy system. It must be determined whether the constraint is time-dependent (or has a different set), and a function must be assigned as the constraint rule. This function needs to have at least two arguments. The first one represents the energy system instance; the second one is the pyomo model instance. For a time-dependent constraint, we would need to add two more arguments, representing period (p) and time-step (t). \n",
"Furthermore, *aristopy* offers the possibility to add additional variables and objective function contributions to the main pyomo model instance. The functions ```add_variable``` and ```add_objective_function_contribution``` are available for this purpose. Please consult the API of the [EnergySystem class](energySystem.rst) for further information."
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"# Use a backdoor to the pyomo model instance and add a user-defined constraint\n",
"# specifying the installed nominal thermal capacity needs to be at least 220 MW.\n",
"def minimal_installed_capacity_heat(ensys, m):\n",
" expr = 0\n",
" for comp in ensys.components.values():\n",
" if comp.name == 'gas_boiler':\n",
" expr += comp.block.CAP\n",
" elif comp.group_name == 'gas_engine':\n",
" expr += 19.5 * comp.block.BI_EX\n",
" return expr >= 220\n",
"es.add_constraint(name='min_cap_heat', has_time_set=False,\n",
" rule=minimal_installed_capacity_heat)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Before we pass the developed model to a solver, we create a copy by using the method ```deepcopy``` of the copy-package. We will use this original model instance later."
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [],
"source": [
"# Create a deepcopy of the original (full scale) model instance for later use\n",
"original_model = copy.deepcopy(es)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Conventional optimization\n",
"\n",
"As in the first example, the model can simply be by calling the method ```optimize``` of the energy system instance. To control the solving process, we add a time limitation in seconds and a requested relative MIP-gap to the solver Gurobi. "
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"# Run the optimization of the original model instance\n",
"es.optimize(solver='gurobi', time_limit=300, optimization_specs='mipgap=1e-3', \n",
" results_file='full_scale_results.json', tee=False)"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'solver_name': 'gurobi',\n",
" 'time_limit': 300,\n",
" 'optimization_specs': 'mipgap=1e-3',\n",
" 'model_build_time': 4,\n",
" 'model_solve_time': 308,\n",
" 'upper_bound': -163634109.14349702,\n",
" 'lower_bound': -163920996.8127097,\n",
" 'sense': 'maximize',\n",
" 'solver_status': 'aborted',\n",
" 'termination_condition': 'maxTimeLimit'}"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# show basic information about the model building and solving process\n",
"es.run_info"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The optimal design solution of the original (full scale) model instance can be accessed with the ```export_component_configuration``` method and is stated in the table below. The objective function contributions of the individual components are shown in the following figure."
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"ar.Plotter('full_scale_results.json').plot_objective(show_plot=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Time series clustering\n",
"\n",
"To reduce the number of variables and to speed up the calculation, we can perform the optimization for only a few typical periods. To do this, we first call the ```cluster``` method and break the time series data into elements of equal length and combine similar items into clusters. In the example below, 20 typical days are formed from the input data, and the optimization is performed with these aggregated data (\"use_clustered_data\" is True)."
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [],
"source": [
"# Cluster time series data to 20 typical days.\n",
"es.cluster(number_of_typical_periods=20,\n",
" number_of_time_steps_per_period=24//hpts)\n",
"\n",
"# Use clustered data to run the optimization\n",
"es.optimize(use_clustered_data=True, tee=False, results_file='clustered_results.json')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The required time to build and solve the model can be reduced significantly by using clustered input data. The optimal design solution found is comparable to the one of the original model."
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'solver_name': 'gurobi',\n",
" 'time_limit': None,\n",
" 'optimization_specs': '',\n",
" 'model_build_time': 0,\n",
" 'model_solve_time': 2,\n",
" 'upper_bound': -170278358.02702567,\n",
" 'lower_bound': -170278358.02702567,\n",
" 'sense': 'maximize',\n",
" 'solver_status': 'ok',\n",
" 'termination_condition': 'optimal'}"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"es.run_info"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
gas_source
\n",
"
heat_sink
\n",
"
elec_sink
\n",
"
heat_bus
\n",
"
gas_boiler
\n",
"
gas_engine_1
\n",
"
gas_engine_2
\n",
"
gas_engine_3
\n",
"
gas_engine_4
\n",
"
gas_engine_5
\n",
"
heat_storage
\n",
"
\n",
" \n",
" \n",
"
\n",
"
BI_EX
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
1
\n",
"
1
\n",
"
1
\n",
"
-0
\n",
"
-0
\n",
"
1
\n",
"
\n",
"
\n",
"
BI_MODULE_EX
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
{1: 1.0, 2: 1.0, 3: 1.0, 4: -0.0, 5: -0.0}
\n",
"
\n",
"
\n",
"
CAP
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
None
\n",
"
161.5
\n",
"
20
\n",
"
20
\n",
"
20
\n",
"
-0
\n",
"
-0
\n",
"
150
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" gas_source heat_sink elec_sink heat_bus gas_boiler gas_engine_1 \\\n",
"BI_EX None None None None None 1 \n",
"BI_MODULE_EX None None None None None None \n",
"CAP None None None None 161.5 20 \n",
"\n",
" gas_engine_2 gas_engine_3 gas_engine_4 gas_engine_5 \\\n",
"BI_EX 1 1 -0 -0 \n",
"BI_MODULE_EX None None None None \n",
"CAP 20 20 -0 -0 \n",
"\n",
" heat_storage \n",
"BI_EX 1 \n",
"BI_MODULE_EX {1: 1.0, 2: 1.0, 3: 1.0, 4: -0.0, 5: -0.0} \n",
"CAP 150 "
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"es.export_component_configuration()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The following figures show the state of charge (SOC) of the heat storage tank and the associated loading and unloading operations for a selected typical period. Besides, the heat flows, entering and leaving the heat bus component are plotted scaled to a full year cycle."
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAt4AAAG+CAYAAACkiAOYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzdd3jUVd6G8fukFxJaEgEBQVFRLCio9A7Sq6jYlZUuCIiCIh0pokgPiyJiQ6QoivQuTcCCCCJY6UloCaRnzvsHrK/rKgTJ5KQ8n+vKlWQmM3NvQb/8cuYcY61FRERERES8y8d1gIiIiIhIfqDBW0REREQkG2jwFhERERHJBhq8RURERESygQZvEREREZFsoMFbRERERCQb+LkO8LaIiAhbpkwZ1xkiIiIikg/s2LEjzlob+Vf35fnBu0yZMmzfvt11hoiIiIjkA8aYX//uPi01ERERERHJBhq8RURERESygQZvEREREZFsoMFbRERERCQbaPAWEREREckGGrxFRERERLKBBm8RERERkWygwVtEREREJBto8BYRERERyQYavEVEREREsoEGbxERERGRbKDBW0REREQkG2jwFhERERHJBhq8RURERESygQZvL7DW8uOPP7rOEBEREcmX9u3bh8fjcZ3xPzR4Z7GUlBRatGjBLbfcwm+//eY6R0RERCRfWbx4MbfddhvPPvus65T/ocE7iwUGBhIaGkpiYiJPPfWU6xwRERGRfGPatGm0bNmSs2fPEhMTQ0ZGhuuk/6LB2wvGjx9PeHg4ixYt4uOPP3adIyIiIpKneTwe+vXrR7du3fB4PAwaNIi33noLX19f12n/RYO3F5QoUYKRI0cC8NRTT3HmzBnHRSIiIiJ5U1JSEvfffz/jxo3Dz8+PmTNnMnToUIwxrtP+hwZvL+natSuVK1fmwIEDDB482HWOiIiISJ4TFxdHgwYN+PDDDwkPD2fJkiU8/vjjrrP+lgZvL/H19SU6OhofHx8mTJjA119/7TpJREREJM/Yt28fVatWZdOmTZQqVYqNGzfSoEED11kXpMHbiypVqkSPHj3IyMigS5cuOW6Bv4iIiEhutHHjRqpWrcr+/fu57bbb2LJlCzfddJPrrIvS4O1lw4cPp0SJEmzdupUZM2a4zhERERHJ1T788EPq16/P8ePHadq0KevXr6dEiRKuszJFg7eXhYeHM3HiRAD69+/P0aNHHReJiIiI5D7WWl5++WXuvfdeUlJS6NKlCx9//DEFChRwnZZpGryzQdu2bWnatCmnT5+mT58+rnNEREREcpX09HS6dev2+6E4Y8eOZerUqfj5+TkuuzQavLOBMYbJkycTHBzM+++/z4oVK1wniYiIiOQKZ86coVWrVkRHRxMYGMjcuXPp169fjtwu8GI0eGeTsmXLMmjQIAC6detGcnKy4yIRERGRnO3w4cPUqlWLzz77jKJFi7Jq1Srat2/vOusf0+Cdjfr27UuFChXYv38/o0aNcp0jIiIikmPt2rWLKlWq8NVXX1GuXDk2b95M9erVXWddFg3e2cjf35/o6GgARo8ezd69ex0XiYiIiOQ8K1eupHr16hw4cICqVauyefNmrr32WtdZl02DdzarUaMGHTt2JDU1la5du2KtdZ0kIiIikmPMmjWLJk2aEB8fT/v27Vm1ahURERGus7KEBm8HxowZQ0REBGvWrOGdd95xnSMiIiLinLWWwYMH8/jjj5Oenk6/fv2YM2cOwcHBrtOyjAZvB4oWLcq4ceOAc+u+T5w44bhIRERExJ3U1FQeffRRhg0bho+PD1OnTmXs2LH4+OStUTVv/afJRR555BFq165NbGws/fv3d50jIiIi4sTJkye5++67efvttwkNDWXRokV07drVdZZXaPB2xBhDdHQ0/v7+zJgxg40bN7pOEhEREclWv/zyC9WrV2ft2rUUK1aM9evX06xZM9dZXqPB26Hy5cvz3HPPAdClSxfS0tIcF4mIiIhkj+3bt1OlShX27NlDhQoV2Lp1K7fffrvrLK/S4O3Y888/zzXXXMOuXbsYP3686xwRERERr1u0aBG1a9fm2LFj1K9fn40bN1K6dGnXWV6nwdux4OBgpk6dCsCQIUP45Zdf3AaJiIiIeNHkyZNp06YNiYmJPPbYY3z22WcULFjQdVa20OCdAzRq1Ij777+fpKQkevToob29RUREJM/xeDz06dOHp556Co/Hw7Bhw5g5cyYBAQGu07KNBu8c4tVXXyU8PJzFixezcOFC1zkiIiIiWSYxMZH27dszfvx4/P39mT17Ni+++CLGGNdp2UqDdw5RvHhxRo0aBUDPnj1JSEhwXCQiIiJy+WJiYqhXrx4LFiygYMGCLF26lIcffth1lhMavHOQzp07c+edd3Lo0CEGDRrkOkdERETksuzdu5eqVauydetWrrrqKjZt2kS9evVcZzmjwTsH8fX1JTo6Gh8fHyZOnMiXX37pOklERETkH9mwYQNVq1blp59+olKlSmzZsoUbb7zRdZZTGrxzmNtuu41evXrh8Xjo0qULGRkZrpNERERELsmcOXNo0KABJ0+epEWLFqxbt45ixYq5znJOg3cONHToUEqWLMm2bduIjo52nSMiIiKSKdZaRo8eTYcOHUhNTaVHjx4sXLiQ0NBQ12k5ggbvHCgsLIyJEycC5w7YOXLkiOMiERERkQtLS0ujc+fODBgwAGMMr776KhMnTsTX19d1Wo6hwTuHat26Nc2bNyc+Pp7evXu7zhERERH5WwkJCbRo0YIZM2YQFBTEvHnz6N27d77bLvBiNHjnUMYYJk+eTEhICB988AHLli1znSQiIiLyPw4ePEjNmjVZtmwZERERrFmzhrZt27rOypE0eOdgV111FUOGDAGgW7duJCUluQ0SERER+YOdO3dSpUoVvvnmG6677jq2bNlClSpVXGflWBq8c7inn36am2++mZ9++omRI0e6zhEREREBYNmyZdSoUYNDhw5Ro0YNNm3axDXXXOM6K0fT4J3D+fv7/76zydixY9mzZ4/jIhEREcnvXn/9dZo1a0ZCQgL3338/K1asoGjRoq6zcjwN3rlAtWrVePLJJ0lLS6Nr165Ya10niYiISD5kreWFF17gySefJCMjgwEDBvDuu+8SFBTkOi1X0OCdS4wePZrIyEjWrVvH7NmzXeeIiIhIPpOSksJDDz3ESy+9hK+vL9OnT+ell17Cx0fjZGbpv6lcokiRIrzyyisAPPPMMxw/ftxxkYiIiOQXJ06coGHDhrz33nsUKFCATz/9lE6dOrnOynU0eOciDz30EPXq1SMuLo5nn33WdY6IiIjkAz/99BPVqlVjw4YNlChRgg0bNtC4cWPXWbmSBu9cxBjD1KlTCQgIYObMmWzYsMF1koiIiORhW7dupUqVKuzdu5dbbrmFrVu3UrFiRddZuZYG71zm+uuvp3///gB06dKF1NRUx0UiIiKSFy1cuJC6desSGxtLo0aN2LBhAyVLlnSdlatp8M6FBgwYQLly5di9ezevvvqq6xwRERHJY1577TXatWtHUlISHTt25NNPPyU8PNx1Vq6nwTsXCgoKYtq0aQAMGzaMn3/+2XGRiIiI5AUZGRn06tWL3r17Y61l5MiRzJgxA39/f9dpeYIG71yqQYMGPPDAAyQlJdG9e3ft7S0iIiKX5ezZs7Rt25aJEycSEBDAu+++y/PPP48xxnVanqHBOxd79dVXKViwIEuWLGH+/Pmuc0RERCSXOnbsGHXr1mXRokUULlyY5cuX88ADD7jOynM0eOdiV1xxBaNHjwagV69exMfHOy4SERGR3GbPnj1UqVKFbdu2UbZsWTZt2kTt2rVdZ+VJGrxzuU6dOlGlShUOHz7Miy++6DpHREREcpG1a9dSrVo1fvnlF+688042b95M+fLlXWflWRq8czkfHx+io6Px9fVl8uTJ7Nixw3WSiIiI5ALvvPMOjRo14tSpU7Ru3Zo1a9ZwxRVXuM7K0zR45wG33norTz/9NB6Ph86dO5ORkeE6SURERHIoay3Dhw/n4YcfJi0tjV69ejFv3jxCQkJcp+V5GrzziCFDhlCqVCl27NjB1KlTXeeIiIhIDpSWlsa//vUvBg0ahDGG1157jddeew1fX1/XafmC08HbGDPTGBNjjNn1h9uKGGNWGGP2nf9c+A/3DTDG7DfG7DXG3O2mOmcqUKAAkyZNAuCFF17g0KFDjotEREQkJzl9+jRNmzZl5syZBAcHs3DhQnr16uU6K19xfcV7FtD4T7f1B1ZZa68FVp3/HmPMjcD9QIXzj5lqjNFfz/6gVatWtGrVioSEBHr37u06R0RERHKIAwcOULNmTVauXElUVBRr166lVatWrrPyHaeDt7V2PXDiTze3At46//VbQOs/3D7HWptirf0Z2A/cmS2hucjEiRMJDQ3lww8/ZMmSJa5zRERExLGvvvqKu+66i2+//Zby5cuzZcsW7rxTI5QLrq94/5UrrLVHAM5/jjp/+5XAgT/83MHzt8kflC5dmqFDhwLQvXt3EhMTHReJiIiIK5999hk1a9bkyJEj1K5dm02bNlG2bFnXWflWThy8/85fnVf6l+ekG2M6GWO2G2O2x8bGejkr5+nZsye33HILP//8MyNGjHCdIyIiIg5Mnz6dli1bcvbsWR588EGWLVtG4cKFL/5A8ZqcOHgfM8YUBzj/Oeb87QeBUn/4uZLA4b96Amvtv621la21lSMjI70amxP5+/szffp0jDG8/PLLfPfdd66TREREJJt4PB6ee+45unTpQkZGBgMHDuTtt98mMDDQdVq+lxMH70XAo+e/fhT4+A+332+MCTTGlAWuBb5w0JcrVKlShU6dOpGenk6XLl3weDyuk0RERMTLkpOT6dChA2PHjsXPz4833niD4cOHY8xfLRyQ7OZ6O8H3gc3A9caYg8aYjsBooKExZh/Q8Pz3WGu/A+YCu4GlQHdrrU6KuYBRo0YRFRXF559/zqxZs1zniIiIiBfFxcXRoEED5s6dS1hYGIsXL+aJJ55wnSV/YKz9y2XSeUblypXt9u3bXWc489577/Hggw9SpEgR9u7dS0REhOskERERyWL79++nadOm7Nu3j5IlS/LZZ59x8803u87Kl4wxO6y1lf/qvpy41ESyUIcOHWjQoAEnTpygX79+rnNEREQki23evJmqVauyb98+KlasyJYtWzR051AavPM4YwxTp04lMDCQWbNmsW7dOtdJIiIikkXmzZtH3bp1iYuLo0mTJqxfv54rr9RuyzmVBu984Nprr2XAgAEAdOnShdTUVMdFIiIicjmstYwbN4727duTkpJC586dWbRoEWFhYa7T5AI0eOcT/fv357rrruP777/n5Zdfdp0jIiIi/1B6ejo9evT4fQnpmDFjmDZtGn5+fo7L5GI0eOcTgYGBTJs2DYARI0bw448/Oi4SERGRS3XmzBlat27N1KlTCQgIYM6cOTz77LPaLjCX0OCdj9SrV4+HHnqI5ORkunfvTl7f0UZERCQv+c+x74sXL6ZIkSKsWrWK++67z3WWXAIN3vnMK6+8QqFChVi2bBlz5851nSMiIiKZ8N1331GlShW+/PJLrrnmGjZv3kyNGjVcZ8kl0uCdz0RFRTFmzBgAnn76aU6fPu24SERERC5k1apVVKtWjd9++40qVaqwefNmrrvuOtdZ8g9o8M6H/vWvf1G1alWOHj3KCy+84DpHRERE/sZbb71F48aNiY+Pp127dqxevZrIyEjXWfIPafDOh3x8fIiOjsbX15epU6eybds210kiIiLyB9ZahgwZwmOPPUZ6ejp9+/Zl7ty5BAcHu06Ty6DBO5+65ZZb6NOnD9ZaOnfuTHp6uuskERERAVJTU3nssccYOnQoPj4+TJ48mXHjxuHjo7Ett9P/gvnY4MGDKV26NF999RVTpkxxnSMiIpLvnTp1isaNGzN79mxCQkL4+OOP6d69u+ssySIavPOx0NBQJk+eDMDAgQM5ePCg4yIREZH869dff6V69eqsWbOGYsWKsX79epo3b+46S7KQBu98rkWLFrRp04YzZ87Qq1cv1zkiIiL50vbt26lSpQq7d+/mxhtvZMuWLVSqVMl1lmQxDd7ChAkTKFCgAAsWLODTTz91nSMiIpKvfPLJJ9SuXZujR49St25dNm7cyFVXXeU6S7xAg7dQqlQphg0bBkCPHj04e/as4yIREZH8YcqUKbRu3ZrExEQeeeQRli5dSqFChVxniZdo8BYAnnrqKSpWrMivv/7K8OHDXeeIiIjkaR6Ph2eeeYYePXrg8XgYMmQIs2bNIiAgwHWaeJEGbwHAz8+P6dOnY4zhlVdeYdeuXa6TRERE8qSkpCTuvfdeXnnlFfz8/Jg1axaDBw/GGOM6TbxMg7f87s4776Rr166kp6fTuXNnPB6P6yQREZE8JTY2lnr16jF//nwKFizI0qVLefTRR11nSTbR4C3/ZeTIkRQrVoxNmzYxc+ZM1zkiIiJ5xg8//ECVKlXYsmULpUuXZuPGjdSvX991lmQjDd7yXwoVKsT48eMBePbZZ4mJiXFcJCIikvt9/vnnVK1alZ9++onbb7+dLVu2UKFCBddZks00eMv/uO+++2jYsCEnT56kX79+rnNERERytQ8++ID69etz4sQJmjVrxrp16yhevLjrLHFAg7f8D2MMU6dOJTAwkNmzZ7NmzRrXSSIiIrmOtZYxY8Zw//33k5qaSrdu3fjoo48oUKCA6zRxRIO3/KVy5coxcOBAALp06UJKSorjIhERkdwjPT2drl270r9/fwDGjRvH5MmT8fPzc1wmLmnwlr/Vr18/rr/+en744QfGjh3rOkdERCRXSEhIoEWLFkyfPp3AwEA+/PBD+vbtq+0CRYO3/L3AwECio6OBc7ud7Nu3z3GRiIhIznbo0CFq1arF0qVLiYiIYPXq1dxzzz2usySH0OAtF1SnTh0eeeQRUlJS6N69O9Za10kiIiI50s6dO6lSpQpff/011157LZs3b6ZatWqusyQH0eAtFzVu3DgKFy7MihUrmDNnjuscERGRHGf58uXUqFGDgwcPUr16dTZt2kS5cuVcZ0kOo8FbLioyMvL3Nd69e/fm1KlTjotERERyjpkzZ9KsWTMSEhK49957WblyJREREa6zJAfS4C2Z8sQTT1C9enWOHTvG888/7zpHRETEOWstgwcPpmPHjqSnp/Pss8/y/vvvExQU5DpNciiT19fsVq5c2W7fvt11Rp6wa9cubrvtNjIyMti8eTN33XWX6yQRkXylTP/FrhPkD+K3L+Lkqn/j4+PDlClT6NKli+skyQGMMTustZX/6j5d8ZZMu+mmm+jbty/WWjp37kx6errrJBEREScSf9zGydWvAzB79mwN3ZIpGrzlkgwaNIgyZcrwzTffMHHiRNc5IiIi2S415mfiFo0F66Fg9Q48+OCDrpMkl9DgLZckJCSEyZMnA+eG8AMHDjguEhERyT4ZZ04SM28YNjWJkBtqU7D6A66TJBfR4C2XrFmzZrRr146zZ8/Ss2dP1zkiIiLZwpOWQsyC4WQkxBJYojwRTXvpNEq5JBq85R+ZMGECBQoU4KOPPmLRokWuc0RERLzKWg/HF48n9cgP+Ba8gsi2AzF+Aa6zJJfR4C3/yJVXXsmIESMA6NGjB2fOnHFcJCIi4j2nN7xL4t7PMQHBRLV7Ed/QQq6TJBfS4C3/WPfu3bn99ts5cOAAQ4cOdZ0jIiLiFWd2reb05g/A+BDZqj8BkWVcJ0kupcFb/jE/Pz+io6MxxjB+/Hh27tzpOklERCRLJR/8juNLz+3iVaRBJ4KvruS4SHIzDd5yWe644w66d+9ORkYGnTt3xuPxuE4SERHJEmknjxC7YCRkpBNWqQVhtzd3nSS5nAZvuWwjRoygePHibNmyhRkzZrjOERERuWye5DPEzBuKJymeoKsrUbjev1wnSR6gwVsuW8GCBXnttdcA6N+/P8eOHXNcJCIi8s/ZjHRiPxpF+omD+EdcRWTL5zA+vq6zJA/Q4C1Zon379jRu3JhTp07Rt29f1zkiIiL/iLWWEyuiSf71G3xCChF1z2B8AkNcZ0ke4ec6QPIGYwxTpkyhQoUKvPvuuzz++OPUr1/fdZbkYB6PhzfeeIORI0eSkpJCVFQUkZGRv3/+49d/vK1QoUI6sEJEvCZh+8ec+WYp+PoT1XYgfgWjXCdJHqLBW7LM1VdfzYsvvsgLL7xA165d2blzJ0FBQa6zJAf6+uuv6dq1K1u2bPn9tqNHj2bqsf7+/kRERPztYP7n2woWLKhBXUQyJXH/Vk6ufgOAiGa9CbyyvOMiyWs0eEuWeuaZZ3jnnXfYs2cPY8aMYfDgwa6TJAeJj49n0KBBTJo0CY/HQ/HixRk/fjzVq1cnNjaW2NhYYmJi/uvzn29LSEjgyJEjHDlyJFOv6e/v//tQnplhPTw8XIO6SD6Ueuwn4ha9DFgK1niQ0BtquU6SPEiDt2SpgIAAoqOjqV27Ni+99BIdOnTguuuuc50ljllrmTt3Lr179+bIkSP4+Pjw9NNPM3ToUMLDwwEoWbJkpp4rOTmZuLi4vx3M/3xbQkIChw8f5vDhw5l6/v8M6hca0P/4tQZ1kdwv/cwJYuYPw6YlE1qhLgWr3e86SfIoDd6S5WrVqsXjjz/Om2++SdeuXVm5cqUGk3xs37599OjRg+XLlwNQpUoVpk2bRsWKFf/R8wUFBVGyZMlLGtQzezU9NjaWM2fOXNKgHhAQkKmr6f/5OiwsTH8eRHIQT1oysfOHk5EQR+CVN1K0cU/9GRWv0eAtXjF27FgWLVrE6tWree+993jwwQddJ0k2S05OZtSoUYwePZrU1FQKFy7MmDFj6NixIz4+2behUlBQEKVKlaJUqVKZ+vmkpKT/GcovNLSfPXuWQ4cOcejQoUw9f0BAwEXfQPrHrzWoi3iPtR6Of/oqqUf34VfwCiLbvoDx83edJXmYBm/xioiICF5++WWeeOIJ+vTpQ9OmTSlcuLDrLMkmS5cupUePHvz4448APP7444wZM4bIyEjHZRcXHBxM6dKlKV26dKZ+/j+DemaupsfGxnL27FkOHjzIwYMHM/X8gYGBmdrt5T9fFyhQQIO6SCadWv82iT9swgSEEHnPYHxDCrpOkjxOg7d4zWOPPcabb77Jhg0bGDBgANHR0a6TxMsOHjxI7969mTdvHgA33XQT06ZNo0aNGo7LvOdSB/XExMRMX02PjY0lMTHxkgb1oKCgC65J//NtoaGhGtQlXzrz7Urit3wIxofI1gMIiMjcn2GRy2Gsta4bvKpy5cp2+/btrjPyrd27d3PrrbeSnp7Opk2bqFq1qusk8YL09HQmTpzI4MGDOXPmDKGhoQwZMoRevXrh769f216Os2fPXtLSl6SkpEt6/qCgoL8c0MuUKcOjjz5KWFiYl/6TyT9Rpv9i1wl5QvJv33LsgxfBk06RRt0Iu63pZT3fL6ObZVGZ5AXGmB3W2sp/dZ+ueItX3XjjjfTr149Ro0bRuXNnduzYoUEsj9m4cSPdunVj586dALRt25bXXnst02uq5cJCQ0MJDQ2lTJkymfr5s2fPZvpqekxMDMnJyfz222/89ttv//Ncb7/9NkuXLtUyMclT0k4cInbhS+BJJ6xyq8seukUuhQZv8bqBAwcyZ84cvv32WyZMmMAzzzzjOkmyQFxcHP379+eNN84dNlG2bFkmT55M06b6l5hLoaGhlC1blrJly170Z621v19R//NgHh0dzRdffEG9evVYvnx5rlifL3IxGUkJxMwfhic5geBr7qBw3SdcJ0k+o6Umki2WLl1KkyZNCAkJYffu3Vx11VWuk+Qf8ng8vPnmmzz33HMcP36cgIAAnnvuOQYMGEBwcLDrPMkiBw4coH79+uzbt48bb7yRlStXUrx4cddZ+Z6WmvxzNiONY3MHk/LbTvwjy1DswbH4BIZkyXNrqYn80YWWmmTfnl6SrzVu3Jj27duTmJhIz549XefIP7Rz505q1qzJv/71L44fP079+vXZuXMnw4YN09Cdx5QqVYr169dToUIFdu/eTe3atTlw4IDrLJF/xFrLieXTSPltJz6hhYi6Z1CWDd0il0KDt2Sb1157jbCwMBYtWsRHH33kOkcuQUJCAn379uX2229n06ZNFCtWjPfff58VK1Zw/fXXu84TLylWrBhr166lYsWK7Nu3j1q1avHTTz+5zhK5ZPFfLOTMzuUYvwCi2r6IX3iU6yTJp7TURLLVpEmT6NmzJyVLlmT37t3aMSGHs9Yyf/58nn76aQ4dOoSPjw89evRg2LBhFCyo/W7zi5MnT9KkSRO2bt3KlVdeyapVq/QXLke01OTSJf6w+dybKbFEtOpPaPm8u72p/DdXS4C01ERyjG7dulGpUiUOHjzIkCFDXOfIBezfv58mTZrQvn17Dh06xJ133sm2bduYMGGChu58pnDhwixfvpyaNWty6NAhateuza5du1xniVxUytH9xH06DrAUqvWIhm5xToO3ZCtfX1+mT5+Oj48PEyZM4Ouvv3adJH+SnJzM0KFDuemmm1i2bBmFChUiOjqaTZs2cfvtt7vOE0fCw8NZsmQJDRo04NixY9SpU4cvv/zSdZbI30pPiCN2/jBsWgqhN9UnvEp710kiGrwl+1WqVIkePXqQkZFB586dycjIcJ0k561YsYKbb76ZIUOGkJKSwqOPPsrevXvp3Lkzvr6+rvPEsdDQUD755BOaNWvG8ePHqVevHlu2bHGdJfI/PKnJxM4fTsaZEwSWrEDRu3vohFbJETR4ixPDhw+nRIkSfPHFF/z73/92nZPvHT58mPvuu49GjRqxf/9+brzxRtatW8esWbOIitKbkOT/BQUFsWDBAtq1a8fp06dp2LAh69atc50l8jtrPcR9Oo7UYz/iV6g4kW2ex/jp4DbJGTR4ixPh4eFMmDABgAEDBnD06FHHRflTeno6r732GuXLl2fu3LmEhIQwZswYvv76a2rVquU6T3KogIAA5syZwwMPPMCZM2do0qQJy5cvd50lAsCptbNI2rcFn8BQou4ZhG+I3pMiOYcGb3GmXbt2NG3alNOnT9OnTx/XOfnO5tOmgqYAACAASURBVM2bqVy5Mr179yYhIYHWrVuzZ88enn32Wfz9dXVILszPz4/Zs2fTsWNHkpKSaNGiBZ988onrLMnnEr5ZTvwXC8DHl4jWA/AvWsp1ksh/0eAtzhhjmDx5MsHBwbz//vu6YpZNjh8/TqdOnahWrRrffPMNZcqU4ZNPPmHhwoWULl3adZ7kIr6+vvz73/+mR48epKam0rZtWz788EPXWZJPJf+6kxPLpwBQpGFXgstUdFwk8r80eItTZcuWZdCgQcC5rQaTkpIcF+Vd/znqvXz58syYMQN/f3+ef/55vvvuO5o3b+46T3IpHx8fJk6cSL9+/UhPT+f+++/nnXfecZ0l+UzaiUPEfvQSeDIIv6MNYRUbu04S+UsavMW5vn37UqFCBX788UdGjRrlOidP+vbbb6lduzZPPPEEcXFx1K1bl507dzJy5EhCQnRsslweYwxjxoxh0KBBeDweHnnkEWbMmOE6S/KJjKR4YuYNwZN8huByd1GozmOuk0T+lgZvcc7f35/o6GgARo8ezffff++4KO84c+YM/fr147bbbuPzzz/niiuu4J133mHVqlWUL1/edZ7kIcYYhg4dyqhRo7DW0qlTJyZNmuQ6S/I4m5FG7MKXSD95BP+oq4lo8QzGR1ufSs7l93d3GGMy8263s9ba6VnYI/lUjRo16NixI2+88QZdu3Zl9erV2nP1MlhrWbhwIb169eLgwYMYY+jevTsjRoygUKFCrvMkD+vfvz8hISH06tWLnj17kpSUxLPPPus6S/Igay3Hl04h5cAufAsUIardIHwCgl1niVzQha549wMKAGEX+OjrrTBjzC/GmG+NMV8bY7afv62IMWaFMWbf+c+FvfX6kv3GjBlDREQEa9eu5e2333adk2v99NNPNG/enHbt2nHw4EEqV67Mtm3bmDx5soZuyRY9e/Zk+vTpGGN47rnnGDJkCNZa11mSx8RvncfZXSsxfoFEthuEX3iE6ySRi/rbK97A29baYRd6sDEmNIt7/qyutTbuD9/3B1ZZa0cbY/qf//45LzdINilatCjjxo3jscceo2/fvjRv3pwiRYq4zso1UlJSGDt2LC+99BLJyckULFiQUaNG0alTJ506KdmuU6dOBAUF8fjjjzN06FCSkpIYPXq0fpMlWSJx7yZOrXsLgIjmfQksVs5xkUjm/O0Vb2vtRX83mJmfyWKtgLfOf/0W0DqbX1+87JFHHqF27drExcXx3HP6O1VmrVy5kltuuYVBgwaRnJzMww8/zN69e+natauGbnHmkUceYc6cOfj5+TF27Fh69eqFx+NxnSW5XMqRfcR9+goAhWo/Rsj11RwXiWTeRd9caYwJNMY8YIx53hgz6D8f2dBmgeXGmB3GmE7nb7vCWnsE4PxnnWWdxxhjiI6Oxt/fn9dff52NGze6TsrRjhw5QocOHWjYsCE//PADN9xwA2vWrGH27NlcccUVrvNEaN++PfPnzycgIIBJkybRuXNnMjIyXGdJLpUeH0fsguHY9BRCb25I+F3tXCeJXJLM7GryMeeuNKcDZ//w4W3VrbW3A02A7saYTJ9fbYzpZIzZbozZHhsb671C8Yry5cv/frW7S5cupKWlOS7KedLT05k0aRLly5dnzpw5BAcHM2rUKL7++mvq1KnjOk/kv7Rs2ZJFixYRFBTE66+/zmOPPUZ6errrLMllPKlJxMwfRsaZEwSWuomid3fT0iXJdTIzeJe01t5nrR1rrX3lPx/eDrPWHj7/OQZYCNwJHDPGFAc4/znmbx77b2ttZWtt5cjISG+nihc8//zzXHPNNezatYtXX33VdU6O8sUXX3DnnXfSs2dP4uPjadGiBbt376Z///4EBAS4zhP5S3fffTdLliwhNDSUd955hw4dOpCamuo6S3IJ68kg7pOXSYv5Cb/CJYhs8zzG1991lsgly8zgvckYc7PXS/7AGBNqjAn7z9dAI2AXsAh49PyPPcq5q/GSBwUHBzNlyrmjf4cOHcovv/ziNigHOHnyJF26dKFKlSp89dVXlC5dmo8//phFixZRpkwZ13kiF1WnTh1WrFhBeHg48+bN45577iE5Odl1luQCJ9e+SdL+L/AJKkDUPYPxDQ53nSTyj/zt4H1+K7+dQA3gS2PMXmPMzj/c7k1XAJ8bY74BvgAWW2uXAqOBhsaYfUDD899LHnX33Xdz3333kZSURI8ePfLtdmTWWt566y2uv/56pk+fjq+vL/3792f37t20bNnSdZ7IJalatSqrV6+mSJEifPLJJ7Rs2ZLExETXWZKDJXy9lIRtH4GPL5Gtn8e/yJWuk0T+MfN3w4wx5qoLPdBa+6tXirJY5cqV7fbt211nyD905MgRypcvT3x8PPPnz6dt27auk7LVd999R9euXdmwYQMAtWvXZurUqdx4442Oy0Quz86dO2nYsCExMTHUqlWLTz/9lLCwMNdZOV6Z/otdJ2SrpF++JmbuILAeijbpSYFbGrlOklzkl9HNnLyuMWaHtbbyX913oaUmvTm3rjrdWvvrnz+8UiryJ8WLF2fUqFHAuUM5EhISHBdlj7Nnz/Lcc89RsWJFNmzYQFRUFLNnz2bNmjUauiVPuOWWW1i3bh0lSpRg/fr1NGrUiFOnTrnOkhwk7fgBYj8aBdZD+F3tNHRLnnChwXs/0AbYeP4UyfeMMd2NMbcZYzKzNlwkS3Tu3Jk77riDQ4cO8eKLL7rO8SprLR999BE33HADY8eOJSMjg65du/L999/z8MMP6x38kqeUL1+e9evXc9VVV7Flyxbq1atHXFzcxR8oeV5G4mli5g3Fppwl+NoqFKr96MUfJJILXOgAncnW2gestWWAqsAC4BrgQ0CXJSTb+Pr6Mn36dHx8fJg0aRJffvml6ySv+Pnnn2nZsiVt2rThwIEDVKpUia1btzJ16lQKFy7sOk/EK6655hrWr19PuXLl+Oqrr6hbty5Hjx51nSUO2fQ0YheOJP3UUQKuuIaI5s+g632SV1zw/8nmnFs4t493K6A2566Ee307QZE/uu222+jZsycejyfPHcCRmprKSy+9RIUKFfj0008JDw9n0qRJbN26lTvuuMN1nojXlS5dmnXr1nHDDTewa9cuateuzcGDB11niQPWWo4vnUjKwd34FihCZLsX8QkIcp0lkmUutKvJCmAL0B1IAV6y1lay1ja21g7NrkCR/xg2bBhXXnkl27dvJzo62nVOlli9ejW33norL7zwAklJSTzwwAPs3buXHj166Kh3yVdKlCjB2rVrufXWW/nhhx+oVauWthHNh+I3z+Xsd2sw/oFEthuEX1iE6ySRLHWhK94/ce7Y9mvPf5QzxuhPgDgTFhbGxIkTgXMH7Bw+fNhx0T939OhRHnroIerXr8/333/P9ddfz6pVq3j33XcpVqyY6zwRJ6Kioli9ejV33HEHP//8MzVr1mTfvn2usySbnN2zgVMb3gYMES36EVisnOskkSx3oTXena21VYDWwFqgEvCOMWaHMeatbOoT+S9t2rShefPmxMfH06dPH9c5lywjI4MpU6ZQvnx53n33XYKCghgxYgTffPMN9erVc50n4lyRIkVYuXIl1atX5+DBg9SqVYvvvvvOdZZ4WcrhvRz/bDwAheo8Tsi1VRwXiXhHZt6tkAIkAknnvy4J3O7NKJG/Y4xh8uTJhISE8MEHH7Bs2TLXSZm2bds27rrrLnr06MHp06dp1qwZu3fv5oUXXiAwMNB1nkiOER4eztKlS6lXrx5Hjx6lTp06fPXVV66zxEvS42OIWTAcm55KgVsaEX5nG9dJIl5zoTXe440xW4EjwDAgDJgOXG+tzdYj5EX+6KqrrmLw4MEAdOvWjaSkJMdFF3bq1Cm6devGXXfdxY4dOyhVqhQLFy7kk08+oWzZsq7zRHKkAgUK8Omnn9KkSRPi4uKoV68eW7dudZ0lWcyTkkjMvGF4zp4i6KpbKNKom7ZNlTztQidX9gQ2AV9Za3PtFhI6uTJvSktL4/bbb2fXrl288MILjBgxIlOPy85T36y1nN29lpOr38CTeAp8fAm/ozUFq92PT0BwtnXkBq5OF5OcLyUlhQ4dOrBw4ULCwsJYvHgxNWvWdJ3lTF46udJ6MohdMIKkH7fhV+RKij38Cr5BBVxnSR6S206u/BzwALcaY27/84dXSkUyyd/fn+nTpwMwduxY9uzZ47jov6XFHeDYnOc5/ukreBJPEViyAsUfm0DhOo9r6Ba5BIGBgXzwwQfcf//9JCQk0LhxY1auXOk6S7LAydVvkPTjNnyCwoi6Z7CGbskX/C5w33bgOyD2/Pd//N2PBfROMHGqWrVqPPnkk8yYMYMuXbqwdu1a57+i9KQlc3rTHOK/+Ag86fgEh1O4bkdCb6rnvE0kt/L39+edd94hKCiIWbNm0bx5c+bPn0+zZvpNSW6V8NVnJOxYBD5+RLZ9Af/CJVwniWSLC13x7guc5tybKt8EWlhr657/0NAtOcLo0aOJjIxk/fr1vPWW2812Evdt5fDr3YjfMg88GRSo2JgST06nwM31NXSLXCZfX1/eeOMNunbtSkpKCm3atGHBggWus+QfSPr5S06sOHcWQ9HGTxFU6ibHRSLZ50LbCY631tYAegClgFXGmLnGmIrZVidyEUWKFOGVV84dpPrMM88QFxeX7Q3pp2OImT+c2AXDyYiPwT/qaoo99DJF7+6Bb3BYtveI5FU+Pj5MmTKFPn36kJaWxr333st7773nOksuQWrcb8R+NBqsh/Aq7Slwc33XSSLZ6qLbCVprfwY+BpYDdwLXeTtK5FI89NBD1K1bl+PHj/Pcc89l2+vajDROb/mQw693JWn/VkxAMIXrd6L4o+MJvLJ8tnWI5CfGGMaNG8fAgQPJyMjgoYceYubMma6zJBMyEk8TO28oNjWRkOuqUajWw66TRLLdhbYTvNoY8/z5LQWHAt8A5a21c7OtTiQTjDFMmzaNgIAAZs6cyYYNG7z+msm/7eTImz05te4tbHoKITfUosST0wmv3BLjo6PeRbzJGMPw4cMZMWIE1lo6duzIlClTXGfJBdj0VGIXjCD99DECil1L0eZ9MCYzR4mI5C0XenPlfmAn5652xwOlgW7/WatqrX3V63UimXT99dfTv39/hg0bRpcuXfjqq68ICAjI8tfJOHuSk2tmcva7NQD4FS5BkYZdCS57W5a/lohc2AsvvEBISAh9+vShR48eJCUl8cwzz7jOkj+x1hK3ZAIph/bgGxZBZLsX8fEPcp0l4sSF/ro5DFjIuS0FC3DuAJ0/fojkKAMGDKBcuXLs3r3793XfWcV6Mkj46jMOz+hybuj29adgjQcp8cQUDd0iDvXu3Ztp06YB0K9fP4YPH87fnU8hbpzeNIfE3esw/kFE3TMIvwJFXCeJOPO3V7yttUOysUPksgUFBTFt2jQaNmzIsGHDuO+++7j66qsv+3lTju7nxPIppB7Zd+51rq5EkQZd8C9c/LKfW0QuX5cuXQgKCqJjx44MGjSIpKQkRo4cqd2EcoCzu9dx+vN3AUNEy34ERF3+P5NFcrMLrfHudLEHZ+ZnRLJTgwYNeOCBB0hOTqZ79+6XdeXLk3yGEyuiOTq7D6lH9uFboCgRrQcQdc8QDd0iOcxjjz3Gu+++i6+vL6NGjaJ379668u1YyqHvifvsNQAK1+tISLm7HBeJuHehNd79jTEX2pvNAL2Af2dtksjleeWVV1i8eDFLly5l3rx5tG/f/pIeb60lcc86Tqx+Hc/ZU2B8CL+jDQWrd8AnMMRL1SJyue6//36CgoK49957mTBhAsnJyUydOhUfH72JL7ulnz5GzIIRkJFGgYqNCavcynWSSI5wocF7HdDiIo9fkYUtIlmiWLFijB49mq5du9KrVy/uvvtuwsPDM/XYtOMHOLFiGsm/7gQg8MobKXJ3NwIiy3ixWESySuvWrVm0aBFt2rRh+vTpJCUl8cYbb+Dnd6F/3UlW8qQkEjNvGJ7EUwRdVZEiDbpo2Y/IeRda4/14doaIZKVOnToxa9Ystm7dysCBA5k4ceIFf96TlszpzXOJ37rg/496r/M4oTfX15ZXIrlM48aNWbx4MS1atGD27NkkJyfzzjvv4O/v7zotz7OeDGI/HkNa3K/4FSlJZOv+GF/9pUfkPzRRSJ7k4+PD9OnT8fX1ZcqUKWzfvv1vfzbxx20cfqM78ZvngiedArfeTYknoylwS0MN3SK5VL169Vi+fDlhYWHMnTuXe+65h5SUFNdZed7JVTNI/nkHPsHhRN0zGJ+gAq6TRHIUTRWSZ91666306tULj8dDly5dyMjI+K/70+NjiFkwgth5Q8k4fQz/qLLnjnpv/BS+wZlbmiIiOVf16tVZtWoVhQsXZtGiRbRq1YrExETXWXlW/I5PSPjyU/D1I7LtC3oTushfuOjgbYzRMXySaw0dOpRSpUqxY8eO30+2sxnpnN46/9xR7/u2nDvqvd6TFH/0NQKvvMFxsYhkpTvuuIM1a9YQGRnJsmXLaNasGWfOnHGdleck/bidk6tmAFC0SS+CSlZwXCSSM2Xmivd+Y8zLxpgbvV4jksUKFCjApEmTABg4cCBnv/+cI7N6cmrtm9i0FEKur0GJf00j/I5WOupdJI+69dZbWbduHcWLF2ft2rU0atSI06dPu87KM1JjfyF20RiwHgpWvY8CFeq6ThLJsTIzeN8C/AC8bozZYozpZIzR7+El12jVqhUtW7YkISGBuI9Hkxb3G36FixPVfiiRrfvjFxbhOlFEvOyGG25g/fr1lC5dms2bN1O/fn2OHz/uOivXyzh7kph5w7CpSYSUr0nBmg+6ThLJ0S46eFtrE6y1M6y11YBngcHAEWPMW8aYcl4vFMkCkyZNomDBgueOeq/+wLmj3q+u5DpLRLJRuXLlWL9+PVdffTU7duygbt26HDt2zHVWruVJSyFmwQgy4mMIKH4dRZs+rTeki1xEptZ4G2NaGmMWAhOAV4CrgU+Az7zcJ5IlSpcuzZ49eyjZdSaFajyA8QtwnSQiDlx11VWsX7+e8uXL8+2331KnTh0OHTrkOivXsdZyfMkEUg/vxTc8kqi2L+LjH+g6SyTHy8xfTfcBrYCXrbW3WWtftdYes9bOA5Z6N08k6xQvXhzf0MKuM0TEsSuvvJJ169Zx88038/3331OrVi1+/fVX11m5yunP3yNxz3pMQDBR7QbhW0D/bBXJjEyt8bbWdrTWbvrzHdbanl5oEhER8aqoqCjWrFlDpUqV+Omnn6hZsyb79+93nZUrnPluDac3vQ/Gh8iWzxIQVdZ1kkiukZnjpF76i6NeTwPbrbUfZ32SiIiI9xUtWpRVq1bRtGlTNm3aRK1atVi5ciU33qhNvP5O8sHdHF8yAYDC9f5F8DV3OC4SyV0yc8U7CKjIuSUn+zi3y0kRoKMx5jUvtomIiHhVwYIFWbZsGXXq1OHIkSPUqVOHb775xnVWjpR26iixC0ZARjoFbmtGWKUWrpNEcp3MDN7lgHrW2knW2klAA+AGoA3QyJtxIiIi3lagQAE+++wz7r77bmJjY6lbty7btm1znZWjeFLOEjtvGJ6keILK3EaRBp34i9+Gi8hFZGbwvhII/cP3oUAJa20GkOKVKhERkWwUHBzMxx9/TKtWrTh58iT169fn888/d52VI1hPBrEfjSbt+G/4Fy1NZOv+OnBM5B/KzOA9FvjaGPOmMWYW8BUwzhgTCqz0ZpyIiEh2CQwM5MMPP+Tee+8lISGBu+++m9WrV7vOcspay4mV00n+5St8QgoSec8gfAJDL/5AEflLFxy8zbnfIy0HqgEfnf+oYa193Vp71lrbLxsaRUREsoW/vz/vvfcejz76KImJiTRr1owlS5a4znImYcciznz1Gfj6E9lmIP6FirlOEsnVLjh4W2st8JG19oi19mNr7UfW2sPZ1CYiIpLtfH19mTlzJp07dyY5OZlWrVqxcOFC11nZLvHHbZxc/QYAEU17EVTyBsdFIrlfZpaabDHGaL8gERHJN3x8fJg2bRq9evUiLS2N9u3bM2fOHNdZ2SY15mfiFo0F66Fg9Q6E3ljHdZJInpCZwbsu54bvH40xO40x3xpjdno7TERExCVjDOPHj2fAgAFkZGTw4IMPMmvWLNdZXpdx5iQx84ZhU5MIuaEWBas/4DpJJM/IzAE6TbxeISIikgMZY3jppZcICQnhxRdf5PHHHycpKYmuXbu6TvMKT1oKMQuGk5EQS0CJ64lo+rS2DRTJQhe94m2t/RUoxbm9vH8FEjPzOBERkbxi4MCBjBs3DoBu3boxfvx4x0VZz1oPxxePJ/XID/iGRxHVdiDGL8B1lkiectEB2hgzGHgOGHD+Jn/gHW9GiYiI5DR9+/Zl8uTJAPTp04eRI0c6Lspapza8S+LezzEBwUTdMwjf0MKuk0TynMxcuW4DtATOApzf1STMm1EiIiI5Uffu3XnjjTcwxjBw4EAGDhzIuQ3Acrczu1YRv/kDMD5EtupPQGQZ10kieVJmBu/U89sKWoDzB+eIiIjkS0888QTvvvsuvr6+jBw5kmeeeSZXD9/JB3ZxfMkkAIo06ETw1ZUcF4nkXZkZvOcaY6YDhYwxT3LutMoZ3s0SERHJuTp06MDcuXPx9/fn1VdfpXv37ng8HtdZlyzt5GFiF74EnnTCKrUg7PbmrpNE8rTMvLlyHDAPmA9cDwyy1k7ydpiIiEhO1rZtWxYuXEhgYCDTpk2jY8eOZGRkuM7KtIzkM8TMG4YnKZ6gqytRuN6/XCeJ5HmZ2U4Qa+0KYIWXW0RERHKVZs2asXjxYlq2bMmsWbNITk5m9uzZ+Pv7u067IJuRTtxHo0g/cRD/iKuIbPkcxsfXdZZInpeZXU3aGmP2GWNOG2PijTEJxpj47IgTERHJ6erXr8/SpUsJCwtjzpw53HfffaSkpLjO+lvWWk6siCb512/wCSlE1D2D8QkMcZ0lki9kZo33WKCltbagtTbcWhtmrQ33dpiIiEhuUbNmTVauXEmhQoVYuHAhbdq0ISkpyXXWX0rY9hFnvlkKvv5EtR2IX8Eo10ki+UZmBu9j1to9Xi8RERHJxe68807WrFlDREQES5YsoXnz5pw9e9Z11n9J3LeVk2tmAhDRrDeBV5Z3XCSSv2Rm8N5ujPnAGNPh/LKTtsaYtl4vExERyWUqVqzI2rVrKVasGKtXr+buu+8mPj5nrM5MPfYTcZ+8DFgK1niQ0BtquU4SyXcyM3iHc+6Y+EZAi/Mf2m9IRETkL1SoUIH169dTqlQpNm7cSIMGDThx4oTTpvSE48TMG4pNSyb0xjoUrHa/0x6R/Oqiu5pYax/PjhAREZG84tprr2X9+vXUq1ePbdu2UbduXVasWEFUVPavp/akJRO7YAQZZ44TeOUNFG3SE2NMtneISOZ2NbnOGLPKGLPr/Pe3GGMGej9NREQk9ypTpgzr16/nuuuuY+fOndSpU4fDhw9na4O1Ho5/+iqpR/fhV/AKItsOxPgFZGuDiPy/zCw1mQEMANIArLU7Af2OSkRE5CJKlizJ+vXruemmm9izZw+1a9fmt99+y7bXP7V+Nok/bMIEhBB5z2B8Qwpm22uLyP/KzOAdYq394k+3pXsjRkREJK+54oorWLNmDbfffjv79++nVq1a/Pjjj15/3TM7VxC/ZR4YHyJbDyAgorTXX1NELiwzg3ecMeYawAIYY+4Bjni1SkREJA+JiIhg1apVVKlShV9//ZVatWrx/fffe+31kn/7luPLpgBQpGEXgsve5rXXEpHMy8zg3R2YDpQ3xhwCnga6eLVKREQkjylUqBDLly+nVq1aHD58mNq1a/Ptt99m+euknThE7MKXwJNOWOVWhN3WNMtfQ0T+mYsO3tban6y1DYBIoLy1toa19lfvp4mIiOQtYWFhLFmyhIYNGxITE0OdOnXYsWNHlj1/RlICMfOH4UlOIPiaOyhc94kse24RuXyZueINgLX2rLU2wZsxIiIieV1ISAiLFi2iRYsWnDhxgnr16rFp06bLfl6bkUbsR6NIP3EI/8gyRLToh/HxzYJiEckqmR68RUREJGsEBQUxb9487rnnHuLj42nUqBFr1679x89nreXE8mmk/LYTn9BCRN0zCJ/AkKwLFpEsocFbRETEgYCAAN5//30eeughzp49S5MmTVi2bNk/eq74LxZwZudyjF8AUW1fxC88+w/qEZGL+9uTK40xbS/0QGvtgqzPERERyT/8/Px46623CA4OZsaMGbRs2ZIPP/yQli1bZvo5En/YzKm1swAo2qwPgSWu91KtiFyuCx0Z3+IC91lAg7eIiMhl8vHxYfr06QQHBzNx4kTatWvHu+++y7333nvRx6Yc3U/cp+MAS6FajxBavob3g0XkH/vbwdta+3h2hoiIiORXxhhee+01goODGTNmDB06dCA5OZlHHnnkbx+TnhBH7Pxh2LQUQm+qR/j/tXfv8VWWd97vPz9iCghqrYBbxQo4HDQBAoSoRSMjjtCKVtRp9VELWkWLtaN9wWyku7W11WE2eXxs6eMBKzha3UWxWB/BLZs9IjhFjlJRQVEOBaUe2uEMSsj1/JGQ4RBCEHKvED7v12u9knWt+/C773uFfLlyres++x8zrFjSF1Fbj3e1iLgYKACa7WxLKd1dX0VJknSkiQj+5V/+hebNm/PTn/6UwYMHs23bNoYOHbrXshWfb+OTZ3/Ojk1/o2nbAk7ofxsRkYOqJR2I/X64MiIeAr4N3AYE8I/AafVclyRJR5yI4K677uJf//VfAbj55pv55S9/udsyKVXw6QtlfP7R+xz15ZNoPWgUcVR+LsqVdIDqMqvJ11JK3wH+M6X0M+Ac4NT6LWvfImJARLwTEe9FxMhc1SFJUn3553/+Z8aOHQvA7bffzujRo6tfWzfjMbYue40mTVvQ5sqfkHf0cbkqU9IBqstQk61V1/94zgAAIABJREFUX7dExMnAX4H29VfSvkVEHvA/gX8A1gDzIuL5lNLbuahHkqT68v3vf59mzZoxdOhQ7rzzTrZs2cLGP61jw9zfQ5M8Wl12J/kn5KwfTNIXUJfg/UJEfBkYAyykckaT39RrVftWAryXUloOEBG/A74JGLwlSY3OjTfeSLNmzRg8eDA///nPq9u/8g/fo3m7ohxWJumL2G/wTint/El/NiJeAJqllNbXb1n7dAqwepfna4CzclSLJEn17tprr6VZs2ZcffXVlJeXc0zvyzimaECuy5L0Bew3eFcN77gYaLdz+YggpXRf/ZZWczk1tKW9FooYCgyFytvyFhcX13dde/nRr6dkvk/tX6vpd+W6BNVg8mvZ/4xKh5O8tufx8wee48FfDKPVjj8T/lsm7VdxccP7OanLUJP/BWwDFgMV9VvOfq1h9w92tgU+3HOhlNI4YBxAcXFxmj9/fjbV7WLyax9lvk/tXy7eC9o/f16k/evctZhVq1blugxJ+1Hb1J51Cd5tU0rdDl05B2Ue0DEi2gMfAFcB/y23JUmSJEn7V5fpBF+MiIvqvZI6SCmVA98HXgKWAE+nlN7KbVWSJEnS/tWlx/s1YHJENAG2UznOOqWUjq3XyvYhpTQVmJqLfUuqH4POPjHXJWgPDv+RpEOvLsH7v1N505zFKaW9PsgoSZIkaf/qMtRkGfCmoVuSJEn64urS470WmBERLwKf7WzM0XSCkiRJ0mGpLsF7RdXjS1UPSZIkSQeo1uBddfOclimlERnVI0mSJDVKtY7xTintAHpmVIskSZLUaNVlqMmiiHgeeAbYvLMxpfT7eqtKkiRJamTqEry/AvwVuGCXtgQYvCVJkqQ62m/wTildn0UhkiRJUmO233m8I6JTRPz/EfFm1fNuEfF/1X9pkiRJUuNRlxvoPALcSeXt4kkpvQFcVZ9FSZIkSY1NXYL30SmluXu0lddHMZIkSVJjVZfg/WlEnE7lByqJiCupvJulJEmSpDqqy6wmtwLjgC4R8QGVd7G8tl6rkiRJkhqZusxqshy4MCJaAE1SShvrvyxJkiSpcdlv8I6IpsAVQDvgqIgAIKV0d71WJkmSJDUidRlq8gdgPbAA+Kx+y5EkSZIap7oE77YppQH1XokkSZLUiNVlVpM/RkTXeq9EkiRJasT22eNddafKiqplro+I5VQONQkgpZS6ZVOiJEmSdPirbajJKUBRVoVIkiRJjVltwXtFSmlVZpVIkiRJjVhtwbtNRPxwXy+mlO6rh3okSZKkRqm24J0HtKRyTLckSZKkg1Bb8F7rTXIkSZKkQ6O26QTt6ZYkSZIOkdqCd7/MqpAkSZIauX0G75TS37IsRJIkSWrM6nLnSkmSJEkHyeAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZcDgLUmSJGXA4C1JkiRlwOAtSZIkZaDBBe+I+GlEfBARi6oe39jltTsj4r2IeCci+ueyTkmSJOlAHJXrAvbhf6SUynZtiIgzgauAAuBkYHpEdEop7chFgZIkSdKBaHA93rX4JvC7lNJnKaUVwHtASY5rkiRJkuqkoQbv70fEGxExPiKOr2o7BVi9yzJrqtokSZKkBi8nwTsipkfEmzU8vgk8CJwOFAFrgf++c7UaNpX2sf2hETE/IuZ/8skn9XIMkiRJ0oHIyRjvlNKFdVkuIh4BXqh6ugY4dZeX2wIf7mP744BxAMXFxTWGc0mSJClLDW6oSUSctMvTQcCbVd8/D1wVEU0joj3QEZibdX2SJEnSF9EQZzX5vyOiiMphJCuBmwFSSm9FxNPA20A5cKszmkiSJOlw0eCCd0rpulpeuwe4J8NyJEmSpEOiwQ01kSRJkhojg7ckSZKUAYO3JEmSlAGDtyRJkpQBg7ckSZKUAYO3JEmSlAGDtyRJkpQBg7ckSZKUAYO3JEmSlAGDtyRJkpQBg7ckSZKUgaNyXYAkqeEZdPaJuS5Bkhode7wlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMHJXrAiRJ0uHh448/Zvjw4SxdupSKiopclyPlTJMmTejSpQtlZWW0adOmzusZvCVJUp0MHz6cv//7v+fRRx8lPz8/1+VIObN9+3aeeOIJhg8fzuOPP17n9RxqIkmS6mTp0qVce+21hm4d8fLz87nuuutYunTpAa1n8JYkSXVSUVFh6Jaq5OfnH/CQK4O3JEk6rNxzzz0UFBTQrVs3ioqKmDNnDp9//jm33347p59+Oh07duSb3/wma9asqV7nL3/5C1dddRWnn346Z555Jt/4xjd49913c3gUOhI5xluSJB02Zs+ezQsvvMDChQtp2rQpn376KZ9//jmjRo1i48aNvPvuu+Tl5TFhwgQuv/xy5syZA8CgQYMYPHgwv/vd7wBYtGgRH330EZ06dcrl4egIY/CWJEmHjbVr19KqVSuaNm0KQKtWrdiyZQsTJkxgxYoV5OXlAXD99dczfvx4/v3f/52IID8/n1tuuaV6O0VFRTmpX0c2h5pIkqQDFhH18tifiy66iNWrV9OpUyeGDRvGK6+8wnvvvcdXv/pVjj322N2WLS4u5q233uLNN9+kV69e9XUqpDozeEuSpMNGy5YtWbBgAePGjaN169Z8+9vf5uWXX64xtKeU6hTmpaw41ESSJB2wlFLO9p2Xl0ffvn3p27cvXbt25eGHH2bVqlVs3LiRY445pnq5hQsXcskllwAwadKkXJUrVbPHW5IkHTbeeecdli1bVv180aJFdO7cmcGDB/PDH/6QHTt2APD444+zZcsWLrjgAi644AI+++wzHnnkker15s2bxyuvvJJ5/Tqy2eMtSZIOG5s2beK2225j3bp1HHXUUfzd3/0d48aN45hjjmH48OF06tSp+nbekydPrh5qMnnyZG6//XZGjx5Ns2bNaNeuHffff3+Oj0ZHGoO3JEk6bPTq1Ys//vGPNb42duxYxo4dW+NrJ598Mk8//XR9libtl0NNJEmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAwYvCVJkqQMGLwlSZKkDBi8JUmSpAw4j7ckSfpC2o2ccki3t3L0xYd0e/vTt29fysrKKC4uPqjtrFu3jqeeeophw4btc5kZM2ZQVlbGCy+8cFD72lOWx7Av9957L6NGjTqo/R9KQ4YMYeDAgVx55ZXVbS1btmTTpk1faHuH8vjs8ZYkSY3KztvGZ2XdunU88MADme7zUDuYY7j33nsPaPmUEhUVFV9oX7lwoMdXG4O3JEk6rFx22WX06tWLgoICxo0bB1T2aP7kJz/hrLPOYvbs2cybN4+vfe1rdO/enZKSEjZu3Fjjtp555hlKSkro1KkTs2bNAiqD+4gRI+jduzfdunXj4YcfBipvV9+vXz969uxJ165d+cMf/gDAyJEjef/99ykqKmLEiBH7rHvDhg0MGjSIM888k1tuuaU6fLZs2bJ6mUmTJjFkyJDq2goLC+nevTulpaX73G5Wx7B27VpKS0spKiqisLCQWbNmMXLkSLZu3UpRURHXXHMNAPfddx+FhYUUFhZy//33A7By5UrOOOMMhg0bRs+ePVm9ejXf+973KC4upqCggLvuuqt6P1OnTqVLly6ce+65/OAHP2DgwIEAbN68mRtuuIHevXvTo0eP6toP1JgxY6rPy677rel9VdPxHQyHmkiSpMPK+PHj+cpXvsLWrVvp3bs3V1xxBZs3b6awsJC7776bzz//nC5dujBx4kR69+7Nhg0baN68eY3bKi8vZ+7cuUydOpWf/exnTJ8+nUcffZTjjjuOefPm8dlnn9GnTx8uuugiTj31VCZPnsyxxx7Lp59+ytlnn82ll17K6NGjefPNN1m0aFGtdc+dO5e3336b0047jQEDBvD73/9+t+EQe7r77rt56aWXOOWUU1i3bt0+l8vqGJ566in69+/Pj370I3bs2MGWLVs477zz+PWvf1293oIFC5gwYQJz5swhpcRZZ53F+eefz/HHH88777zDhAkTqnvW77nnHr7yla+wY8cO+vXrxxtvvEGnTp24+eabmTlzJu3bt+fqq6+u3v8999zDBRdcwPjx41m3bh0lJSVceOGFtGjRYq9aR4wYwS9+8Yu92qdNm8ayZcuYO3cuKSUuvfRSZs6cSWlpaY3vq9GjR+92fAfL4C1Jkg4rv/rVr5g8eTIAq1evZtmyZeTl5XHFFVcA8M4773DSSSfRu3dvAI499th9buvyyy8HoFevXqxcuRKoDGdvvPEGkyZNAmD9+vUsW7aMtm3bMmrUKGbOnEmTJk344IMP+Oijj+pcd0lJCR06dADg6quv5tVXX601ePfp04chQ4bwrW99q7rOXB5D7969ueGGG9i+fTuXXXYZRUVFey3z6quvMmjQoOowfPnllzNr1iwuvfRSTjvtNM4+++zqZZ9++mnGjRtHeXk5a9eu5e2336aiooIOHTrQvn376vO0s/d52rRpPP/885SVlQGwbds2/vznP3PGGWfsVceYMWP2GuO9cxvTpk2jR48eQOVfAJYtW0ZpaWmN76sTTjihTuemrgzekiTpsDFjxgymT5/O7NmzOfroo+nbty/btm2jWbNm5OXlAZVjiCNir3Wvv/56Xn/9dU4++WSmTp0KQNOmTQHIy8ujvLy8ev2xY8fSv3//3dZ/7LHH+OSTT1iwYAH5+fm0a9eObdu21bn2PWva+XzX9l2399BDDzFnzhymTJlCUVERixYtYvjw4Tk7htLSUmbOnMmUKVO47rrrGDFiBN/5znd2WyaltM/1d+2ZXrFiBWVlZcybN4/jjz+eIUOGsG3btlrXTynx7LPP0rlz593aa7qutW3jzjvv5Oabb96tfV/vq0PNMd6SJOmwsX79eo4//niOPvpoli5dymuvvbbXMl26dOHDDz9k3rx5AGzcuJHy8nImTJjAokWL9hvO+vfvz4MPPsj27dsBePfdd9m8eTPr16+nTZs25Ofn8/LLL7Nq1SoAjjnmmH2OId/V3LlzWbFiBRUVFUycOJFzzz0XgBNPPJElS5ZQUVFR3eMK8P7773PWWWdx991306pVK1avXp3TY1i1ahVt2rThpptu4rvf/S4LFy4EID8/v3o/paWlPPfcc2zZsoXNmzczefJkzjvvvL22tWHDBlq0aMFxxx3HRx99xIsvvghUXrvly5dX99xPnDhxt2MaO3ZsdTh//fXXAep8TnZuY/z48dUznHzwwQd8/PHHtb6vdj2+g2WPtyRJ+kKynv4PYMCAATz00EN069aNzp077zZ0YacvfelLTJw4kdtuu42tW7fSvHlzpk+fvtuHGGtz4403snLlSnr27ElKidatW/Pcc89xzTXXcMkll1BcXExRURFdunQB4IQTTqBPnz4UFhby9a9/nTFjxtS43XPOOYeRI0eyePFiSktLGTRoEACjR49m4MCBnHrqqRQWFlaHwhEjRrBs2TJSSvTr14/u3bvX+TzVxzHMmDGDMWPGkJ+fT8uWLXn88ccBGDp0KN26daNnz548+eSTDBkyhJKSkuo6evToUR2kd+revTs9evSgoKCADh060KdPHwCaN2/OAw88wIABA2jVqlX1dgB+/OMfc/vtt9OtWzdSSrRr1+6Ap2e86KKLWLJkCeeccw5QOQTlt7/9ba3vqz2P72BEbV36jUFxcXGaP39+5vud/Frdx3wpO4POPjHXJUjSYau4uJhc/E7VkWXTpk20bNmSlBK33norHTt25I477sh1WTWq6WciIhaklGqcWN2hJpIkSWowHnnkEYqKiigoKGD9+vV7jcc+nDnURJIk6RBZvHgx11133W5tTZs2Zc6cOTmq6MDl+hjuuOOOBtvDfbAM3pIkSYdI165dD9mcz7nSGI6hoXKoiSRJkpSBnATviPjHiHgrIioioniP1+6MiPci4p2I6L9Le6+IWFz12q+ipgk6JUmSpAYqVz3ebwKXAzN3bYyIM4GrgAJgAPBARORVvfwgMBToWPUYkFm1kiRJ0kHKSfBOKS1JKb1Tw0vfBH6XUvospbQCeA8oiYiTgGNTSrNT5fyHjwOXZViyJEmSdFAa2ocrTwF2vQXVmqq27VXf79kuSZJy5FDfsyLrey307duXsrIyiotrnHK5ztatW8dTTz3FsGHDDnjde++9l1GjRh3U/g+lLM/JjBkzKCsrO+Cb4OzPkCFDGDhwIFdeeWV1W8uWLatvTHSgDuU1qrfgHRHTgf+jhpd+lFL6w75Wq6Et1dK+r30PpXJYCl/96lf3U2n98EYtkiTlxo4dO8jLy9v/gofIunXreOCBBzIJ3iklUko0adKw58c4mHPS0BzK4F1vVy2ldGFKqbCGx75CN1T2ZJ+6y/O2wIdV7W1raN/XvsellIpTSsWtW7c+mMOQJEkNzGWXXUavXr0oKChg3LhxQGWP5k9+8hPOOussZs+ezbx58/ja175G9+7dKSkpYePGjTVu65lnnqGkpIROnToxa9YsoDK4jxgxgt69e9OtWzcefvhhoPKOiv369aNnz5507dqVP/yhMtKMHDmS999/n6KiIkaMGFHjftauXUtpaSlFRUUUFhYya9YsRo4cydatWykqKuKaa64B4L777qOwsJDCwkLuv/9+AFauXMkZZ5zBsGHD6NmzJ6tXr+Z73/sexcXFFBQUcNddd1XvZ+rUqXTp0oVzzz2XH/zgBwwcOBCAzZs3c8MNN9C7d2969OhRXXsuzwnAhg0bGDRoEGeeeSa33HILFRUV1ddzp0mTJjFkyJDq2goLC+nevTulpaX73G5txowZU30cu567mt5XNV2jg7Lzf065eAAzgOJdnhcAfwKaAu2B5UBe1WvzgLOp7P1+EfhGXfbRq1evJEmSDt6ev1N/P/svh/RRV3/9619TSilt2bIlFRQUpE8//TQBaeLEiSmllD777LPUvn37NHfu3JRSSuvXr0/bt2/fazvnn39++uEPf5hSSmnKlCmpX79+KaWUHn744fTzn/88pZTStm3bUq9evdLy5cvT9u3b0/r161NKKX3yySfp9NNPTxUVFWnFihWpoKCg1prLysrSL37xi5RSSuXl5WnDhg0ppZRatGhRvcz8+fNTYWFh2rRpU9q4cWM688wz08KFC9OKFStSRKTZs2fvdQ7Ky8vT+eefn/70pz+lrVu3prZt26bly5enlFK66qqr0sUXX5xSSunOO+9MTzzxREoppf/8z/9MHTt2TJs2bcrpOXn55ZdT06ZN0/vvv5/Ky8vThRdemJ555pm9zsszzzyTBg8enFJKqbCwMK1Zs6b6OGoyePDg1K5du9S9e/fqx87tvfTSS+mmm25KFRUVaceOHeniiy9Or7zyym7ndNf31Z617KmmnAnMT/vIpTkZ4x0Rg4CxQGtgSkQsSin1Tym9FRFPA28D5cCtKaUdVat9D3gMaE5l8H4x+8olSVKu/epXv2Ly5MkArF69mmXLlpGXl8cVV1wBwDvvvMNJJ51E7969ATj22GP3ua3LL78cgF69erFy5UoApk2bxhtvvMGkSZMAWL9+PcuWLaNt27aMGjWKmTNn0qRJEz744AM++qhu49x79+7NDTfcwPbt27nssssoKiraa5lXX32VQYMG0aJFi+raZs2axaWXXsppp53G2WefXb3s008/zbhx4ygvL2ft2rW8/fbbVFRU0KFDB9q3bw/A1VdfXd1zO23aNJ5//nnKysoA2LZtG3/+858544wzcnZOAEpKSujQoUN1va+++upuY7P31KdPH4YMGcK3vvWt6jprMmbMmL3GeO88jmnTptGjRw+gssd+2bJllJaW1vi+OuGEE+p8LHWRk+CdUpoMTN7Ha/cA99TQPh8orOfSJElSAzZjxgymT5/O7NmzOfroo+nbty/btm2jWbNm1eO6U0rUdLuP66+/ntdff52TTz6ZqVOnApW3QgfIy8ujvLy8ev2xY8fSv3//3dZ/7LHH+OSTT1iwYAH5+fm0a9eObdu21anu0tJSZs6cyZQpU7juuusYMWIE3/nOd3ZbprKztGY7wzjAihUrKCsrY968eRx//PEMGTKEbdu21bp+Solnn32Wzp0779aey3MC7HWddj7ftX3X7T300EPMmTOHKVOmUFRUxKJFixg+fPhex1Dbebjzzju5+eabd2vf1/vqUGvYI/MlSZJ2sX79eo4//niOPvpoli5dymuvvbbXMl26dOHDDz9k3rx5AGzcuJHy8nImTJjAokWL9hvO+vfvz4MPPsj27dsBePfdd9m8eTPr16+nTZs25Ofn8/LLL7Nq1SoAjjnmmH2OId9p1apVtGnThptuuonvfve7LFy4EID8/Pzq/ZSWlvLcc8+xZcsWNm/ezOTJkznvvPP22taGDRto0aIFxx13HB999BEvvvhi9XEvX768upd64sSJux3T2LFjq8P566+/DpDTcwIwd+5cVqxYQUVFBRMnTuTcc88F4MQTT2TJkiVUVFRU90IDvP/++5x11lncfffdtGrVitWrV9f5GHYex/jx46tnOPnggw/4+OOPa31f7XqNDlZDm05QkiQdJnIxg9eAAQN46KGH6NatG507d95t+MVOX/rSl5g4cSK33XYbW7dupXnz5kyfPn23D+zV5sYbb2TlypX07NmTlBKtW7fmueee45prruGSSy6huLiYoqIiunTpAsAJJ5xAnz59KCws5Otf/zpjxozZa5szZsxgzJgx5Ofn07JlSx5//HEAhg4dSrdu3ejZsydPPvkkQ4YMoaSkpLqOHj16VAfpnbp3706PHj0oKCigQ4cO9OnTB4DmzZvzwAMPMGDAAFq1alW9HYAf//jH3H777XTr1o2UEu3atTugafzq45wAnHPOOYwcOZLFixdTWlrKoEGDABg9ejQDBw7k1FNPpbCwsDoojxgxgmXLlpFSol+/fnTv3r3OxwBw0UUXsWTJEs455xygcgjKb3/721rfV3teo4MRtf1ZojEoLi5O8+fPz3UZkiQd9oqLi/F3asO2adMmWrZsSUqJW2+9lY4dO3LHHXfkuqxGq6afiYhYkFKqcSJ0h5pIkiQ1Eo888ghFRUUUFBSwfv36vcYyK7ccaiJJknSILF68mOuuu263tqZNmzJnzpxM9n/HHXc0uB7uXJ+ThsTgLUmSdIh07dqVRYsW5bqMBsVz8l8caiJJkuqkSZMmh2x2B+lwt337dpo0ObAobfCWJEl10qVLF5544gnDt45427dv54knnqiexaWuHGoiSZLqpKysjOHDh/PQQw9RUVGR63KknGnSpAldunSpvhNoXRm8JUlSnbRp06Z6/mlJB86hJpIkSVIGDN6SJElSBgzekiRJUgYa/S3jI+ITYFUOdt0K+DQH+1XtvC4Nk9el4fGaNExel4bJ69Iw5eq6nJZSal3TC40+eOdKRMxPKRXnug7tzuvSMHldGh6vScPkdWmYvC4NU0O8Lg41kSRJkjJg8JYkSZIyYPCuP+NyXYBq5HVpmLwuDY/XpGHyujRMXpeGqcFdF8d4S5IkSRmwx1uSJEnKgMG7HkTEgIh4JyLei4iRua5HEBHjI+LjiHgz17WoUkScGhEvR8SSiHgrIv4p1zUJIqJZRMyNiD9VXZef5bomVYqIvIh4PSJeyHUt+i8RsTIiFkfEooiYn+t6BBHx5YiYFBFLq37HnJPrmnZyqMkhFhF5wLvAPwBrgHnA1Smlt3Na2BEuIkqBTcDjKaXCXNcjiIiTgJNSSgsj4hhgAXCZPyu5FREBtEgpbYqIfOBV4J9SSq/luLQjXkT8ECgGjk0pDcx1PaoUESuB4pSS83g3EBHxb8CslNJvIuJLwNEppXW5rgvs8a4PJcB7KaXlKaXPgd8B38xxTUe8lNJM4G+5rkP/JaW0NqW0sOr7jcAS4JTcVqVUaVPV0/yqhz00ORYRbYGLgd/kuhapIYuIY4FS4FGAlNLnDSV0g8G7PpwCrN7l+RoME1KtIqId0AOYk9tKBNVDGhYBHwP/X0rJ65J79wP/DFTkuhDtJQHTImJBRAzNdTGiA/AJMKFqaNZvIqJFrovayeB96EUNbfYWSfsQES2BZ4HbU0obcl2PIKW0I6VUBLQFSiLC4Vk5FBEDgY9TSgtyXYtq1Cel1BP4OnBr1dBG5c5RQE/gwZRSD2Az0GA+b2fwPvTWAKfu8rwt8GGOapEatKoxxM8CT6aUfp/rerS7qj/PzgAG5LiUI10f4NKqscS/Ay6IiN/mtiTtlFL6sOrrx8BkKoecKnfWAGt2+UvdJCqDeINg8D705gEdI6J91YD+q4Dnc1yT1OBUfYjvUWBJSum+XNejShHROiK+XPV9c+BCYGluqzqypZTuTCm1TSm1o/J3yr+nlK7NcVkCIqJF1YfDqRrOcBHg7Fk5lFL6C7A6IjpXNfUDGsyH9o/KdQGNTUqpPCK+D7wE5AHjU0pv5bisI15E/D9AX6BVRKwB7kopPZrbqo54fYDrgMWXkWiJAAADR0lEQVRV44kBRqWUpuawJsFJwL9VzdDUBHg6peT0dVLNTgQmV/YjcBTwVErp/81tSQJuA56s6gBdDlyf43qqOZ2gJEmSlAGHmkiSJEkZMHhLkiRJGTB4S5IkSRkweEuSJEkZMHhLkiRJGTB4S1IDEhEnRMSiqsdfIuKDXZ7/MaMavhwRw+p5Hz+tOra79/H6kxHxt4i4sj7rkKQsOZ2gJDVQEfFTYFNKqSzj/bYDXkgp1dut4utybBHxWFUdk+qrDknKkj3eknSYiIhNVV/7RsQrEfF0RLwbEaMj4pqImBsRiyPi9KrlWkfEsxExr+rRp4ZtFlSttygi3oiIjsBo4PSqtjFVy42o2sYbEfGzqrZ2EbE0Iv6tqn1SRBxd9droiHi7qj3T/zhIUkPlnSsl6fDUHTgD+BuVd2b7TUqpJCL+icq7tt0O/BL4HymlVyPiq1TeUfeMPbZzC/DLlNLOu7zlASOBwpRSEUBEXAR0BEqAAJ6PiFLgz0Bn4Lsppf+IiPHAsKqvg4AuKaW08xb0knSks8dbkg5P81JKa1NKnwHvA9Oq2hcD7aq+vxD4dUQsAp4Hjo2IY/bYzmxgVET8n8BpKaWtNezroqrH68BCoAuVQRxgdUrpP6q+/y1wLrAB2Ab8JiIuB7Yc1JFKUiNhj7ckHZ4+2+X7il2eV/Bf/7Y3Ac7ZR5gGIKX0VETMAS4GXoqIG6nsQd9VAP+SUnp4t8bKseB7flAopZTKI6IE6AdcBXwfuKC2g4mIl4ATgfkppRtrW1aSDlf2eEtS4zWNytALQEQU7blARHQAlqeUfkVlr3g3YCOwa8/4S8ANEdGyap1TIqJN1WtfjYhzqr6/Gni1arnjUkpTqRzystd+95RS6p9SKjJ0S2rM7PGWpMbrB8D/jIg3qPz3fiaVY7p39W3g2ojYDvwFuDul9LeI+I+IeBN4MaU0IiLOAGZHBMAm4FpgB7AEGBwRDwPLgAeB44A/REQzKnvL76jvA5Wkw4HTCUqSvpCDmXbQ6QQlHYkcaiJJyoVNwNDabqADnE/lhzQlqVGwx1uSJEnKgD3ekiRJUgYM3pIkSVIGDN6SJElSBgzekiRJUgYM3pIkSVIGDN6SJElSBv43ETnWhHDCBQ0AAAAASUVORK5CYII=\n",
"text/plain": [
"
"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plotter = ar.Plotter('clustered_results.json')\n",
"plotter.plot_operation('heat_storage', 'Heat', ylabel='Thermal energy [MWh]',\n",
" file_name='heat_storage', plot_single_period_with_index=7,\n",
" show_plot=True)\n",
"plotter.plot_operation('heat_bus', 'Heat', ylabel='Thermal energy [MWh]',\n",
" file_name='heat_bus', scale_to_hourly_resolution=True, \n",
" show_plot=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Two-stage model run\n",
"\n",
"The last part of this example shows the possibility of dividing the original problem into two levels to reduce computation time and enable the fast generation of several high-quality design solution candidates.\n",
"The developers of *aristopy* have published a comparable method in [this paper](https://doi.org/10.1115/IMECE2019-11519). \n",
"In the first optimization stage, a simplified problem is considered, and a design candidate is generated. The simplification is achieved by the aggregation of time series data and the relaxation of binary operation variables. In the second step, the design candidate's objective function value is determined for the original problem by considering the full-scale model.\n",
"\n",
"Building a large optimization model in Pyomo can take a considerable time. Therefore Pyomo offers the possibility to create so-called persistent solvers (only for Gurobi and CPLEX). This solver instance is created only once, and the model is added via its persistent interface. Later, this model can be solved various times without the time-consuming building phase. The flag \"declare_persistent\" must be set to True during model declaration to exploit these persistent models in *aristopy*. Minor changes can be added to the persistent solver at a later stage. More information about the persistent solver interface of Pyomo can be found [here](https://pyomo.readthedocs.io/en/stable/advanced_topics/persistent_solvers.html)."
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [],
"source": [
"# Construct (declare) a persistent model of the original (full scale) instance\n",
"original_model.declare_model(use_clustered_data=False, declare_persistent=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In the following steps, the cluster function is called for the energy system model (not the original/full-scale model!), and the model is declared with aggregated data. Then the binary operating variables (ON/OFF) are relaxed so that they can take on values between zero and one. The existence variables remain binary."
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [],
"source": [
"# Cluster time series data, declare the model and relax the integrality of the\n",
"# binary operation variables\n",
"es.cluster(number_of_typical_periods=20,\n",
" number_of_time_steps_per_period=24//hpts)\n",
"es.declare_model(use_clustered_data=True)\n",
"es.relax_integrality(include_time_dependent=True,\n",
" include_existence=False, include_modules=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Subsequently, several optimization problems are iteratively solved, and the results are stored in a pandas DataFrame for further comparison. In the first step, the problem is solved with clustered data (model is not declared again!). The configuration (existence and capacity of the components) is imported into the full-scale problem. Now, this original problem is solved with a fixed design as a pure dispatch model. \n",
"In the last step, an integer cut constraint is introduced in the aggregated model. Thus, the current design solution (the combination of binary variables) is excluded from the solution space, and another design solution must be generated in the next model run."
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {},
"outputs": [],
"source": [
"# Create a number of integer-cuts\n",
"nbr_of_icc = 5\n",
"\n",
"# Store some results in a DataFrame to compare solutions\n",
"results = pd.DataFrame(index=range(nbr_of_icc), columns=[\n",
" 'Obj. value [Mio.€]', 'Nbr. of engines [-]', \n",
" 'Capacity storage [MWh]', 'Capacity boiler [MW]', \n",
" '1st stage time [s]', '2nd stage time [s]'])\n",
"\n",
"for icc in range(nbr_of_icc):\n",
"\n",
" # Run optimization with clustered data\n",
" es.optimize(use_clustered_data=True, declare_model=False,\n",
" tee=False, results_file=f'1st_stage_results_{icc}.json')\n",
"\n",
" # Export the configuration to a DataFrame and store it in an Excel-File\n",
" config = es.export_component_configuration(\n",
" write_to_excel=True, file_name=f'1st_stage_config_{icc}.xlsx')\n",
"\n",
" # Import the optimal configuration and fix design variables\n",
" original_model.import_component_configuration(config)\n",
"\n",
" # Run the already declared persistent, full scale model\n",
" original_model.optimize(declare_model=False, tee=False,\n",
" results_file=f'2nd_stage_results_{icc}.json')\n",
"\n",
" # Reset variables to their original (unfixed) values\n",
" original_model.reset_component_variables()\n",
"\n",
" # Add an integer-cut constraint for the gas engines to the aggregated model\n",
" # => Exclude previously found optimal design from the solution space\n",
" es.add_design_integer_cut_constraint(which_instances=['gas_engine'])\n",
" \n",
" \n",
" results['Obj. value [Mio.€]'].loc[icc] = original_model.run_info['lower_bound'] / 1e6 \n",
" results['Nbr. of engines [-]'].loc[icc] = sum(comp.block.BI_EX.value \n",
" for comp in es.components.values() \n",
" if comp.group_name == 'gas_engine')\n",
" results['Capacity storage [MWh]'].loc[icc] = config['heat_storage'].loc['CAP']\n",
" results['Capacity boiler [MW]'].loc[icc] = config['gas_boiler'].loc['CAP']\n",
" results['1st stage time [s]'].loc[icc] = es.run_info['model_solve_time']\n",
" results['2nd stage time [s]'].loc[icc] = original_model.run_info['model_solve_time']"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The following summary of the results shows that a single iteration can be performed much faster by splitting it into two optimization levels. The detected solutions are comparable to the globally optimal solution of the original problem (see above). In addition, further design candidates are evaluated, and a deeper understanding of the present model is gained."
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"