diff --git a/docs/control/system_level_control/system_level_control_base.md b/docs/control/system_level_control/system_level_control_base.md index 7319cad36..9ae137c93 100644 --- a/docs/control/system_level_control/system_level_control_base.md +++ b/docs/control/system_level_control/system_level_control_base.md @@ -41,3 +41,13 @@ Helper functions for cost-aware controllers. :undoc-members: :show-inheritance: ``` + + +## Heterogeneous Commodity Hybrid System + +```{eval-rst} +.. autoclass:: h2integrate.control.control_strategies.system_level.system_level_control_base.HCHSConfig + :members: + :undoc-members: + :show-inheritance: +``` diff --git a/examples/35_system_level_control/complex_multi_commodity/driver_config.yaml b/examples/35_system_level_control/complex_multi_commodity/driver_config.yaml new file mode 100644 index 000000000..e6b823fec --- /dev/null +++ b/examples/35_system_level_control/complex_multi_commodity/driver_config.yaml @@ -0,0 +1,4 @@ +name: driver_config +description: This analysis runs a hybrid plant to match the first example in H2Integrate +general: + folder_output: outputs diff --git a/examples/35_system_level_control/complex_multi_commodity/plant_config.yaml b/examples/35_system_level_control/complex_multi_commodity/plant_config.yaml new file mode 100644 index 000000000..5702ddf9e --- /dev/null +++ b/examples/35_system_level_control/complex_multi_commodity/plant_config.yaml @@ -0,0 +1,111 @@ +name: plant_config +description: This plant is located in MN, USA... +sites: + site: + latitude: 32.31714 + longitude: -100.18 + resources: + wind_resource: + resource_model: WTKNLRDeveloperAPIWindResource + resource_parameters: + resource_year: 2013 + solar_resource: + resource_model: GOESAggregatedSolarAPI + resource_parameters: + resource_year: 2013 +# array of arrays containing left-to-right technology +# interconnections; can support bidirectional connections +# with the reverse definition. +# this will naturally grow as we mature the interconnected tech +technology_interconnections: + - [wind, combiner, electricity, cable] + - [solar, combiner, electricity, cable] + - [combiner, battery, electricity, cable] + # combine the battery output with the wind and solar generation + - [battery, elec_combiner, electricity, cable] + - [combiner, elec_combiner, electricity, cable] + # connect the electricity supply to the electrolyzer + - [elec_combiner, electrolyzer, electricity, cable] + - [electrolyzer, h2_storage, hydrogen, pipe] + # combine the h2 from the electrolyzer and the h2_storage + - [electrolyzer, h2_combiner, hydrogen, pipe] + - [h2_storage, h2_combiner, hydrogen, pipe] + # subtract the hydrogen supplied from the hydrogen demand + # - [h2_combiner, h2_load_demand, hydrogen, pipe] + # connect the hydrogen supplied as the demand to the ammonia model + - [h2_combiner, haber_bosch, hydrogen, pipe] + - [n2_feedstock, haber_bosch, nitrogen, pipe] + - [electricity_feedstock, haber_bosch, electricity, cable] + - [haber_bosch, nh3_load_demand, ammonia, pipe] + # etc +tech_to_dispatch_connections: + - [combiner, battery] + - [battery, battery] +resource_to_tech_connections: + # connect the wind resource to the wind technology + - [site.wind_resource, wind, wind_resource_data] + - [site.solar_resource, solar, solar_resource_data] +system_level_control: + control_strategy: DemandFollowingControl + demand_component: nh3_load_demand + control_parameters: + use_average_conversion_factor: true + solver_options: + solver_name: gauss_seidel + max_iter: 20 + convergence_tolerance: 1.0e-6 +plant: + plant_life: 30 +finance_parameters: + finance_groups: + finance_model: ProFastLCO + model_inputs: + params: + analysis_start_year: 2032 + installation_time: 36 # months + inflation_rate: 0.0 # 0 for real analysis + discount_rate: 0.06 # nominal return based on 2024 ATB baseline workbook for land-based wind + debt_equity_ratio: 0.724 # 2024 ATB uses 72.4% debt for land-based wind + property_tax_and_insurance: 0.025 # percent of CAPEX estimated based on https://www.nlr.gov/docs/fy25osti/91775.pdf https://www.house.mn.gov/hrd/issinfo/clsrates.aspx + total_income_tax_rate: 0.2574 # 0.257 tax rate in 2024 atb baseline workbook, value here is based on federal (21%) and state in MN (9.8) + capital_gains_tax_rate: 0.15 # H2FAST default + sales_tax_rate: 0.0 # average combined state and local sales tax https://taxfoundation.org/location/texas/ + debt_interest_rate: 0.07 # based on 2024 ATB nominal interest rate for land-based wind + debt_type: Revolving debt # can be "Revolving debt" or "One time loan". Revolving debt is H2FAST default and leads to much lower LCOH + loan_period_if_used: 0 # H2FAST default, not used for revolving debt + cash_onhand_months: 1 # H2FAST default + admin_expense: 0.00 # percent of sales H2FAST default + capital_items: + depr_type: MACRS # can be "MACRS" or "Straight line" + depr_period: 7 # 5 years - for clean energy facilities as specified by the IRS MACRS schedule https://www.irs.gov/publications/p946#en_US_2020_publink1000107507 + refurb: [0.] + cost_adjustment_parameters: + cost_year_adjustment_inflation: 0.025 # used to adjust modeled costs to target_dollar_year + target_dollar_year: 2022 + finance_subgroups: + h2: + commodity: hydrogen + commodity_stream: electrolyzer + technologies: [wind, solar, battery, electrolyzer, h2_storage] + nh3_produced: + commodity: ammonia + commodity_stream: haber_bosch + technologies: + - wind + - solar + - battery + - electrolyzer + - h2_storage + - haber_bosch + - n2_feedstock + nh3_delivered: + commodity: ammonia + commodity_stream: nh3_load_demand + technologies: + - wind + - solar + - battery + - electrolyzer + - h2_storage + - haber_bosch + - n2_feedstock diff --git a/examples/35_system_level_control/complex_multi_commodity/run_complex_multicommod.py b/examples/35_system_level_control/complex_multi_commodity/run_complex_multicommod.py new file mode 100644 index 000000000..af19119d8 --- /dev/null +++ b/examples/35_system_level_control/complex_multi_commodity/run_complex_multicommod.py @@ -0,0 +1,16 @@ +import os + +from h2integrate import EXAMPLE_DIR +from h2integrate.core.h2integrate_model import H2IntegrateModel + + +os.chdir(EXAMPLE_DIR / "35_system_level_control" / "complex_multi_commodity") + +################################## +# Create an H2I model with a fixed electricity load demand +h2i = H2IntegrateModel("top_level_config.yaml") + +h2i.setup() + +# Run the model +h2i.run() diff --git a/examples/35_system_level_control/complex_multi_commodity/tech_config.yaml b/examples/35_system_level_control/complex_multi_commodity/tech_config.yaml new file mode 100644 index 000000000..9c26c5058 --- /dev/null +++ b/examples/35_system_level_control/complex_multi_commodity/tech_config.yaml @@ -0,0 +1,279 @@ +name: technology_config +description: This hybrid plant produces ammonia +technologies: + wind: + performance_model: + model: PYSAMWindPlantPerformanceModel + cost_model: + model: ATBWindPlantCostModel + model_inputs: + performance_parameters: + num_turbines: 148 # number of turbines in the farm + turbine_rating_kw: 6000 + hub_height: 115 + rotor_diameter: 170 + create_model_from: default + config_name: WindPowerSingleOwner + pysam_options: !include pysam_options_6MW.yaml + + # hub_height: 115.0 # turbine hub-height + # operational_losses: 10.49 # percentage of non-wake losses + # floris_wake_config: !include "floris_v4_default_template.yaml" #floris wake model file + # floris_turbine_config: !include "floris_turbine_NREL_6MW_170.yaml" #turbine model file formatted for floris + # resource_data_averaging_method: average #"weighted_average", "average" or "nearest" + # operation_model: cosine-loss # turbine operation model + # default_turbulence_intensity: 0.06 + # enable_caching: true # whether to use cached results + # cache_dir: cache # directory to save or load cached data + layout: + layout_mode: basicgrid + layout_options: + row_D_spacing: 7.0 + turbine_D_spacing: 7.0 + rotation_angle_deg: 0.0 + row_phase_offset: 0.0 + layout_shape: square + cost_parameters: + capex_per_kW: 1380.0 + opex_per_kW_per_year: 29.0 + cost_year: 2019 + solar: + performance_model: + model: PYSAMSolarPlantPerformanceModel + cost_model: + model: ATBResComPVCostModel + model_inputs: + shared_parameters: + pv_capacity_kWdc: 400000 # 400 MWdc + performance_parameters: + dc_ac_ratio: 1.3 + create_model_from: default + config_name: PVWattsSingleOwner + tilt_angle_func: lat-func + pysam_options: + SystemDesign: + inv_eff: 96.0 + module_type: 0 # 19% efficiency + losses: 14.08 # dc losses + Lifetime: + dc_degradation: [0] + cost_parameters: + capex_per_kWdc: 1323 + opex_per_kWdc_per_year: 18 + cost_year: 2019 + combiner: + performance_model: + model: GenericCombinerPerformanceModel + dispatch_rule_set: + model: PyomoDispatchGenericConverter + model_inputs: + shared_parameters: + commodity: electricity + commodity_rate_units: kW + battery: + dispatch_rule_set: + model: PyomoRuleStorageBaseclass + control_strategy: + model: HeuristicLoadFollowingStorageController + performance_model: + model: PySAMBatteryPerformanceModel + cost_model: + model: ATBBatteryCostModel + model_inputs: + shared_parameters: + commodity: electricity + commodity_rate_units: kW + max_charge_rate: 96.0 # kW + max_capacity: 96.0 # kWh + init_soc_fraction: 0.9 + max_soc_fraction: 1.0 + min_soc_fraction: 0.2 + performance_parameters: + chemistry: LFPGraphite + demand_profile: 640000 # 640 MW + cost_parameters: + cost_year: 2019 + energy_capex: 310 # $/kWh from 2024 ATB year 2025 + power_capex: 311 # $/kW from 2024 ATB year 2025 + opex_fraction: 0.025 + control_parameters: + n_control_window_hours: 24 + system_commodity_interface_limit: 1e12 + elec_combiner: + performance_model: + model: GenericCombinerPerformanceModel + model_inputs: + performance_parameters: + commodity: electricity + commodity_rate_units: kW + electrolyzer: + performance_model: + model: ECOElectrolyzerPerformanceModel + cost_model: + model: SingliticoCostModel + model_inputs: + shared_parameters: + location: onshore + electrolyzer_capex: 1295 # $/kW overnight installed capital costs for a 1 MW system in 2022 USD/kW (DOE hydrogen program record 24005 Clean Hydrogen Production Cost Scenarios with PEM Electrolyzer Technology 05/20/24) (https://www.hydrogen.energy.gov/docs/hydrogenprogramlibraries/pdfs/24005-clean-hydrogen-production-cost-pem-electrolyzer.pdf?sfvrsn=8cb10889_1) + performance_parameters: + size_mode: normal + n_clusters: 16 + cluster_rating_MW: 40 + eol_eff_percent_loss: 10 # eol defined as x% change in efficiency from bol + uptime_hours_until_eol: 80000 # number of 'on' hours until electrolyzer reaches eol + include_degradation_penalty: true # include degradation + turndown_ratio: 0.1 # turndown_ratio = minimum_cluster_power/cluster_rating_MW + financial_parameters: + capital_items: + depr_period: 7 # based on PEM Electrolysis H2A Production Case Study Documentation estimate of 7 years. also see https://www.irs.gov/publications/p946#en_US_2020_publink1000107507 + replacement_cost_percent: 0.15 # percent of capex - H2A default case + h2_storage: + performance_model: + model: StoragePerformanceModel + control_strategy: + model: DemandOpenLoopStorageController + cost_model: + model: GenericStorageCostModel + model_inputs: + shared_parameters: + commodity: hydrogen + commodity_rate_units: kg/h + commodity_amount_units: kg + demand_profile: 9306.754158 # 50 kg/h + # performance_parameters: + min_soc_fraction: 0.0 + max_soc_fraction: 1.0 + charge_efficiency: 1.0 + discharge_efficiency: 1.0 + max_capacity: 1500.0 + init_soc_fraction: 0.0 + max_charge_rate: 500.0 + cost_parameters: + capacity_capex: 200.0 + charge_capex: 240.0 + opex_fraction: 0.05 + cost_year: 2020 + # # since the storage is being auto-sized by the performance model, + # # we set the sizing mode to 'auto' rather than defining the capacities + # # in the input file + # sizing_mode: auto # set as "auto" so storage capacity doesnt have to be defined + h2_combiner: + performance_model: + model: GenericCombinerPerformanceModel + model_inputs: + performance_parameters: + commodity: hydrogen + commodity_rate_units: kg/h + h2_load_demand: + performance_model: + model: GenericDemandComponent + model_inputs: + performance_parameters: + commodity: hydrogen + commodity_rate_units: kg/h + demand_profile: 9306.754158 # 50 kg/h + n2_feedstock: + performance_model: + model: FeedstockPerformanceModel + cost_model: + model: FeedstockCostModel + model_inputs: + shared_parameters: + commodity: nitrogen + commodity_rate_units: t/h + performance_parameters: + rated_capacity: 50.0 # metric tonnes of N2/hour + cost_parameters: + cost_year: 2022 + price: 5.0 + annual_cost: 0. + start_up_cost: 0.0 + electricity_feedstock: + performance_model: + model: FeedstockPerformanceModel + cost_model: + model: FeedstockCostModel + model_inputs: + shared_parameters: + commodity: electricity + commodity_rate_units: MW + performance_parameters: + rated_capacity: 29.0 # MW of electricity + cost_parameters: + cost_year: 2022 + price: 0.0 + annual_cost: 0. + start_up_cost: 0.0 + haber_bosch: + performance_model: + model: AmmoniaSynLoopPerformanceModel + cost_model: + model: AmmoniaSynLoopCostModel + model_inputs: # See converters/ammonia/Ammonia cost breakdown-ANL source.xlsx + shared_parameters: + production_capacity: 52777.6 + catalyst_consumption_rate: 0.000091295354067341 + catalyst_replacement_interval: 3 + performance_parameters: + size_mode: normal + capacity_factor: 0.9 + energy_demand: 0.530645243 + heat_output: 0.8299956 + feed_gas_t: 25.8 + feed_gas_p: 20 + feed_gas_x_n2: 0.25 + feed_gas_x_h2: 0.75 + feed_gas_mass_ratio: 1.13 + purge_gas_t: 7.5 + purge_gas_p: 275 + purge_gas_x_n2: 0.26 + purge_gas_x_h2: 0.68 + purge_gas_x_ar: 0.02 + purge_gas_x_nh3: 0.04 + purge_gas_mass_ratio: 0.07 + # --- Dynamic operating constraints (optional) --- + # Turndown ratio: minimum production as a fraction of rated capacity. + turndown_ratio: 0.2 + # Per-hour ramp limits as a fraction of rated capacity. + ramp_up_rate_fraction: 0.5 + ramp_down_rate_fraction: 0.5 + # Cold start: triggered after a long off-period; introduces a multi-hour delay. + include_cold_start: true + off_hours_cold_start: 6 + cold_start_delay_hours: 4 + # Warm start: triggered after any short off-period; introduces a sub-hour delay. + include_warm_start: true + off_hours_warm_start: 0.5 + warm_start_delay_hours: 0.5 + cost_parameters: + baseline_capacity: 52777.6 + base_cost_year: 2016 + capex_scaling_exponent: 0.6 + labor_scaling_exponent: 0.25 + asu_capex_base: 236920646 # See ASU + HB capex-NETL source.xlsx + synloop_capex_base: 302460908 # See ASU + HB capex-NETL source.xlsx + heat_capex_base: 7069100 + cool_capex_base: 4799200 + other_eqpt_capex_base: 0 + land_capex_base: 4112701.84103543 + deprec_noneq_capex_rate: 0.42 + labor_rate_base: 57 + num_workers_base: 50 + hours_yr: 2080 + gen_admin: 0.2 + prop_tax_ins: 0.02 + maint_rep: 0.005 + oxygen_byproduct_rate: 0.29405077250145 + water_consumption_rate: 0.049236824 + rebuild_cost_base: 0 + cooling_water_cost_base: 0.000113349938601175 + catalyst_cost_base: 23.19977341 + oxygen_price_base: 0.0285210891617726 + nh3_load_demand: + performance_model: + model: GenericDemandComponent + model_inputs: + performance_parameters: + commodity: ammonia + commodity_rate_units: kg/h + demand_profile: 47499.84 # kg/h diff --git a/examples/35_system_level_control/complex_multi_commodity/top_level_config.yaml b/examples/35_system_level_control/complex_multi_commodity/top_level_config.yaml new file mode 100644 index 000000000..e09f3dcda --- /dev/null +++ b/examples/35_system_level_control/complex_multi_commodity/top_level_config.yaml @@ -0,0 +1,4 @@ +name: H2Integrate_config +driver_config: driver_config.yaml +plant_config: plant_config.yaml +technology_config: tech_config.yaml diff --git a/examples/35_system_level_control/nh3_with_storage/driver_config.yaml b/examples/35_system_level_control/nh3_with_storage/driver_config.yaml new file mode 100644 index 000000000..e6b823fec --- /dev/null +++ b/examples/35_system_level_control/nh3_with_storage/driver_config.yaml @@ -0,0 +1,4 @@ +name: driver_config +description: This analysis runs a hybrid plant to match the first example in H2Integrate +general: + folder_output: outputs diff --git a/examples/35_system_level_control/nh3_with_storage/plant_config.yaml b/examples/35_system_level_control/nh3_with_storage/plant_config.yaml new file mode 100644 index 000000000..50e1650ee --- /dev/null +++ b/examples/35_system_level_control/nh3_with_storage/plant_config.yaml @@ -0,0 +1,129 @@ +name: plant_config +description: This plant is located in MN, USA... +sites: + site: + latitude: 32.31714 + longitude: -100.18 + resources: + wind_resource: + resource_model: WTKNLRDeveloperAPIWindResource + resource_parameters: + resource_year: 2013 + solar_resource: + resource_model: GOESAggregatedSolarAPI + resource_parameters: + resource_year: 2013 +# array of arrays containing left-to-right technology +# interconnections; can support bidirectional connections +# with the reverse definition. +# this will naturally grow as we mature the interconnected tech +technology_interconnections: + - [wind, combiner, electricity, cable] + - [solar, combiner, electricity, cable] + - [combiner, battery, electricity, cable] + # combine the battery output with the wind and solar generation + - [battery, elec_combiner, electricity, cable] + - [combiner, elec_combiner, electricity, cable] + # connect the electricity supply to the electrolyzer + - [elec_combiner, electrolyzer, electricity, cable] + - [electrolyzer, h2_storage, hydrogen, pipe] + # combine the h2 from the electrolyzer and the h2_storage + - [electrolyzer, h2_combiner, hydrogen, pipe] + - [h2_storage, h2_combiner, hydrogen, pipe] + # subtract the hydrogen supplied from the hydrogen demand + # - [h2_combiner, h2_load_demand, hydrogen, pipe] + # connect the hydrogen supplied as the demand to the ammonia model + - [h2_combiner, haber_bosch, hydrogen, pipe] + - [n2_feedstock, haber_bosch, nitrogen, pipe] + - [electricity_feedstock, haber_bosch, electricity, cable] + # ammonia system + - [haber_bosch, nh3_storage, ammonia, pipe] + - [nh3_storage, nh3_combiner, ammonia, pipe] + - [haber_bosch, nh3_combiner, ammonia, pipe] + - [nh3_combiner, nh3_load_demand, ammonia, pipe] + + # etc +tech_to_dispatch_connections: + - [combiner, battery] + - [battery, battery] +resource_to_tech_connections: + # connect the wind resource to the wind technology + - [site.wind_resource, wind, wind_resource_data] + - [site.solar_resource, solar, solar_resource_data] +system_level_control: + control_strategy: DemandFollowingControl + demand_component: nh3_load_demand + control_parameters: + use_average_conversion_factor: true + solver_options: + solver_name: gauss_seidel + max_iter: 20 + convergence_tolerance: 1.0e-6 +plant: + plant_life: 30 +finance_parameters: + finance_groups: + finance_model: ProFastLCO + model_inputs: + params: + analysis_start_year: 2032 + installation_time: 36 # months + inflation_rate: 0.0 # 0 for real analysis + discount_rate: 0.06 # nominal return based on 2024 ATB baseline workbook for land-based wind + debt_equity_ratio: 0.724 # 2024 ATB uses 72.4% debt for land-based wind + property_tax_and_insurance: 0.025 # percent of CAPEX estimated based on https://www.nlr.gov/docs/fy25osti/91775.pdf https://www.house.mn.gov/hrd/issinfo/clsrates.aspx + total_income_tax_rate: 0.2574 # 0.257 tax rate in 2024 atb baseline workbook, value here is based on federal (21%) and state in MN (9.8) + capital_gains_tax_rate: 0.15 # H2FAST default + sales_tax_rate: 0.0 # average combined state and local sales tax https://taxfoundation.org/location/texas/ + debt_interest_rate: 0.07 # based on 2024 ATB nominal interest rate for land-based wind + debt_type: Revolving debt # can be "Revolving debt" or "One time loan". Revolving debt is H2FAST default and leads to much lower LCOH + loan_period_if_used: 0 # H2FAST default, not used for revolving debt + cash_onhand_months: 1 # H2FAST default + admin_expense: 0.00 # percent of sales H2FAST default + capital_items: + depr_type: MACRS # can be "MACRS" or "Straight line" + depr_period: 7 # 5 years - for clean energy facilities as specified by the IRS MACRS schedule https://www.irs.gov/publications/p946#en_US_2020_publink1000107507 + refurb: [0.] + cost_adjustment_parameters: + cost_year_adjustment_inflation: 0.025 # used to adjust modeled costs to target_dollar_year + target_dollar_year: 2022 + finance_subgroups: + h2: + commodity: hydrogen + commodity_stream: electrolyzer + technologies: [wind, solar, battery, electrolyzer, h2_storage] + nh3_produced: + commodity: ammonia + commodity_stream: haber_bosch + technologies: + - wind + - solar + - battery + - electrolyzer + - h2_storage + - haber_bosch + - n2_feedstock + ammonia_available: + commodity: ammonia + commodity_stream: nh3_combiner + technologies: + - wind + - solar + - battery + - electrolyzer + - h2_storage + - haber_bosch + - n2_feedstock + - nh3_storage + nh3_delivered: + commodity: ammonia + commodity_stream: nh3_load_demand + technologies: + - wind + - solar + - battery + - electrolyzer + - h2_storage + - haber_bosch + - n2_feedstock + - nh3_storage diff --git a/examples/35_system_level_control/nh3_with_storage/run_nh3_with_storage.py b/examples/35_system_level_control/nh3_with_storage/run_nh3_with_storage.py new file mode 100644 index 000000000..eb0e4940d --- /dev/null +++ b/examples/35_system_level_control/nh3_with_storage/run_nh3_with_storage.py @@ -0,0 +1,16 @@ +import os + +from h2integrate import EXAMPLE_DIR +from h2integrate.core.h2integrate_model import H2IntegrateModel + + +os.chdir(EXAMPLE_DIR / "35_system_level_control" / "nh3_with_storage") + +################################## +# Create an H2I model with a fixed electricity load demand +h2i = H2IntegrateModel("top_level_config.yaml") + +h2i.setup() + +# Run the model +h2i.run() diff --git a/examples/35_system_level_control/nh3_with_storage/tech_config.yaml b/examples/35_system_level_control/nh3_with_storage/tech_config.yaml new file mode 100644 index 000000000..f985503c0 --- /dev/null +++ b/examples/35_system_level_control/nh3_with_storage/tech_config.yaml @@ -0,0 +1,303 @@ +name: technology_config +description: This hybrid plant produces ammonia +technologies: + wind: + performance_model: + model: PYSAMWindPlantPerformanceModel + cost_model: + model: ATBWindPlantCostModel + model_inputs: + performance_parameters: + num_turbines: 148 # number of turbines in the farm + turbine_rating_kw: 6000 + hub_height: 115 + rotor_diameter: 170 + create_model_from: default + config_name: WindPowerSingleOwner + pysam_options: !include pysam_options_6MW.yaml + layout: + layout_mode: basicgrid + layout_options: + row_D_spacing: 7.0 + turbine_D_spacing: 7.0 + rotation_angle_deg: 0.0 + row_phase_offset: 0.0 + layout_shape: square + cost_parameters: + capex_per_kW: 1380.0 + opex_per_kW_per_year: 29.0 + cost_year: 2019 + solar: + performance_model: + model: PYSAMSolarPlantPerformanceModel + cost_model: + model: ATBResComPVCostModel + model_inputs: + shared_parameters: + pv_capacity_kWdc: 400000 # 400 MWdc + performance_parameters: + dc_ac_ratio: 1.3 + create_model_from: default + config_name: PVWattsSingleOwner + tilt_angle_func: lat-func + pysam_options: + SystemDesign: + inv_eff: 96.0 + module_type: 0 # 19% efficiency + losses: 14.08 # dc losses + Lifetime: + dc_degradation: [0] + cost_parameters: + capex_per_kWdc: 1323 + opex_per_kWdc_per_year: 18 + cost_year: 2019 + combiner: + performance_model: + model: GenericCombinerPerformanceModel + dispatch_rule_set: + model: PyomoDispatchGenericConverter + model_inputs: + shared_parameters: + commodity: electricity + commodity_rate_units: kW + battery: + dispatch_rule_set: + model: PyomoRuleStorageBaseclass + control_strategy: + model: HeuristicLoadFollowingStorageController + performance_model: + model: PySAMBatteryPerformanceModel + cost_model: + model: ATBBatteryCostModel + model_inputs: + shared_parameters: + commodity: electricity + commodity_rate_units: kW + max_charge_rate: 96.0 # kW + max_capacity: 96.0 # kWh + init_soc_fraction: 0.9 + max_soc_fraction: 1.0 + min_soc_fraction: 0.2 + performance_parameters: + chemistry: LFPGraphite + demand_profile: 640000 # 640 MW + cost_parameters: + cost_year: 2019 + energy_capex: 310 # $/kWh from 2024 ATB year 2025 + power_capex: 311 # $/kW from 2024 ATB year 2025 + opex_fraction: 0.025 + control_parameters: + n_control_window_hours: 24 + system_commodity_interface_limit: 1e12 + elec_combiner: + performance_model: + model: GenericCombinerPerformanceModel + model_inputs: + performance_parameters: + commodity: electricity + commodity_rate_units: kW + electrolyzer: + performance_model: + model: ECOElectrolyzerPerformanceModel + cost_model: + model: SingliticoCostModel + model_inputs: + shared_parameters: + location: onshore + electrolyzer_capex: 1295 # $/kW overnight installed capital costs for a 1 MW system in 2022 USD/kW (DOE hydrogen program record 24005 Clean Hydrogen Production Cost Scenarios with PEM Electrolyzer Technology 05/20/24) (https://www.hydrogen.energy.gov/docs/hydrogenprogramlibraries/pdfs/24005-clean-hydrogen-production-cost-pem-electrolyzer.pdf?sfvrsn=8cb10889_1) + performance_parameters: + size_mode: normal + n_clusters: 16 + cluster_rating_MW: 40 + eol_eff_percent_loss: 10 # eol defined as x% change in efficiency from bol + uptime_hours_until_eol: 80000 # number of 'on' hours until electrolyzer reaches eol + include_degradation_penalty: true # include degradation + turndown_ratio: 0.1 # turndown_ratio = minimum_cluster_power/cluster_rating_MW + financial_parameters: + capital_items: + depr_period: 7 # based on PEM Electrolysis H2A Production Case Study Documentation estimate of 7 years. also see https://www.irs.gov/publications/p946#en_US_2020_publink1000107507 + replacement_cost_percent: 0.15 # percent of capex - H2A default case + h2_storage: + performance_model: + model: StoragePerformanceModel + control_strategy: + model: DemandOpenLoopStorageController + cost_model: + model: GenericStorageCostModel + model_inputs: + shared_parameters: + commodity: hydrogen + commodity_rate_units: kg/h + commodity_amount_units: kg + demand_profile: 9306.754158 # 50 kg/h + # performance_parameters: + min_soc_fraction: 0.0 + max_soc_fraction: 1.0 + charge_efficiency: 1.0 + discharge_efficiency: 1.0 + max_capacity: 1500.0 + init_soc_fraction: 0.0 + max_charge_rate: 500.0 + cost_parameters: + capacity_capex: 200.0 + charge_capex: 240.0 + opex_fraction: 0.05 + cost_year: 2020 + # # since the storage is being auto-sized by the performance model, + # # we set the sizing mode to 'auto' rather than defining the capacities + # # in the input file + # sizing_mode: auto # set as "auto" so storage capacity doesnt have to be defined + h2_combiner: + performance_model: + model: GenericCombinerPerformanceModel + model_inputs: + performance_parameters: + commodity: hydrogen + commodity_rate_units: kg/h + h2_load_demand: + performance_model: + model: GenericDemandComponent + model_inputs: + performance_parameters: + commodity: hydrogen + commodity_rate_units: kg/h + demand_profile: 9306.754158 # 50 kg/h + n2_feedstock: + performance_model: + model: FeedstockPerformanceModel + cost_model: + model: FeedstockCostModel + model_inputs: + shared_parameters: + commodity: nitrogen + commodity_rate_units: t/h + performance_parameters: + rated_capacity: 50.0 # metric tonnes of N2/hour + cost_parameters: + cost_year: 2022 + price: 5.0 + annual_cost: 0. + start_up_cost: 0.0 + electricity_feedstock: + performance_model: + model: FeedstockPerformanceModel + cost_model: + model: FeedstockCostModel + model_inputs: + shared_parameters: + commodity: electricity + commodity_rate_units: MW + performance_parameters: + rated_capacity: 29.0 # MW of electricity + cost_parameters: + cost_year: 2022 + price: 0.0 + annual_cost: 0. + start_up_cost: 0.0 + haber_bosch: + performance_model: + model: AmmoniaSynLoopPerformanceModel + cost_model: + model: AmmoniaSynLoopCostModel + model_inputs: # See converters/ammonia/Ammonia cost breakdown-ANL source.xlsx + shared_parameters: + production_capacity: 52777.6 + catalyst_consumption_rate: 0.000091295354067341 + catalyst_replacement_interval: 3 + performance_parameters: + size_mode: normal + capacity_factor: 0.9 + energy_demand: 0.530645243 + heat_output: 0.8299956 + feed_gas_t: 25.8 + feed_gas_p: 20 + feed_gas_x_n2: 0.25 + feed_gas_x_h2: 0.75 + feed_gas_mass_ratio: 1.13 + purge_gas_t: 7.5 + purge_gas_p: 275 + purge_gas_x_n2: 0.26 + purge_gas_x_h2: 0.68 + purge_gas_x_ar: 0.02 + purge_gas_x_nh3: 0.04 + purge_gas_mass_ratio: 0.07 + # --- Dynamic operating constraints (optional) --- + # Turndown ratio: minimum production as a fraction of rated capacity. + turndown_ratio: 0.2 + # Per-hour ramp limits as a fraction of rated capacity. + ramp_up_rate_fraction: 0.5 + ramp_down_rate_fraction: 0.5 + # Cold start: triggered after a long off-period; introduces a multi-hour delay. + include_cold_start: true + off_hours_cold_start: 6 + cold_start_delay_hours: 4 + # Warm start: triggered after any short off-period; introduces a sub-hour delay. + include_warm_start: true + off_hours_warm_start: 0.5 + warm_start_delay_hours: 0.5 + cost_parameters: + baseline_capacity: 52777.6 + base_cost_year: 2016 + capex_scaling_exponent: 0.6 + labor_scaling_exponent: 0.25 + asu_capex_base: 236920646 # See ASU + HB capex-NETL source.xlsx + synloop_capex_base: 302460908 # See ASU + HB capex-NETL source.xlsx + heat_capex_base: 7069100 + cool_capex_base: 4799200 + other_eqpt_capex_base: 0 + land_capex_base: 4112701.84103543 + deprec_noneq_capex_rate: 0.42 + labor_rate_base: 57 + num_workers_base: 50 + hours_yr: 2080 + gen_admin: 0.2 + prop_tax_ins: 0.02 + maint_rep: 0.005 + oxygen_byproduct_rate: 0.29405077250145 + water_consumption_rate: 0.049236824 + rebuild_cost_base: 0 + cooling_water_cost_base: 0.000113349938601175 + catalyst_cost_base: 23.19977341 + oxygen_price_base: 0.0285210891617726 + nh3_load_demand: + performance_model: + model: GenericDemandComponent + model_inputs: + performance_parameters: + commodity: ammonia + commodity_rate_units: kg/h + demand_profile: 47499.84 # kg/h + # below is only used in V2 + nh3_combiner: + performance_model: + model: GenericCombinerPerformanceModel + model_inputs: + performance_parameters: + commodity: ammonia + commodity_rate_units: kg/h + nh3_storage: + performance_model: + model: StoragePerformanceModel + control_strategy: + model: DemandOpenLoopStorageController + cost_model: + model: GenericStorageCostModel + model_inputs: + shared_parameters: + commodity: ammonia + commodity_rate_units: kg/h + commodity_amount_units: kg + demand_profile: 47499.84 # 50 kg/h + # performance_parameters: + min_soc_fraction: 0.0 + max_soc_fraction: 1.0 + charge_efficiency: 1.0 + discharge_efficiency: 1.0 + max_capacity: 15000.0 + init_soc_fraction: 0.0 + max_charge_rate: 5000.0 + cost_parameters: + capacity_capex: 200.0 + charge_capex: 240.0 + opex_fraction: 0.05 + cost_year: 2020 diff --git a/examples/35_system_level_control/nh3_with_storage/top_level_config.yaml b/examples/35_system_level_control/nh3_with_storage/top_level_config.yaml new file mode 100644 index 000000000..5f8defd20 --- /dev/null +++ b/examples/35_system_level_control/nh3_with_storage/top_level_config.yaml @@ -0,0 +1,4 @@ +name: H2Integrate_config +driver_config: driver_config.yaml +plant_config: plant_config_v2.yaml +technology_config: tech_config_v2.yaml diff --git a/h2integrate/control/control_strategies/system_level/demand_following_control.py b/h2integrate/control/control_strategies/system_level/demand_following_control.py index bb83712d0..03d6ebde8 100644 --- a/h2integrate/control/control_strategies/system_level/demand_following_control.py +++ b/h2integrate/control/control_strategies/system_level/demand_following_control.py @@ -1,10 +1,19 @@ +import warnings + import numpy as np +from attrs import field, define +from h2integrate.core.utilities import BaseConfig from h2integrate.control.control_strategies.system_level.system_level_control_base import ( SystemLevelControlBase, ) +@define(kw_only=True) +class DemandFollowingControlConfig(BaseConfig): + use_average_conversion_factor: bool = field(default=False) + + class DemandFollowingControl(SystemLevelControlBase): """Demand-following system-level controller. @@ -31,55 +40,175 @@ class DemandFollowingControl(SystemLevelControlBase): (each receives ``remaining_demand / n_dispatchable``). """ - def compute(self, inputs, outputs): - commodity = self.commodity - demand = inputs[self.demand_input_name].copy() + def setup(self): + super().setup() + + self.config = DemandFollowingControlConfig.from_dict( + self.options["plant_config"]["system_level_control"].get("control_parameters", {}), + strict=False, + ) + + self.tech_demands_set = [] + + def get_setpoints_for_commodity_subset( + self, inputs, outputs, commodity, commodity_demand, tech_subset: list | set | None = None + ): + # TODO: rename this method + if tech_subset is None: + tech_subset = set(self.input_techs) + + fixed_tech_subset = set(self.fixed_techs) & set(tech_subset) + flexible_tech_subest = set(self.flexible_techs) & set(tech_subset) + storage_tech_subset = set(self.storage_techs) & set(tech_subset) + dispatchable_tech_subset = set(self.dispatchable_techs) & set(tech_subset) # 1. Fixed techs: always produce, subtract from demand - for fixed_tech in self.fixed_techs: + for fixed_tech in fixed_tech_subset: commodity_from_tech = self._get_commodity_for_tech(fixed_tech) for tech_commodity in commodity_from_tech: if tech_commodity == commodity: - demand = self._subtract_fixed(fixed_tech, demand, commodity, inputs) + commodity_demand = self._subtract_fixed( + fixed_tech, commodity_demand, commodity, inputs + ) + self.tech_demands_set.append((fixed_tech, tech_commodity)) # 2. Flexible techs: operate at full production - for flexible_tech in self.flexible_techs: + for flexible_tech in flexible_tech_subest: commodity_from_tech = self._get_commodity_for_tech(flexible_tech) for tech_commodity in commodity_from_tech: if tech_commodity == commodity: - demand = self._subtract_flexible( - flexible_tech, demand, commodity, inputs, outputs + commodity_demand = self._subtract_flexible( + flexible_tech, commodity_demand, commodity, inputs, outputs ) + self.tech_demands_set.append((flexible_tech, tech_commodity)) else: if f"{flexible_tech}_rated_{tech_commodity}_production" in inputs: # set the per-tech set-point as the rated production outputs[f"{flexible_tech}_{tech_commodity}_set_point"] = inputs[ f"{flexible_tech}_rated_{tech_commodity}_production" ] * np.ones(self.n_timesteps) + self.tech_demands_set.append((flexible_tech, tech_commodity)) # 3. Storage dispatch # number of storage components that produce the demanded commodity n_storage = len( - [s for s in self.storage_techs if commodity in self._get_commodity_for_tech(s)] + [s for s in storage_tech_subset if commodity in self._get_commodity_for_tech(s)] ) - for storage_tech in self.storage_techs: + for storage_tech in storage_tech_subset: commodity_from_tech = self._get_commodity_for_tech(storage_tech) if commodity in commodity_from_tech: - demand = self._dispatch_storage( - storage_tech, demand / n_storage, commodity, inputs, outputs + commodity_demand = self._dispatch_storage( + storage_tech, commodity_demand / n_storage, commodity, inputs, outputs ) + self.tech_demands_set.append((storage_tech, commodity)) # 4. Dispatchable techs - remaining_demand = np.maximum(demand, 0.0) + remaining_demand = np.maximum(commodity_demand, 0.0) # calculate the number of dispatchable technologies that # produce the demanded commodity n_dispatchable = len( - [s for s in self.dispatchable_techs if commodity in self._get_commodity_for_tech(s)] + [s for s in dispatchable_tech_subset if commodity in self._get_commodity_for_tech(s)] ) - for dispatchable_tech in self.dispatchable_techs: + for dispatchable_tech in dispatchable_tech_subset: commodity_from_tech = self._get_commodity_for_tech(dispatchable_tech) if commodity in commodity_from_tech: outputs[f"{dispatchable_tech}_{commodity}_set_point"] = ( remaining_demand / n_dispatchable ) + self.tech_demands_set.append((dispatchable_tech, commodity)) + + return outputs + + def get_conversion_factors(self, converters, converter_upstreams, inputs): + conversion_factors = {} + for converter_info in list(converters): + input_cmod, tech, output_cmod = converter_info + tech_ancestors = converter_upstreams[(input_cmod, tech)] + conversion_ratio = self.get_converter_conversion_ratio( + inputs, input_cmod, output_cmod, tech, list(tech_ancestors) + ) + + has_nan = np.isnan(conversion_ratio).any() + has_inf = np.isinf(conversion_ratio).any() + is_zero = np.all(conversion_ratio == 0.0) + if has_inf or has_nan or is_zero: + # not all values are finite + if is_zero: + bad_indices = list(np.arange(0, len(conversion_ratio), 1)) + else: + inf_indices = np.argwhere(~np.isfinite(conversion_ratio)).flatten() + nan_indices = np.argwhere(~np.isnan(conversion_ratio)).flatten() + bad_indices = list(set(inf_indices) | set(nan_indices)) + + capacity_ratio = self.get_converter_capacity_ratio( + inputs, + input_cmod, + output_cmod, + tech, + list(tech_ancestors), + ) + conversion_ratio[bad_indices] = capacity_ratio + + if self.config.use_average_conversion_factor: + conversion_ratio = conversion_ratio.mean() + + conversion_factors[converter_info] = conversion_ratio + return conversion_factors + + def compute(self, inputs, outputs): + if not self.multi_commodity_system: + self.get_setpoints_for_commodity_subset( + inputs, outputs, self.commodity, inputs[self.demand_input_name].copy() + ) + return + + converter_conversion_factors = self.get_conversion_factors( + self.rename_me_config.converters, self.rename_me_config.converter_upstreams, inputs + ) + + conversion_factor_of_1 = ( + 1.0 if self.config.use_average_conversion_factor else np.ones(self.n_timesteps) + ) + + non_converter_conversion_factors = dict( + zip( + self.rename_me_config.non_converter_conversion_factor_keys, + [conversion_factor_of_1] + * len(self.rename_me_config.non_converter_conversion_factor_keys), + ) + ) + conversion_factors = non_converter_conversion_factors | converter_conversion_factors + + self.tech_demands_set = [] + + demand_techs = self.rename_me_config.converter_upstreams[(self.commodity, self.demand_tech)] + + outputs = self.get_setpoints_for_commodity_subset( + inputs, + outputs, + self.commodity, + inputs[self.demand_input_name].copy(), + tech_subset=demand_techs, + ) + + conversion_factors_tracker = {} + for recipe_name, recipe in self.rename_me_config.conversion_recipes.items(): + commodity_to_demand = recipe_name[1] + techs_to_demand = self.get_techs_to_demand_from_recipe(recipe_name) + conversion_factor = self.get_conversion_from_recipe(conversion_factors, recipe) + demand = inputs[self.demand_input_name].copy() * conversion_factor + outputs = self.get_setpoints_for_commodity_subset( + inputs, + outputs, + commodity_to_demand, + demand, + tech_subset=techs_to_demand, + ) + conversion_factors_tracker[recipe_name] = conversion_factor + unset_techs_cmods = self.techs_to_commodities - set(self.tech_demands_set) + unset_techs = [k for k in list(unset_techs_cmods) if k[0] not in self.feedstock_comps] + if unset_techs: + warnings.warn( + f"Commands not set for these technologies: {unset_techs}", UserWarning, stacklevel=3 + ) diff --git a/h2integrate/control/control_strategies/system_level/system_level_control_base.py b/h2integrate/control/control_strategies/system_level/system_level_control_base.py index 362b551ed..78c3da77d 100644 --- a/h2integrate/control/control_strategies/system_level/system_level_control_base.py +++ b/h2integrate/control/control_strategies/system_level/system_level_control_base.py @@ -1,3 +1,6 @@ +import warnings +import itertools + import numpy as np import networkx as nx import openmdao.api as om @@ -81,6 +84,149 @@ def _get_buy_price_default_and_shape(tech_config, tech_name, n_timesteps, plant_ return 0.0, n_timesteps +class HCHSConfig: + """Configuration class for a Heterogeneous Commodity Hybrid System. + The inputs to this configuration class are made in the + ``_post_setup_multi_commodity()`` method of ``SystemLevelControlBase``. + + Attributes: + converter_upstreams (dict): describes the technologies that provide + each of the input commodities for converter technologies + converters (set[tuple]): set of tuples describing the commodity conversions of + all technologies that convert one commodity into another + grouped_techs (dict): groups of technologies based on the commodities they produce. + This is built from ``converter_upstreams`` + simple_graph (nx.DiGraph): directional graph representation of ``grouped_techs`` + converter_tech_names (set[str]): names of the converter technologies in the system + conversion_recipes (dict): instructions on how to convert the demanded commodity + into a demand profile for each group of technologies + non_converter_conversion_factor_keys (set[tuple]): equivalent of ``converters`` + but for all other technologies in the system that are not in ``converters`` + + Examples: + Below highlights what these attributes look like if we have the following system: + + >>> technology_interconnections = [ + ... ["wind", "elec_combiner", "electricity", "cable"], + ... ["solar", "elec_combiner", "electricity", "cable"], + ... ["elec_combiner", "electrolyzer", "electricity", "cable"], + ... ["electrolyzer", "haber_bosch", "hydrogen", "pipe"], + ... ["electricity_feedstock", "haber_bosch", "electricity", "cable"], + ... ["haber_bosch", "nh3_storage", "ammonia", "pipe"], + ... ["haber_bosch", "nh3_combiner", "ammonia", "pipe"], + ... ["nh3_storage", "nh3_combiner", "ammonia", "pipe"], + ... ["nh3_combiner", "nh3_load_demand", "ammonia", "pipe"], + ... ] + + >>> converters # tuples formatted as (input_commodity, tech_name, output_commodity) + { + # (input_commodity, tech_name, output_commodity) + ("electricity", "electrolyzer", "hydrogen"), + ("electricity", "haber_bosch", "ammonia"), + ("hydrogen", "haber_bosch", "ammonia") + } + + >>> converter_upstreams # keys formatted as (input_commodity, tech) + { + # (input_commodity, tech): [upstream technologies providing input_commodity to tech] + ("electricity", "electrolyzer"): ["wind", "solar", "elec_combiner"], + ("electricity", "haber_bosch"): ["electricity_feedstock"], + ("hydrogen", "haber_bosch"): ["electrolyzer"], + ("ammonia", "nh3_load_demand"): ["nh3_combiner", "nh3_storage", "haber_bosch"] + } + + >>> converter_tech_names # set of strings + {"electrolyzer", "haber_bosch"} + + >>> non_converter_conversion_factor_keys + { + # (output_commodity, tech, output_commodity) + ("ammonia", "nh3_combiner", "ammonia"), + ("ammonia", "nh3_storage", "ammonia"), + ("electricity", "elec_combiner", "electricity"), + ("electricity", "wind", "electricity"), + ("electricity", "solar", "electricity"), + ("electricity", "electricity_feedstock", "electricity"), + } + + >>> grouped_techs + { + # group_name: [technologies in group] + "electricity-0": ["solar", "wind", "elec_combiner"], + "electricity-1": ["electricity_feedstock"], + "hydrogen-2": ["electrolyzer"], + "ammonia-3": ["nh3_combiner", "nh3_storage", "haber_bosch"] + } + + >>> list(conversion_recipes.keys(()) + [ + # (input_commodity, output_commodity, group_name) + ('ammonia', 'electricity', 'ammonia-3'), + ('ammonia', 'hydrogen', 'ammonia-3'), + ('hydrogen', 'electricity', 'hydrogen-2') + ] + + Recipe to calculate electricity demand for ammonia plant + + >>> conversion_recipes[("ammonia", "electricity", "ammonia-3")] + [ + [ + ('ammonia', 'nh3_combiner', 'ammonia'), + ('ammonia', 'nh3_storage', 'ammonia'), + ('electricity', 'haber_bosch', 'ammonia') + ] + ] + + Recipe to calculate hydrogen demand for ammonia plant + + >>> conversion_recipes[("ammonia", "hydrogen", "ammonia-3")] + [ + [ + ('ammonia', 'nh3_combiner', 'ammonia'), + ('ammonia', 'nh3_storage', 'ammonia'), + ('hydrogen', 'haber_bosch', 'ammonia') + ] + ] + + Recipe to calculate electricity demand for hydrogen used in the ammonia plant + + >>> conversion_recipes[("hydrogen", "electricity", "hydrogen-2")] + [ + # recipe of hydrogen to ammonia + [ + ('ammonia', 'nh3_combiner', 'ammonia'), + ('ammonia', 'nh3_storage', 'ammonia'), + ('hydrogen', 'haber_bosch', 'ammonia') + ], + # recipe of electricity to hydrogen + [ + ('electricity', 'electrolyzer', 'hydrogen'), + ('hydrogen', 'h2_combiner', 'hydrogen'), + ('hydrogen', 'h2_storage', 'hydrogen') + ] + ] + + """ + + def __init__( + self, + converter_upstreams: dict, + converters: set, + grouped_techs: dict, + simple_graph: nx.DiGraph, + converter_tech_names: set, + conversion_recipes: dict, + non_converter_keys: set, + ): + self.converter_upstreams = converter_upstreams + self.converters = converters + self.grouped_techs = grouped_techs + self.simple_graph = simple_graph + self.converter_tech_names = converter_tech_names + self.conversion_recipes = conversion_recipes + self.non_converter_conversion_factor_keys = non_converter_keys + + class SystemLevelControlBase(om.ExplicitComponent): """Base class for system-level controllers. @@ -202,6 +348,8 @@ def setup(self): self._setup_tech_category("storage", self.storage_techs) self._setup_feedstock_category(self.feedstock_comps) + self._post_setup_multi_commodity() + def _setup_commodity( self, tech_name, @@ -560,6 +708,8 @@ def _dispatch_storage(self, storage_tech, remaining_demand, commodity, inputs, o # Storage tech has its own sub-controller: emit a combined demand # signal (always positive) equal to the commodity flowing into # storage from upstream techs plus any remaining demand. + # TODO: possibly replace self.get_upstream_techs_for_commodity with + # get_successors_for_tech_with_input_cmod upstream_techs = self.get_upstream_techs_for_commodity(storage_tech, commodity) commodity_into_storage = np.zeros(self.n_timesteps) for tech_name in upstream_techs: @@ -821,6 +971,142 @@ def _feedstock_marginal_cost(self, inputs, marginal_cost_data): return np.full(self.n_timesteps, marginal_cost_scalar) + def _post_setup_multi_commodity(self): + """This method creates sets the attribute ``rename_me_config``, which is a + ``HCHSConfig`` object. This method is only used in + heterogeneous commodity hybrid system (HCHS). Below is a summary of what this method does: + + 1. Find the converter technologies and the technologies upstream of them. + This is done by calling the method ``_find_converter_techs()``. + + 2. Use the ``converter_upstreams`` made in Step 1 to group together technologies + with the same output commodity and the same downstream converter. A single + technology may exist in multiple groups if it has multiple commodities + connected to another component. + + 3. Create ``simple_graph`` - a directional graph representation of the + grouped technologies from Step 2. + + 4. Create keys for a conversion recipe based on the technologies that + are do not convert one commodity to another. + + 5. Create recipes to convert from the demand to the demand for each group + of technologies/commodities, using ``_make_conversion_factor_recipes()`` + + """ + if not self.multi_commodity_system: + return + # converter upstreams now has values of lists intead of sets + converters, converter_upstreams = self._find_converter_techs() + + demand_group_techs = self.get_successors_for_tech_with_input_cmod( + self.demand_tech, self.commodity + ) + + converter_upstreams[(self.commodity, self.demand_tech)] = demand_group_techs + # converter_info.add((self.commodity, self.demand_tech, self.commodity)) + # conversion fator recipes requires simple_graph, converters, demand_tech, grouped_techs + # grouped_techs[f"{self.commodity}-{len(converter_upstreams)+1}"] = demand_group_techs + # alt_grouped_techs[(self.commodity, f"{len(converter_upstreams)+1}")] = demand_group_techs + + # dictionary with keys as a group name mapping to a list of technologies in that groups + grouped_techs = {} + groups_to_commodities = {} + # Not doing this as 1-liners just in case any ordering could change (unlikely) + # for i, (key, value) in enumerate(converter_upstreams.items()): + # grouped_techs[f"{key[0]}-{i}"] = value + # groups_to_commodities[f"{key[0]}-{i}"] = key[0] + + for i, ((commodity, _), upstream_commodity_techs) in enumerate(converter_upstreams.items()): + group_name = f"{commodity}-{i}" + grouped_techs[group_name] = upstream_commodity_techs + groups_to_commodities[group_name] = commodity + + reversed_grouped_techs = {} + reversed_commodity_groups = {} + for group_name, techs_in_group in grouped_techs.items(): + group_commodity = groups_to_commodities[group_name] + + for tech_name in list(techs_in_group): + if (tech_name, group_commodity) in reversed_commodity_groups: + msg = ( + f"The tech/commodity pair {tech_name}/{group_commodity} " + "should not be duplicated. This may be due to a splitter " + "in the system which is not currently supported." + ) + # this error will get raised if using a splitter. + raise ValueError(msg) + + reversed_commodity_groups[(tech_name, group_commodity)] = group_name + + if tech_name in reversed_grouped_techs: + reversed_grouped_techs[tech_name] = reversed_grouped_techs[tech_name] + [ + group_name + ] + else: + reversed_grouped_techs[tech_name] = [group_name] + + simple_graph = nx.DiGraph() + for e in list(self.technology_graph.edges(data="commodity")): + source_tech, dest_tech, commodities = e # source_tech, dest_tech, commodity + if commodities is None or len(commodities) == 0: + # skip if no commodity is passed + continue + + for commod in commodities: + if (source_tech, commod) not in reversed_commodity_groups: + msg = ( + f"The technology/commodity pair {source_tech}/{commod} " "is not in a group" + ) + raise ValueError(msg) + # groups containing ``dest_tech`` + destination_groups = reversed_grouped_techs.get(dest_tech, [dest_tech]) + # group containing ``source_tech`` that output ``commod`` + source_group = reversed_commodity_groups[(source_tech, commod)] + for dest_group in destination_groups: + if dest_group == source_group: + # skip if in same group + continue + if not simple_graph.has_edge(source_group, dest_group): + # does not have edge + simple_graph.add_edge(source_group, dest_group, commodity=commod) + + else: + msg = ( + f"The edge for ({source_group}, {dest_group}, {commod}) " + "should not already exist." + ) + warnings.warn(msg, UserWarning, stacklevel=3) + + non_converter_keys = set() + converter_tech_names = {c[1] for c in converters} + + for converter_info, upstream_techs in converter_upstreams.items(): + input_cmod, _ = converter_info + non_converter_keys |= { + (input_cmod, t, input_cmod) for t in upstream_techs if t not in converter_tech_names + } + + non_converter_keys |= { + (self.commodity, t, self.commodity) + for t in demand_group_techs + if t not in converter_tech_names + } + + conversion_recipes = self._make_conversion_factor_recipes( + converters, simple_graph, grouped_techs + ) + + self.rename_me_config = HCHSConfig( + converter_upstreams, + converters, + grouped_techs, + simple_graph, + converter_tech_names, + conversion_recipes, + non_converter_keys, + ) + def get_upstream_techs_for_commodity( self, tech_name: str, commodity: str, include_feedstock_sources=True ): @@ -844,8 +1130,9 @@ def get_upstream_techs_for_commodity( if include_feedstock_sources: input_techs = self.input_techs | set(self.feedstock_comps) else: - input_techs = self.input_techs.copy() + input_techs = set(self.input_techs) + # TODO: refactor to call get_successors_for_tech_with_input_cmod # All graph ancestors of tech_name (any depth) ancestors = nx.ancestors(self.technology_graph, tech_name) @@ -860,79 +1147,454 @@ def get_upstream_techs_for_commodity( # Intersect with controller-managed techs return list(ancestors_with_commodity & input_techs) - def find_converter_techs(self, include_feedstock_sources=True): + def get_successors_for_tech_with_input_cmod(self, tech, input_commodity): + """Find technologies upstream of ``tech`` that produce ``input_commodity`` + for ``tech``. + + Args: + tech (str): Technology whose upstream suppliers are sought. + commodity (str): Commodity of interest that is an input commodity to ``tech`` + (e.g. ``"electricity"``). + + Returns: + list[str]: Controller-managed technologies upstream of ``tech`` + that produce ``commodity``. + """ + in_flows = dict(self.technology_graph.in_degree) + if in_flows[tech] < 1: + # Tech does not have any input commodiites + return [] + + successor_techs_with_commod = set() + upstream_techs = set(self.technology_graph.predecessors(tech)) + for upstream_tech in upstream_techs: + produces_cmod = False + if ( + commod := self.technology_graph.edges[upstream_tech, tech].get("commodity") + ) is not None: + if isinstance(commod, str) and commod == input_commodity: + # this if-statement is outdated and could be removed + successor_techs_with_commod.add(upstream_tech) + produces_cmod = True + if isinstance(commod, list) and input_commodity in commod: + successor_techs_with_commod.add(upstream_tech) + produces_cmod = True + if in_flows[upstream_tech] > 1 and produces_cmod: + # if only use >1, then it wouldn't catch splitters + # use `in_flows[upstream_tech] >= 1` to properly handle splitters + new_techs = self.get_successors_for_tech_with_input_cmod( + upstream_tech, input_commodity + ) + if new_techs: + successor_techs_with_commod |= set(new_techs) + + return list(successor_techs_with_commod) + + def _find_converter_techs(self): """Identify technologies that transform one commodity into another. A "converter" is a tech whose output commodities differ from the commodities produced by its upstream ancestors (e.g. an electrolyzer: electricity → hydrogen). + Returns: + tuple[set, dict]: 2-element tuple containing ``converters`` and ``converter_upstreams``. + + **converters** *(set[tuple])*: Set of tuples formatted as + ``(input_commodity, tech_name, output_commodity)`` tuples. + + >>> converters # tuples formatted as (input_commodity, tech_name, output_commodity) + { + # (input_commodity, tech_name, output_commodity) + ("electricity", "electrolyzer", "hydrogen"), + ("electricity", "haber_bosch", "ammonia"), + ("hydrogen", "haber_bosch", "ammonia") + } + + **converter_upstreams** *(dict[tuple[str,str], list[str]])*: Keys are set of + ``(input_commodity, tech_name)`` and the values are a set of + upstream technologies that output the `input_commodity` to `tech_name`. An + example of this variable is shown below. + + >>> converter_upstreams # keys formatted as (input_commodity, tech) + { + # (input_commodity, tech) : [techs that provide input_commodity to tech] + ("electricity", "electrolyzer"): ["wind", "solar", "elec_combiner"], + ("electricity", "haber_bosch"): ["electricity_feedstock"], + ("hydrogen", "haber_bosch"): ["electrolyzer"], + ("ammonia", "nh3_load_demand"): ["nh3_combiner", "nh3_storage", "haber_bosch"] + } + + """ + in_flows = dict(self.technology_graph.in_degree) + out_flows = dict(self.technology_graph.out_degree) + + non_converter_techs = [ + k for k in list(self.technology_graph.nodes) if in_flows[k] < 1 or out_flows[k] < 1 + ] + likely_converter_techs = ( + set(self.technology_graph.nodes) - set(non_converter_techs) - set(self.storage_techs) + ) & set(self.input_techs) + + converter_info = set() + converter_upstreams = {} + for converter in list(likely_converter_techs): + # predecessors are upstream and directly connected to the converter + predecessor_techs = set(self.technology_graph.predecessors(converter)) + # succesor techs are directly downstream of the converter + successor_techs = set(self.technology_graph.successors(converter)) + + input_commods = set() + for upstream_tech in predecessor_techs: + if ( + cmod := self.technology_graph.edges[upstream_tech, converter].get("commodity") + ) is not None: + if isinstance(cmod, str): + input_commods.add(cmod) + else: + input_commods |= set(cmod) + + output_commods = set() + for downstream_tech in list(successor_techs): + if ( + cmod := self.technology_graph.edges[converter, downstream_tech].get("commodity") + ) is not None: + if isinstance(cmod, str): + output_commods.add(cmod) + else: + output_commods |= set(cmod) + # A converter has commodities that appear only on one side: + # upstream-only commodities are consumed, output-only are produced. + consumed = input_commods - output_commods + produced = output_commods - input_commods + if consumed and produced: + # If both sides have unique commodities, this tech is a converter + + for input_commod in input_commods: + upstream_techs_with_commod = self.get_successors_for_tech_with_input_cmod( + converter, input_commod + ) + converter_upstreams[(input_commod, converter)] = upstream_techs_with_commod + for output_commod in output_commods: + converter_info.add((input_commod, converter, output_commod)) + + return converter_info, converter_upstreams + + def get_converter_capacity_ratio( + self, inputs, in_cmod, out_cmod, converter_tech, tech_ancestors + ): + """Get capacity ratio of ``in_cmod/out_cmod`` for technology ``converter_tech`` + Args: - include_feedstock_sources (bool, optional): If True, include feedstock techs - in the set of candidate technologies. Defaults to True. + inputs (dict): OpenMDAO inputs + in_cmod (str): commodity input to the ``converter_tech`` + out_cmod (str): commodity output from the ``converter_tech`` + converter_tech (str): name of the converter technologies + tech_ancestors (list[str] | set[str] | tuple[str]): upstream technologies + that produce ``in_cmod`` to the ``converter_tech`` Returns: - set[tuple[str, str, str]]: Set of ``(input_commodity, tech_name, output_commodity)`` - tuples for each detected conversion. Returns ``None`` for single-commodity systems. + float | np.ndarray: capacity ratio of ``in_cmod/out_cmod``. """ - if include_feedstock_sources: - input_techs = self.input_techs | set(self.feedstock_comps) - else: - input_techs = self.input_techs.copy() + rated_name_fmt = "{tech}_rated_{commod}_production" + feedstock_name_fmt = "{tech}_{commod}_out" + in_names = [rated_name_fmt.format(tech=t, commod=in_cmod) for t in list(tech_ancestors)] + in_feedstock_names = [ + feedstock_name_fmt.format(tech=t, commod=in_cmod) + for t in list(tech_ancestors) + if t in self.feedstock_comps + ] - # Single-commodity systems have no special handling by definition - if not self.multi_commodity_system: - return + total_in_cmod_capac = [inputs[n] for n in in_names if n in inputs] + avg_feedstock_capac = [inputs[n].mean() for n in in_feedstock_names if n in inputs] - converter_techs = set() - node_order = list(self.technology_graph.nodes()) - edges = list(self.technology_graph.edges(data="commodity")) + total_input_capac = np.array(total_in_cmod_capac).sum() + total_feedstock_capac = np.array(avg_feedstock_capac).sum() - # Track the most recently discovered converter so we can scope - # upstream searches for chained converters (A→B→C where B and C - # both convert). Without this, C would see A's commodity as upstream - # input even though B already consumed it. - last_converter = None + total_commodity_in_capacity = total_input_capac + total_feedstock_capac - for source_tech, _, _ in edges: - if source_tech not in input_techs: - continue + total_output_capac = inputs[rated_name_fmt.format(tech=converter_tech, commod=out_cmod)] + return total_commodity_in_capacity / total_output_capac[0] + + def get_converter_conversion_ratio( + self, inputs, in_cmod, out_cmod, converter_tech, tech_ancestors + ): + """Get conversion ratio of ``in_cmod/out_cmod`` for technology ``converter_tech`` + + Args: + inputs (dict): OpenMDAO inputs + in_cmod (str): commodity input to the ``converter_tech`` + out_cmod (str): commodity output from the ``converter_tech`` + converter_tech (str): name of the converter technologies + tech_ancestors (list[str] | set[str] | tuple[str]): upstream technologies + that produce ``in_cmod`` to the ``converter_tech`` + + Returns: + np.ndarray: conversion ratio of ``in_cmod/out_cmod``. + """ + input_name_fmt = "{tech}_{commod}_out" + in_names = [input_name_fmt.format(tech=t, commod=in_cmod) for t in list(tech_ancestors)] + total_in_cmod = [inputs[n] for n in in_names if n in inputs] + total_input = np.array(total_in_cmod).sum(axis=0) + total_output = inputs[input_name_fmt.format(tech=converter_tech, commod=out_cmod)] + + conversion_factor = total_input / np.abs(total_output) + return conversion_factor - # Get the commodities produced by this tech (the "output" side of the conversion) - output_commodities = set(self._get_commodity_for_tech(source_tech)) + def _make_recipe_from_grouped_path( + self, simple_graph, grouped_techs, converter_tech_names, path + ): + compounding_conversion_factor_recipes = {} + + reverse_path = path[::-1] + commodity_conversions = [ + simple_graph.edges[p0, p1].get("commodity", None) + for p0, p1 in zip(reverse_path[1:], reverse_path[:-1]) + ] + commodity_nodes = list(itertools.pairwise(commodity_conversions)) + techs = reverse_path[1:] - # Find controlled ancestors of this tech - all_ancestors = nx.ancestors(self.technology_graph, source_tech) & input_techs + commodity_graph = nx.DiGraph() # nodes are commodities + for i, commod_node in enumerate(commodity_nodes): + # ammonia, hydrogen + down_cmod, up_cmod = commod_node + commodity_graph.add_edge(down_cmod, up_cmod, tech=techs[i]) - if last_converter is not None: - # Only consider ancestors that appear after the last converter - # in topological order, preventing double-counting across - # chained converters. - converter_idx = node_order.index(last_converter) - nodes_after_converter = set(node_order[converter_idx + 1 :]) - ancestors = all_ancestors & nodes_after_converter + commodity_edges = commodity_graph.edges(data="tech") + + path_recipe = [] + + for edge in commodity_edges: + # in_cmod is demand of next tech + out_cmod, in_cmod, tech = edge + if tech in grouped_techs: + techs_in_group = list(grouped_techs[tech]) + + recipe = [] + for t in techs_in_group: + if t in converter_tech_names: + recipe.append((in_cmod, t, out_cmod)) + else: + recipe.append((out_cmod, t, out_cmod)) + # TODO: add check if any other non-converter techs have a non-1 conversion factor else: - ancestors = all_ancestors + recipe = [(in_cmod, tech, out_cmod)] - # Keep only ancestors actually connected (reachable) to this tech - connected_ancestors = [ - t for t in ancestors if nx.has_path(self.technology_graph, t, source_tech) + path_recipe.append(recipe) + compounding_conversion_factor_recipes[(out_cmod, in_cmod, tech)] = path_recipe.copy() + return compounding_conversion_factor_recipes + + def _make_conversion_factor_recipes( + self, converters, simple_graph, grouped_techs, use_complex_keys=False + ): + """Make recipes to for compounding conversion factor calculations. + + Args: + converters (set[tuple]): Set of tuples formatted as + ``(input_commodity, tech_name, output_commodity)``. An + example of this variable is shown below: + + >>> converters + { + # (input_commodity, tech_name, output_commodity) + ("electricity", "electrolyzer", "hydrogen"), + ("electricity", "haber_bosch", "ammonia"), + ("hydrogen", "haber_bosch", "ammonia") + } + + simple_graph (nx.DiGraph): graph representing the connections + of the technology groups in ``grouped_techs`` + grouped_techs (dict): dictionary with keys as the group name + and values of the technologies within that group. + + >>> grouped_techs + { + "electricity-0": ["solar", "wind", "elec_combiner"], + "electricity-1": ["electricity_feedstock"], + "hydrogen-2": ["electrolyzer"], + "ammonia-3": ["nh3_combiner", "nh3_storage", "haber_bosch"] + } + + use_complex_keys (bool, optional): If True, use key names formatted as + ``(output_commodity, input_commodity, (i, converter_tech_group))``. + Defaults to False. + + Returns: + dict[tuple[str], list[list]]: recipes to calculate the conversion ratio from + the demand commodity to all upstream subsystems. Keys are the recipe name, which + are tuples ``(output_commodity, input_commodity, converter_tech_group)``. + + Values are embedded lists. Each list defines the technologies in a + step of the conversion. Each element of a list is a tuple formatted as + ``(input_commodity, technology, output_commodity)``. An example is shown below. + + >>> conversion_recipes[("hydrogen", "electricity", "hydrogen-2")] + [ + [ + ('ammonia', 'nh3_combiner', 'ammonia'), + ('ammonia', 'nh3_storage', 'ammonia'), + ('hydrogen', 'haber_bosch', 'ammonia') + ], + [ + ('electricity', 'electrolyzer', 'hydrogen'), + ('hydrogen', 'h2_combiner', 'hydrogen'), + ('hydrogen', 'h2_storage', 'hydrogen') + ] ] - # Gather all commodities produced by connected ancestors - input_commodities = set() - for ancestor in connected_ancestors: - input_commodities.update(self._get_commodity_for_tech(ancestor)) + Note: + For more complext system architecturs, the conversion recipe keys + may be formatted as + ``(output_commodity, input_commodity, (i, converter_tech_group))`` + where ``i`` is a unique number to distinguish the recipe. - # A converter has commodities that appear only on one side: - # upstream-only commodities are consumed, output-only are produced. - consumed = input_commodities - output_commodities - produced = output_commodities - input_commodities + """ + if not self.multi_commodity_system: + return {} + + converter_tech_names = {v[1] for v in list(converters)} + + # 6. Get the compounding conversion factors + in_degs = dict(simple_graph.in_degree) + starting_techs = {k for k, v in in_degs.items() if v == 0} + needs_complex_keys = False + compounding_conversion_factor_recipes = {} + cnt = 0 + for starting_tech in list(starting_techs): + paths = list(nx.all_simple_paths(simple_graph, starting_tech, self.demand_tech)) + for path in paths: + res = self._make_recipe_from_grouped_path( + simple_graph, grouped_techs, converter_tech_names, path + ) - # If both sides have unique commodities, this tech is a converter - if consumed and produced: - for in_comm in consumed: - for out_comm in produced: - converter_techs.add((in_comm, source_tech, out_comm)) - last_converter = source_tech + if duplicate_recipe := set(res) & set(compounding_conversion_factor_recipes): + mismatched_recipes = [] + + for recipe_name in list(duplicate_recipe): + recipe_1 = compounding_conversion_factor_recipes[recipe_name] + recipe_2 = res[recipe_name] + if len(recipe_1) != len(recipe_2): + mismatched_recipes.append(recipe_name) + + continue + # have the same length of recipes + for i in range(len(recipe_1)): + if len(recipe_1[i]) != len(recipe_2[i]): + mismatched_recipes.append(recipe_name) + + continue + if set(recipe_1[i]) != set(recipe_2[i]): + mismatched_recipes.append(recipe_name) + + if mismatched_recipes: + needs_complex_keys = True + + # if needs_complex_keys: + # break + + if use_complex_keys: + new_res = { + (k[0], k[1], (cnt + i, k[2])): v for i, (k, v) in enumerate(res.items()) + } + compounding_conversion_factor_recipes |= new_res + cnt += len(res) + else: + compounding_conversion_factor_recipes |= res + + if needs_complex_keys and use_complex_keys: + warnings.warn( + "Duplicate recipes still exist with complex keys", UserWarning, stacklevel=3 + ) + + if needs_complex_keys and not use_complex_keys: + compounding_conversion_factor_recipes = self._make_conversion_factor_recipes( + converters, simple_graph, grouped_techs, use_complex_keys=True + ) + return compounding_conversion_factor_recipes + + return compounding_conversion_factor_recipes - return converter_techs + def get_techs_to_demand_from_recipe(self, recipe_name): + """Get a list of technologies that are in a subsystem that + outputs ``input_commodity`` to the ``tech_group_name``. + + Args: + recipe_name (tuple): name of recipe formatted as a tuple of + ``(input_commodity, output_commodity, tech_group_name)`` or + ``(input_commodity, output_commodity, (i,tech_group_name))``. + This should be a key from the dictionary returned from + ``_make_conversion_factor_recipes()``. + + Raises: + ValueError: there are multiple techs + + Returns: + list[str]: list of technologies that output the ``input_commodity`` + and are connected upstream of ``tech_group_name`` + """ + _, input_cmod, tech_group = recipe_name + if isinstance(tech_group, tuple): + _, tech_group = tech_group + + techs_to_demand = [ + s + for s in list(self.rename_me_config.simple_graph.predecessors(tech_group)) + if self.rename_me_config.simple_graph.edges[s, tech_group].get("commodity", "") + == input_cmod + ] + if len(techs_to_demand) != 1: + raise ValueError("Unexpected situation!") + if techs_to_demand[0] in self.rename_me_config.grouped_techs: + techs_in_group = list(self.rename_me_config.grouped_techs[techs_to_demand[0]]) + else: + techs_in_group = techs_to_demand[0] + return techs_in_group + + def get_conversion_from_recipe(self, conversion_factors, recipe): + """Get the conversion factor from a recipe. + + Args: + conversion_factors (dict): dictionary with keys of 3 element tuples + formatted as ``(input_commodity, tech, output_commodity)``. + Values are an array or float of the conversion factor + ``input_commodity/output_commodity``. An example is shown below: + + >>> conversion_factors + { + ('electricity','electrolyzer','hydrogen'): 55.5, + ('water','electrolyzer','hydrogen'): 40.0, + ('electricity','electrolyzer','oxygen'): 60.0 + } + + recipe (list[list[tuples]]): embedded list of conversions, + a value from from the ``conversion_recipes`` attribute. + This should be a value from the dictionary returned from + ``_make_conversion_factor_recipes()``. + + >>> recipe + [ + [ + ('ammonia', 'nh3_combiner', 'ammonia'), + ('ammonia', 'nh3_storage', 'ammonia'), + ('hydrogen', 'haber_bosch', 'ammonia') + ], + [ + ('electricity', 'electrolyzer', 'hydrogen'), + ('hydrogen', 'h2_combiner', 'hydrogen'), + ('hydrogen', 'h2_storage', 'hydrogen') + ] + ] + + + Returns: + float | np.ndarray: conversion factor created from the recipe. + + + """ + path_conversion = 1.0 + # TODO: update to handle more complex systems + # (or maybe do it external to this method) + for path in recipe: + for tech_conversion in path: + path_conversion *= conversion_factors.get(tech_conversion, 1.0) + + return path_conversion diff --git a/h2integrate/control/control_strategies/system_level/test/test_slc_baseclass.py b/h2integrate/control/control_strategies/system_level/test/test_slc_baseclass.py new file mode 100644 index 000000000..c2eb27397 --- /dev/null +++ b/h2integrate/control/control_strategies/system_level/test/test_slc_baseclass.py @@ -0,0 +1,669 @@ +import numpy as np +import pytest +import openmdao.api as om + +from h2integrate import EXAMPLE_DIR, H2IntegrateModel +from h2integrate.core.inputs.validation import load_tech_yaml, load_plant_yaml +from h2integrate.control.control_strategies.system_level.system_level_control_base import ( + SystemLevelControlBase, +) + + +def make_tech_classifiers(tech_list): + fixed_techs = [] + flexible_techs = ["wind", "solar", "desalination"] + dispatchable_techs = ["electrolyzer", "haber_bosch", "natural_gas_plant", "grid_buy", "grid"] + storage_techs = ["battery", "h2_storage", "nh3_storage"] + feedstock_techs = ["ng_feedstock", "n2_feedstock", "electricity_feedstock", "ocean"] + classifiers = {k: "flexible" for k in flexible_techs} + classifiers |= {k: "dispatchable" for k in dispatchable_techs} + classifiers |= {k: "storage" for k in storage_techs} + classifiers |= {k: "feedstock" for k in feedstock_techs} + classifiers |= {k: "fixed" for k in fixed_techs} + + classifiers |= {k: "connector" for k in tech_list if "combiner" in k} + classifiers |= {k: "splitter" for k in tech_list if "splitter" in k} + classifiers |= {k: "feedstock" for k in tech_list if "feedstock" in k} + classifiers |= {k: "demand" for k in tech_list if "demand" in k} + + classified_techs = list(set(tech_list) & set(classifiers)) + tech_control_classifiers = {k: classifiers[k] for k in classified_techs} + return tech_control_classifiers + + +def make_slc_topology(plant_config, tech_config): + model = object.__new__(H2IntegrateModel) + model.slc = True + # plant_config["system_level_control"].pop("demand_component") + model.plant_config = plant_config + + tech_control_classifiers = make_tech_classifiers(list(tech_config["technologies"])) + model.tech_control_classifiers = tech_control_classifiers + model.technology_config = tech_config + model.technology_graph = model.create_technology_graph( + plant_config.get("technology_interconnections", {}) + ) + slc_topology = model._classify_slc_technologies() + return slc_topology + + +def make_and_setup_slc_baseclass(plant_config, tech_config) -> SystemLevelControlBase: + slc_config = make_slc_topology(plant_config, tech_config) + slc = object.__new__(SystemLevelControlBase) + # run the start of setup() + slc.n_timesteps = plant_config["plant"]["simulation"]["n_timesteps"] + slc.commodity = slc_config["demand_commodity"] + slc.commodity_rate_units = slc_config.get("demand_commodity_rate_units", None) + slc.demand_tech = slc_config["demand_tech"] + slc.storage_techs_to_control = slc_config.get("storage_techs_to_control", {}) + slc.technology_graph = slc_config["technology_graph"] + slc.fixed_techs = [k for k, v in slc_config["tech_control_classifiers"].items() if v == "fixed"] + slc.flexible_techs = [ + k for k, v in slc_config["tech_control_classifiers"].items() if v == "flexible" + ] + slc.dispatchable_techs = [ + k for k, v in slc_config["tech_control_classifiers"].items() if v == "dispatchable" + ] + slc.storage_techs = [ + k for k, v in slc_config["tech_control_classifiers"].items() if v == "storage" + ] + slc.feedstock_comps = [ + k for k, v in slc_config["tech_control_classifiers"].items() if v == "feedstock" + ] + + slc.input_techs = set( + slc.fixed_techs + slc.flexible_techs + slc.dispatchable_techs + slc.storage_techs + ) + + slc.demand_input_name = f"{slc.commodity}_demand" + + slc.techs_to_commodities = slc_config["tech_to_commodity"] + + slc.multi_commodity_system = ( + True if len({e[-1] for e in slc.techs_to_commodities}) > 1 else False + ) + return slc + + +# Test methods in _post_setup_multi_commodity +# _find_converter_techs() +# _make_conversion_factor_recipes() + + +@pytest.mark.unit +def test_find_converter_techs_fake_system(subtests): + # Test methods in _post_setup_multi_commodity + # _find_converter_techs() + tech_connections = [ + ["ocean", "desalination", "salt_water", ""], + ["desalination", "electrolyzer", "fresh_water", ""], + ["wind", "elec_combiner", "electricity", ""], + ["solar", "elec_combiner", "electricity", ""], + ["elec_combiner", "battery", "electricity", ""], + ["battery", "elec_combiner_2", "electricity", ""], + ["elec_combiner", "elec_combiner_2", "electricity", ""], + ["elec_combiner_2", "electrolyzer", "electricity", ""], + ["electrolyzer", "h2_storage", "hydrogen", ""], + ["electrolyzer", "h2_combiner", "hydrogen", ""], + ["electrolyzer", "haber_bosch", "oxygen", ""], + ["electrolyzer", "haber_bosch", "heat", ""], + ["h2_storage", "h2_combiner", "hydrogen", ""], + ["h2_combiner", "haber_bosch", "hydrogen", ""], + ["grid", "haber_bosch", "electricity", ""], + ["n2_feedstock", "haber_bosch", "nitrogen", ""], + ["haber_bosch", "nh3_storage", "ammonia", ""], + ["haber_bosch", "nh3_combiner", "ammonia", ""], + ["nh3_storage", "nh3_combiner", "ammonia", ""], + ["nh3_combiner", "nh3_load_demand", "ammonia", ""], + ] + + example_folder = EXAMPLE_DIR / "35_system_level_control" / "nh3_with_storage" + plant_config = load_plant_yaml(example_folder / "plant_config.yaml") + tech_config = load_tech_yaml(example_folder / "tech_config.yaml") + + plant_config["technology_interconnections"] = tech_connections + extra_tech_config_keys = { + k[0]: {} for k in tech_connections if k[0] not in tech_config["technologies"] + } + extra_tech_config_keys |= { + k[1]: {} for k in tech_connections if k[1] not in tech_config["technologies"] + } + tech_config_fake = tech_config["technologies"] | extra_tech_config_keys + + slc = make_and_setup_slc_baseclass(plant_config, {"technologies": tech_config_fake}) + + converters, converter_upstreams = slc._find_converter_techs() + pem_output_cmod = ["hydrogen", "heat", "oxygen"] + electrolyzer_conversions = [("electricity", "electrolyzer", c) for c in pem_output_cmod] + electrolyzer_conversions += [("fresh_water", "electrolyzer", c) for c in pem_output_cmod] + with subtests.test("Converter elements for electrolyzer"): + assert all(k in converters for k in electrolyzer_conversions) + + hb_conversions = [(c, "haber_bosch", "ammonia") for c in pem_output_cmod] + with subtests.test("Converter elements for haber_bosch"): + assert all(k in converters for k in hb_conversions) + + expected_converter_upstreams = { + ("salt_water", "desalination"): ["ocean"], + ("hydrogen", "haber_bosch"): ["h2_combiner", "electrolyzer", "h2_storage"], + ("heat", "haber_bosch"): ["electrolyzer"], + ("oxygen", "haber_bosch"): ["electrolyzer"], + ("nitrogen", "haber_bosch"): ["n2_feedstock"], + ("electricity", "haber_bosch"): ["grid"], + ("fresh_water", "electrolyzer"): ["desalination"], + ("electricity", "electrolyzer"): [ + "elec_combiner", + "elec_combiner_2", + "battery", + "solar", + "wind", + ], + } + + mismatched_upstreams = [ + k for k, v in expected_converter_upstreams.items() if set(converter_upstreams[k]) != set(v) + ] + with subtests.test("Converter upstreams"): + assert len(mismatched_upstreams) == 0 + + slc._post_setup_multi_commodity() + + with subtests.test("12 distinct conversions"): + assert len(slc.rename_me_config.converters) == len(converters) + assert len(converters) == 12 + + with subtests.test("9 groups"): + assert len(slc.rename_me_config.grouped_techs) == 9 + + with subtests.test("17 conversion recipes"): + assert len(slc.rename_me_config.conversion_recipes) == 17 + + conversion_to_group_name = {} + for input_cmod, tech, output_cmod in list(slc.rename_me_config.converters): + group_name = [ + k[2][1] + for k in slc.rename_me_config.conversion_recipes + if k[0] == output_cmod and k[1] == input_cmod + ] + conversion_to_group_name[(input_cmod, tech, output_cmod)] = list(set(group_name))[0] + + grouped_tech_in_degrees = dict(slc.rename_me_config.simple_graph.in_degree()) + grouped_tech_out_degrees = dict(slc.rename_me_config.simple_graph.out_degree()) + electrolyzer_in_degs = [ + grouped_tech_in_degrees[conversion_to_group_name[con]] for con in electrolyzer_conversions + ] + with subtests.test("electrolyzer in degrees"): + assert all(k == 2 for k in electrolyzer_in_degs) + hb_in_degs = [grouped_tech_in_degrees[conversion_to_group_name[con]] for con in hb_conversions] + with subtests.test("haber bosch in degrees"): + assert all(k == 5 for k in hb_in_degs) + multi_output_grouped_techs = [k for k, v in grouped_tech_out_degrees.items() if v > 1] + electrolyzer_groups = {conversion_to_group_name[con] for con in electrolyzer_conversions} + with subtests.test("3 electrolyzer groups"): + assert len(electrolyzer_groups) == 3 + with subtests.test("out degrees are <1 except for electrolyzer upstream"): + for source_group in multi_output_grouped_techs: + assert all( + slc.rename_me_config.simple_graph.has_edge(source_group, dest_group) + for dest_group in electrolyzer_groups + ) + assert grouped_tech_out_degrees[source_group] == 3 + + recipe_names_long = [ + k + for k in slc.rename_me_config.conversion_recipes + if k[0] == "fresh_water" and k[1] == "salt_water" + ] + r0_partial = [("ammonia", "nh3_storage", "ammonia"), ("ammonia", "nh3_combiner", "ammonia")] + r2 = [("salt_water", "desalination", "fresh_water")] + with subtests.test("desalination recipes"): + for recipe_name in recipe_names_long: + recipe = slc.rename_me_config.conversion_recipes[recipe_name] + + cmod_diff = list(set(recipe[0]) - set(r0_partial)) + cmod = cmod_diff[0][0] + assert cmod in pem_output_cmod + assert ("fresh_water", "electrolyzer", cmod) in recipe[1] + assert recipe[2] == r2 + + +@pytest.mark.unit +def test_find_converter_techs_nh3_system(subtests): + # Test methods in _post_setup_multi_commodity + # _find_converter_techs() + example_folder = EXAMPLE_DIR / "35_system_level_control" / "nh3_with_storage" + plant_config = load_plant_yaml(example_folder / "plant_config.yaml") + tech_config = load_tech_yaml(example_folder / "tech_config.yaml") + slc = make_and_setup_slc_baseclass(plant_config, tech_config) + + # Test _find_converter_techs() + converters, converter_upstreams = slc._find_converter_techs() + + expected_converters = { + ("nitrogen", "haber_bosch", "ammonia"), + ("electricity", "haber_bosch", "ammonia"), + ("hydrogen", "haber_bosch", "ammonia"), + ("electricity", "electrolyzer", "hydrogen"), + } + + expected_converter_upstreams = { + ("electricity", "electrolyzer"): {"solar", "battery", "wind"}, + ("hydrogen", "haber_bosch"): {"electrolyzer", "h2_storage"}, + ("electricity", "haber_bosch"): {"electricity_feedstock"}, + ("nitrogen", "haber_bosch"): {"n2_feedstock"}, + } + + with subtests.test("converters"): + assert converters == expected_converters + with subtests.test("converter_upstreams"): + for k, v in expected_converter_upstreams.items(): + input_tech_upstreams = set(converter_upstreams.get(k)) & ( + set(slc.input_techs) | set(slc.feedstock_comps) + ) + assert input_tech_upstreams == v + + # Test _find_demand_tech_group() + # non_converter_input_techs_in_group, demand_group = slc._find_demand_tech_group( + # converters, converter_upstreams + # ) + + # with subtests.test("non main techs in demand group"): + # assert non_converter_input_techs_in_group == ["nh3_storage"] + + # expected_demand_group = {"ammonia-5": {"nh3_combiner", "haber_bosch", "nh3_storage"}} + # with subtests.test("demand_group"): + # assert demand_group == expected_demand_group + + +@pytest.mark.unit +def test_multi_commodity_post_setup_nh3_system(subtests): + # Test methods in _post_setup_multi_commodity + # _find_converter_techs() + # _make_conversion_factor_recipes() + example_folder = EXAMPLE_DIR / "35_system_level_control" / "nh3_with_storage" + plant_config = load_plant_yaml(example_folder / "plant_config.yaml") + tech_config = load_tech_yaml(example_folder / "tech_config.yaml") + slc_config = make_slc_topology(plant_config, tech_config) + + prob = om.Problem() + + feedstock_techs = [ + k for k, v in slc_config["tech_control_classifiers"].items() if v == "feedstock" + ] + feedstock_subsystem_names = [] + for fi, feedstock_tech in enumerate(feedstock_techs): + feedstock_commodity = [ + e[-1] for e in slc_config["tech_to_commodity"] if e[0] == feedstock_tech + ] + feedstock_comp = prob.model.add_subsystem(f"IVC{fi}", om.Group()) + feedstock_comp.add_subsystem( + "feedstock", + om.IndepVarComp( + name=f"{feedstock_tech}_{feedstock_commodity[0]}_out", + val=np.full(plant_config["plant"]["simulation"]["n_timesteps"], 1e9), + units="MMBtu/h", + ), + ) + + feedstock_subsystem_names.append( + f"IVC{fi}.feedstock.{feedstock_tech}_{feedstock_commodity[0]}_out" + ) + + slc = SystemLevelControlBase( + plant_config=plant_config, + tech_config=tech_config, + driver_config={}, + slc_topology=slc_config, + ) + prob.model.add_subsystem("slc", slc) + + for feedstock_name in feedstock_subsystem_names: + connection_destination = feedstock_name.split(".")[-1] + prob.model.connect(feedstock_name, f"slc.{connection_destination}") + + prob.setup() + + # Check converters + expected_converters = { + ("nitrogen", "haber_bosch", "ammonia"), + ("electricity", "haber_bosch", "ammonia"), + ("hydrogen", "haber_bosch", "ammonia"), + ("electricity", "electrolyzer", "hydrogen"), + } + + converters = prob.model.slc.rename_me_config.converters + + with subtests.test("converters"): + assert converters == expected_converters + + # Check converter_upstreams + expected_converter_upstreams = { + ("electricity", "electrolyzer"): {"solar", "battery", "wind"}, + ("hydrogen", "haber_bosch"): {"electrolyzer", "h2_storage"}, + ("electricity", "haber_bosch"): {"electricity_feedstock"}, + ("nitrogen", "haber_bosch"): {"n2_feedstock"}, + # ("ammonia", "nh3_load_demand"): {"haber_bosch", "nh3_storage"}, + } + + converter_upstreams = prob.model.slc.rename_me_config.converter_upstreams + # with subtests.test("converter upstreams"): + for k, v in expected_converter_upstreams.items(): + with subtests.test(f"converter upstreams {k}"): + input_tech_upstreams = set(converter_upstreams.get(k)) & ( + set(slc.input_techs) | set(slc.feedstock_comps) + ) + assert input_tech_upstreams == v + + # Check simple_graph + simple_graph = prob.model.slc.rename_me_config.simple_graph + edges = list(simple_graph.edges(data="commodity")) + with subtests.test("number of edges"): + assert len(edges) == 5 + non_numbered_edges = [(k[0].split("-")[0], k[1].split("-")[0], k[2]) for k in edges] + expected_edges = [ + ("electricity", "hydrogen", "electricity"), + ("hydrogen", "ammonia", "hydrogen"), + ("ammonia", "nh3_load_demand", "ammonia"), + ("nitrogen", "ammonia", "nitrogen"), + ("electricity", "ammonia", "electricity"), + ] + + with subtests.test("simple_graph edges"): + # assert not bool(set(edges) ^ set(expected_edges)) + assert set(non_numbered_edges) == set(expected_edges) + + # Check grouped_techs + grouped_techs = prob.model.slc.rename_me_config.grouped_techs + expected_groups = [ + {"solar", "battery", "wind", "combiner", "elec_combiner"}, + {"electrolyzer", "h2_storage", "h2_combiner"}, + {"electricity_feedstock"}, + {"n2_feedstock"}, + {"nh3_combiner", "haber_bosch", "nh3_storage"}, + ] + failed_groups = [] + for group, techs_in_group in grouped_techs.items(): + if not any(g == set(techs_in_group) for g in expected_groups): + failed_groups.append(group) + with subtests.test("Grouped technologies is correct"): + assert len(failed_groups) == 0 + + # Check conversion_recipes + conversion_recipes_list = prob.model.slc.rename_me_config.conversion_recipes + conversion_recipes = {} + for k, v in conversion_recipes_list.items(): + v_as_set = [set(vi) for vi in v] + conversion_recipes[k] = v_as_set + demand_group_general = [ + ("ammonia", "nh3_storage", "ammonia"), + ("ammonia", "nh3_combiner", "ammonia"), + ] + + n2_nh3_recipe = [("nitrogen", "haber_bosch", "ammonia"), *demand_group_general] + n2_nh3_recipe_name = [k for k in conversion_recipes if k[0] == "ammonia" and k[1] == "nitrogen"] + with subtests.test("Nitrogen to Ammonia Recipe"): + assert conversion_recipes[n2_nh3_recipe_name[0]] == [set(n2_nh3_recipe)] + + electricity_nh3_recipe = [("electricity", "haber_bosch", "ammonia"), *demand_group_general] + electricity_nh3_recipe_name = [ + k for k in conversion_recipes if k[0] == "ammonia" and k[1] == "electricity" + ] + with subtests.test("Electricity to Ammonia Recipe"): + assert conversion_recipes[electricity_nh3_recipe_name[0]] == [set(electricity_nh3_recipe)] + + h2_nh3_recipe = [("hydrogen", "haber_bosch", "ammonia"), *demand_group_general] + h2_nh3_recipe_name = [k for k in conversion_recipes if k[0] == "ammonia" and k[1] == "hydrogen"] + + with subtests.test("Hydrogen to Ammonia Recipe"): + assert conversion_recipes[h2_nh3_recipe_name[0]] == [set(h2_nh3_recipe)] + + h2_elec_subrecipe = { + ("hydrogen", "h2_storage", "hydrogen"), + ("electricity", "electrolyzer", "hydrogen"), + ("hydrogen", "h2_combiner", "hydrogen"), + } + h2_elec_recipe = [set(h2_nh3_recipe), h2_elec_subrecipe] + + h2_elec_recipe_name = [ + k for k in conversion_recipes if k[0] == "hydrogen" and k[1] == "electricity" + ] + with subtests.test("Electricity for Hydrogen Recipe"): + assert conversion_recipes[h2_elec_recipe_name[0]] == h2_elec_recipe + + with subtests.test("4 recipes"): + assert len(conversion_recipes) == 4 + + # Check non_converter_conversion_factor_keys + non_converter_keys = prob.model.slc.rename_me_config.non_converter_conversion_factor_keys + non_converter_techs = [k[1] for k in non_converter_keys] + expected_non_converter_techs = [ + "nh3_storage", + "battery", + "n2_feedstock", + "wind", + "electricity_feedstock", + "solar", + "h2_storage", + "elec_combiner", + "combiner", + "h2_combiner", + "nh3_combiner", + ] + with subtests.test("Non converter techs"): + assert set(non_converter_techs) == set(expected_non_converter_techs) + with subtests.test("wind key"): + assert ("electricity", "wind", "electricity") in non_converter_keys + with subtests.test("n2_feedstock key"): + assert ("nitrogen", "n2_feedstock", "nitrogen") in non_converter_keys + with subtests.test("h2_combiner key"): + assert ("hydrogen", "h2_combiner", "hydrogen") in non_converter_keys + with subtests.test("nh3_storage key"): + assert ("ammonia", "nh3_storage", "ammonia") in non_converter_keys + + converter_tech_names = prob.model.slc.rename_me_config.converter_tech_names + with subtests.test("Converter tech names"): + assert converter_tech_names == {"haber_bosch", "electrolyzer"} + + +# Test methods used by Demand Following +# `get_converter_capacity_ratio` +# `get_converter_conversion_ratio` +# `get_conversion_from_recipe` +# `get_techs_to_demand_from_recipe` + + +@pytest.mark.unit +def test_multi_commodity_conversion_factor_nh3_system(subtests): + # Test methods available in SLC baseclass that are not used directly within SLC baseclass + + # --- Same setup as ``test_multi_commodity_post_setup_nh3_system`` --- + + example_folder = EXAMPLE_DIR / "35_system_level_control" / "nh3_with_storage" + plant_config = load_plant_yaml(example_folder / "plant_config.yaml") + tech_config = load_tech_yaml(example_folder / "tech_config.yaml") + slc_config = make_slc_topology(plant_config, tech_config) + + prob = om.Problem() + + feedstock_techs = [ + k for k, v in slc_config["tech_control_classifiers"].items() if v == "feedstock" + ] + feedstock_subsystem_names = [] + for fi, feedstock_tech in enumerate(feedstock_techs): + feedstock_commodity = [ + e[-1] for e in slc_config["tech_to_commodity"] if e[0] == feedstock_tech + ] + feedstock_comp = prob.model.add_subsystem(f"IVC{fi}", om.Group()) + feedstock_comp.add_subsystem( + "feedstock", + om.IndepVarComp( + name=f"{feedstock_tech}_{feedstock_commodity[0]}_out", + val=np.full(plant_config["plant"]["simulation"]["n_timesteps"], 1e9), + units="MMBtu/h", + ), + ) + + feedstock_subsystem_names.append( + f"IVC{fi}.feedstock.{feedstock_tech}_{feedstock_commodity[0]}_out" + ) + + slc = SystemLevelControlBase( + plant_config=plant_config, + tech_config=tech_config, + driver_config={}, + slc_topology=slc_config, + ) + prob.model.add_subsystem("slc", slc) + + for feedstock_name in feedstock_subsystem_names: + connection_destination = feedstock_name.split(".")[-1] + prob.model.connect(feedstock_name, f"slc.{connection_destination}") + + prob.setup() + # --------------------------- End of setup --------------------------- + h2_storage_profile = np.tile( + np.concatenate([np.arange(-5.0, 6.0, 1), np.arange(6.0, -5, -1)]), 399 + )[:8760] + fake_inputs = { + "wind_rated_electricity_production": np.array([30.0]), + "wind_electricity_out": np.full(8760, 15.0), + "solar_rated_electricity_production": np.array([25.0]), + "solar_electricity_out": np.full(8760, 20.0), + "battery_rated_electricity_production": np.array([16.0]), + "battery_electricity_out": np.zeros(8760), + "electrolyzer_rated_hydrogen_production": np.array([71.0]), + "electrolyzer_hydrogen_out": np.full(8760, 39.0), + "h2_storage_rated_hydrogen_production": np.array([14.0]), + "h2_storage_hydrogen_out": h2_storage_profile, + "haber_bosch_rated_ammonia_production": np.array([50.0]), + "haber_bosch_ammonia_out": np.full(8760, 40), + "nh3_storage_rated_ammonia_production": np.array([4.0]), + "nh3_storage_ammonia_out": np.tile(np.array([-1, 1]), 4380), + "n2_feedstock_nitrogen_out": np.full(8760, 2.5), + "electricity_feedstock_electricity_out": np.full(8760, 13.0), + } + + # Test `get_converter_capacity_ratio` and `get_converter_conversion_ratio` + # Electricity to hydrogen + elec_per_h2_ratio = prob.model.slc.get_converter_conversion_ratio( + fake_inputs, "electricity", "hydrogen", "electrolyzer", ["battery", "wind", "solar"] + ) + elec_per_h2_capac_ratio = prob.model.slc.get_converter_capacity_ratio( + fake_inputs, "electricity", "hydrogen", "electrolyzer", ["battery", "wind", "solar"] + ) + elec_capac = 30.0 + 25.0 + 16.0 + elec_gen = 15.0 + 20.0 + with subtests.test("Electricity/Hydrogen conversion ratio"): + assert pytest.approx(elec_gen / 39.0, rel=1e-6) == elec_per_h2_ratio.mean() + with subtests.test("Electricity/Hydrogen capacity ratio"): + assert pytest.approx(elec_capac / 71.0, rel=1e-6) == elec_per_h2_capac_ratio + + # Hydrogen to ammonia + h2_per_nh3_ratio = prob.model.slc.get_converter_conversion_ratio( + fake_inputs, "hydrogen", "ammonia", "haber_bosch", ["electrolyzer", "h2_storage"] + ) + h2_per_nh3_capac_ratio = prob.model.slc.get_converter_capacity_ratio( + fake_inputs, "hydrogen", "ammonia", "haber_bosch", ["electrolyzer", "h2_storage"] + ) + h2_capac = 71.0 + 14.0 + h2_gen = h2_storage_profile + np.full(8760, 39.0) + with subtests.test("Hydrogen/Ammonia conversion ratio"): + assert pytest.approx((h2_gen / 40).mean(), rel=1e-6) == h2_per_nh3_ratio.mean() + with subtests.test("Hydrogen/Ammonia capacity ratio"): + assert pytest.approx(h2_capac / 50.0, rel=1e-6) == h2_per_nh3_capac_ratio + + # Nitrogen to ammonia + n2_per_nh3_ratio = prob.model.slc.get_converter_conversion_ratio( + fake_inputs, "nitrogen", "ammonia", "haber_bosch", ["n2_feedstock"] + ) + n2_per_nh3_capac_ratio = prob.model.slc.get_converter_capacity_ratio( + fake_inputs, "nitrogen", "ammonia", "haber_bosch", ["n2_feedstock"] + ) + with subtests.test("Nitrogen/Ammonia conversion ratio"): + assert pytest.approx(2.5 / 40, rel=1e-6) == n2_per_nh3_ratio.mean() + with subtests.test("Nitrogen/Ammonia capacity ratio"): + assert pytest.approx(2.5 / 50.0, rel=1e-6) == n2_per_nh3_capac_ratio + + # Electricity to ammonia + elec_per_nh3_ratio = prob.model.slc.get_converter_conversion_ratio( + fake_inputs, "electricity", "ammonia", "haber_bosch", ["electricity_feedstock"] + ) + elec_per_nh3_capac_ratio = prob.model.slc.get_converter_capacity_ratio( + fake_inputs, "electricity", "ammonia", "haber_bosch", ["electricity_feedstock"] + ) + with subtests.test("Electricity/Ammonia conversion ratio"): + assert pytest.approx(13.0 / 40, rel=1e-6) == elec_per_nh3_ratio.mean() + with subtests.test("Electricity/Ammonia capacity ratio"): + assert pytest.approx(13.0 / 50.0, rel=1e-6) == elec_per_nh3_capac_ratio + + # Test `get_conversion_from_recipe` and `get_techs_to_demand_from_recipe` + conversion_factors = { + ("electricity", "electrolyzer", "hydrogen"): elec_gen / 39.0, + ("hydrogen", "haber_bosch", "ammonia"): (h2_gen / 40).mean(), + ("nitrogen", "haber_bosch", "ammonia"): 2.5 / 40, + ("electricity", "haber_bosch", "ammonia"): 13.0 / 40, + } + non_converter_keys = prob.model.slc.rename_me_config.non_converter_conversion_factor_keys + non_converter_factor = 1.0 + non_converter_conversion_factors = dict( + zip(non_converter_keys, [non_converter_factor] * len(non_converter_keys)) + ) + all_conversion_factors = conversion_factors | non_converter_conversion_factors + + conversion_recipes = prob.model.slc.rename_me_config.conversion_recipes + + # Nitrogen/Ammonia + n2_recipe_name = [k for k in conversion_recipes if k[0] == "ammonia" and k[1] == "nitrogen"][0] + # n2_recipe_name = ("ammonia", "nitrogen", "ammonia-5") + with subtests.test("Nitrogen/Ammonia Conversion Factor"): + conversion_factor = prob.model.slc.get_conversion_from_recipe( + all_conversion_factors, conversion_recipes[n2_recipe_name] + ) + assert pytest.approx(2.5 / 40.0, rel=1e-6) == conversion_factor + with subtests.test("Nitrogen/Ammonia Techs"): + techs_to_demand = prob.model.slc.get_techs_to_demand_from_recipe(n2_recipe_name) + assert ["n2_feedstock"] == techs_to_demand + + # Electricity/Ammonia + elec_recipe_name = [ + k for k in conversion_recipes if k[0] == "ammonia" and k[1] == "electricity" + ][0] + # elec_recipe_name = ("ammonia", "electricity", "ammonia-5") + with subtests.test("Electricity/Ammonia Conversion Factor"): + conversion_factor = prob.model.slc.get_conversion_from_recipe( + all_conversion_factors, conversion_recipes[elec_recipe_name] + ) + assert pytest.approx(13.0 / 40.0, rel=1e-6) == conversion_factor + with subtests.test("Electricity/Ammonia Techs"): + techs_to_demand = prob.model.slc.get_techs_to_demand_from_recipe(elec_recipe_name) + assert ["electricity_feedstock"] == techs_to_demand + + # Hydrogen/Ammonia + h2_recipe_name = [k for k in conversion_recipes if k[0] == "ammonia" and k[1] == "hydrogen"][0] + # h2_recipe_name = ("ammonia", "hydrogen", "ammonia-5") + with subtests.test("Hydrogen/Ammonia Conversion Factor"): + conversion_factor = prob.model.slc.get_conversion_from_recipe( + all_conversion_factors, conversion_recipes[h2_recipe_name] + ) + assert pytest.approx((h2_gen / 40).mean(), rel=1e-6) == conversion_factor + + with subtests.test("Hydrogen/Ammonia Techs"): + techs_to_demand = prob.model.slc.get_techs_to_demand_from_recipe(h2_recipe_name) + expected_techs = ["h2_storage", "electrolyzer", "h2_combiner"] + assert set(expected_techs) == set(techs_to_demand) + + # Electricity/Hydrogen/Ammonia + eh2_recipe_name = [ + k for k in conversion_recipes if k[0] == "hydrogen" and k[1] == "electricity" + ][0] + # eh2_recipe_name = ("hydrogen", "electricity", "hydrogen-1") + with subtests.test("Electricity/Hydrogen/Ammonia Conversion Factor"): + conversion_factor = prob.model.slc.get_conversion_from_recipe( + all_conversion_factors, conversion_recipes[eh2_recipe_name] + ) + expected_conversion_factor = (h2_gen / 40).mean() * (elec_gen / 39.0) + assert pytest.approx(expected_conversion_factor, rel=1e-6) == conversion_factor + with subtests.test("Electricity/Hydrogen/Ammonia Techs"): + techs_to_demand = prob.model.slc.get_techs_to_demand_from_recipe(eh2_recipe_name) + expected_techs = ["battery", "wind", "solar", "combiner", "elec_combiner"] + assert set(expected_techs) == set(techs_to_demand) diff --git a/h2integrate/control/control_strategies/system_level/test/test_slc_controllers.py b/h2integrate/control/control_strategies/system_level/test/test_slc_controllers.py index cfd46a14c..545dd259f 100644 --- a/h2integrate/control/control_strategies/system_level/test/test_slc_controllers.py +++ b/h2integrate/control/control_strategies/system_level/test/test_slc_controllers.py @@ -2,9 +2,9 @@ import numpy as np import pytest -import networkx as nx import openmdao.api as om +from h2integrate import H2IntegrateModel from h2integrate.control.control_strategies.system_level.demand_following_control import ( DemandFollowingControl, ) @@ -38,14 +38,8 @@ def _build_plant_config( def _build_technology_graph(technology_interconnections): - technology_graph = nx.DiGraph() - for connection in technology_interconnections: - source = connection[0] - destination = connection[1] - if len(connection) == 4: - technology_graph.add_edge(source, destination, commodity=connection[2]) - else: - technology_graph.add_edge(source, destination) + model = object.__new__(H2IntegrateModel) + technology_graph = model.create_technology_graph(technology_interconnections) return technology_graph @@ -68,9 +62,10 @@ def _build_slc_topology( demand_commodity_rate_units: str = "kW", storage_techs_with_control: list = [], ): - sources_to_commodities = { - (e[0], e[-1]) for e in technology_graph.edges(data="commodity") if e[-1] is not None - } + sources_to_commodities = set() + for source, _, commodities in technology_graph.edges(data="commodity"): + if commodities is not None: + sources_to_commodities.update((source, commodity) for commodity in commodities) tech_to_commodities = { (e[0], e[-1]) for e in sources_to_commodities if e[0] in tech_control_classifiers diff --git a/h2integrate/control/control_strategies/system_level/test/test_slc_examples.py b/h2integrate/control/control_strategies/system_level/test/test_slc_examples.py index e97e63bb2..ed42c9853 100644 --- a/h2integrate/control/control_strategies/system_level/test/test_slc_examples.py +++ b/h2integrate/control/control_strategies/system_level/test/test_slc_examples.py @@ -185,9 +185,9 @@ def test_slc_yes_hydrogen(subtests, temp_copy_of_example): with subtests.test("LCOH"): assert ( pytest.approx( - model.prob.get_val("finance_subgroup_hydrogen.LCOH", units="USD/kg"), rel=1e-6 + model.prob.get_val("finance_subgroup_hydrogen.LCOH", units="USD/kg")[0], rel=1e-6 ) - == 14.878096642042243 + == 14.46645200483752 ) @@ -403,3 +403,231 @@ def test_slc_upstream_demand(subtests, temp_copy_of_example): # check that no hydrogen systems are in model.prob.get_val(slc_h2s_output_var, units="kg/h") assert f"Variable '{slc_h2s_output_var}' not found. " in str(excinfo.value) + + +@pytest.mark.integration +@pytest.mark.parametrize( + "example_folder,resource_example_folder", + [("35_system_level_control/complex_multi_commodity", None)], +) +def test_slc_complex_multi_commodity_v1(subtests, temp_copy_of_example): + ex_folder = temp_copy_of_example + + h2i = H2IntegrateModel(ex_folder / "top_level_config.yaml") + + h2i.setup() + + h2i.run() + + with subtests.test("LCOH"): + assert ( + pytest.approx(4.064419131023322, rel=1e-6) + == h2i.model.get_val("finance_subgroup_h2.LCOH", units="USD/kg")[0] + ) + + with subtests.test("LCOA - Produced"): + assert ( + pytest.approx(1.306352207437524, rel=1e-6) + == h2i.model.get_val("finance_subgroup_nh3_produced.LCOA", units="USD/kg")[0] + ) + + with subtests.test("LCOA - Delivered"): + assert ( + pytest.approx(1.404495041524232, rel=1e-6) + == h2i.model.get_val("finance_subgroup_nh3_delivered.LCOA", units="USD/kg")[0] + ) + + with subtests.test("Unmet Ammonia Demand"): + assert ( + pytest.approx(102882.4504724315, rel=1e-6) + == h2i.model.get_val("nh3_load_demand.unmet_ammonia_demand_out", units="t/h").sum() + ) + + with subtests.test("Ammonia Demand Capacity Factor"): + assert ( + pytest.approx(75.27450203676736, rel=1e-6) + == h2i.model.get_val("nh3_load_demand.capacity_factor", units="percent")[0] + ) + + with subtests.test("Wind electricity set point (flexible tech)"): + wind_set_point = h2i.model.get_val( + "system_level_controller.wind_electricity_set_point", units="kW" + ) + wind_capacity = h2i.model.get_val( + "system_level_controller.wind_rated_electricity_production", units="kW" + ) + assert np.all(wind_set_point == wind_capacity) + with subtests.test("Solar electricity set point (flexible tech)"): + solar_set_point = h2i.model.get_val( + "system_level_controller.solar_electricity_set_point", units="kW" + ) + solar_capacity = h2i.model.get_val( + "system_level_controller.solar_rated_electricity_production", units="kW" + ) + assert np.all(solar_set_point == solar_capacity) + with subtests.test("Battery electricity set point (storage)"): + np.testing.assert_allclose( + h2i.model.get_val("system_level_controller.battery_electricity_set_point", units="kW"), + np.full(8760, 1119087.470593775), + rtol=1e-6, + atol=1e-6, + ) + with subtests.test("Electrolyzer hydrogen set point (dispatchable)"): + pem_set_point = h2i.model.get_val( + "system_level_controller.electrolyzer_hydrogen_set_point", units="kg/h" + ) + np.testing.assert_allclose( + pem_set_point, np.full(8760, 11744.704419585461), rtol=1e-6, atol=1e-6 + ) + with subtests.test("H2 Storage hydrogen set point (storage)"): + h2s_set_point_max = h2i.model.get_val( + "system_level_controller.h2_storage_hydrogen_set_point", units="kg/h" + ).max() + h2s_set_point_min = h2i.model.get_val( + "system_level_controller.h2_storage_hydrogen_set_point", units="kg/h" + ).min() + assert pytest.approx(11744.704419585461, rel=1e-6) == h2s_set_point_min + assert pytest.approx(23489.408839170923, rel=1e-6) == h2s_set_point_max + + with subtests.test("Haber Bosch: ammonia set point"): + assert np.all( + h2i.model.get_val("system_level_controller.haber_bosch_ammonia_set_point", units="kg/h") + == 47499.84 + ) + with subtests.test("Haber Bosch: hydrogen consumption"): + assert ( + pytest.approx(67359635.13205291, rel=1e-6) + == h2i.model.get_val("haber_bosch.hydrogen_consumed", units="kg/h").sum() + ) + with subtests.test("Haber Bosch: nitrogen consumption"): + assert ( + pytest.approx(312006087.72971725, rel=1e-6) + == h2i.model.get_val("haber_bosch.nitrogen_consumed", units="kg/h").sum() + ) + with subtests.test("Haber Bosch: electricity consumption"): + assert ( + pytest.approx(178149.21787066793, rel=1e-6) + == h2i.model.get_val("haber_bosch.electricity_consumed", units="MW").sum() + ) + + +@pytest.mark.integration +@pytest.mark.parametrize( + "example_folder,resource_example_folder", [("35_system_level_control/nh3_with_storage", None)] +) +def test_slc_complex_nh3_with_storage(subtests, temp_copy_of_example): + ex_folder = temp_copy_of_example + + h2i = H2IntegrateModel(ex_folder / "top_level_config.yaml") + + h2i.setup() + + h2i.run() + + with subtests.test("LCOH"): + assert pytest.approx(3.8867863862476097, rel=1e-6) == h2i.model.get_val( + "finance_subgroup_h2.LCOH", units="USD/kg" + ) + + with subtests.test("LCOA - Produced"): + assert pytest.approx(1.2054326039224676, rel=1e-6) == h2i.model.get_val( + "finance_subgroup_nh3_produced.LCOA", units="USD/kg" + ) + + with subtests.test("LCOA - Available"): + assert pytest.approx(1.207172852022115, rel=1e-6) == h2i.model.get_val( + "finance_subgroup_ammonia_available.LCOA", units="USD/kg" + ) + + with subtests.test("LCOA - Delivered"): + assert pytest.approx(1.3022897725787497, rel=1e-6) == h2i.model.get_val( + "finance_subgroup_nh3_delivered.LCOA", units="USD/kg" + ) + + with subtests.test("Unmet Ammonia Demand"): + assert ( + pytest.approx(79364.20686173553, rel=1e-6) + == h2i.model.get_val("nh3_load_demand.unmet_ammonia_demand_out", units="t/h").sum() + ) + + with subtests.test("Ammonia Demand Capacity Factor"): + assert pytest.approx(81.00662271288984, rel=1e-6) == h2i.model.get_val( + "nh3_load_demand.capacity_factor", units="percent" + ) + + with subtests.test("Wind electricity set point (flexible tech)"): + wind_set_point = h2i.model.get_val( + "system_level_controller.wind_electricity_set_point", units="kW" + ) + wind_capacity = h2i.model.get_val( + "system_level_controller.wind_rated_electricity_production", units="kW" + ) + assert np.all(wind_set_point == wind_capacity) + with subtests.test("Solar electricity set point (flexible tech)"): + solar_set_point = h2i.model.get_val( + "system_level_controller.solar_electricity_set_point", units="kW" + ) + solar_capacity = h2i.model.get_val( + "system_level_controller.solar_rated_electricity_production", units="kW" + ) + assert np.all(solar_set_point == solar_capacity) + with subtests.test("Battery electricity set point (storage)"): + np.testing.assert_allclose( + h2i.model.get_val("system_level_controller.battery_electricity_set_point", units="kW"), + np.full(8760, 1119087.470593775), + rtol=1e-6, + atol=1e-6, + ) + with subtests.test("Electrolyzer hydrogen set point (dispatchable)"): + pem_set_point = h2i.model.get_val( + "system_level_controller.electrolyzer_hydrogen_set_point", units="kg/h" + ) + np.testing.assert_allclose( + pem_set_point, np.full(8760, 11744.704419585461), rtol=1e-6, atol=1e-6 + ) + with subtests.test("H2 Storage hydrogen set point (storage)"): + h2s_set_point_max = h2i.model.get_val( + "system_level_controller.h2_storage_hydrogen_set_point", units="kg/h" + ).max() + h2s_set_point_min = h2i.model.get_val( + "system_level_controller.h2_storage_hydrogen_set_point", units="kg/h" + ).min() + assert pytest.approx(11744.704419585461, rel=1e-6) == h2s_set_point_min + assert pytest.approx(23489.408839170923, rel=1e-6) == h2s_set_point_max + + with subtests.test("Haber Bosch: ammonia set point"): + assert np.all( + h2i.model.get_val("system_level_controller.haber_bosch_ammonia_set_point", units="kg/h") + == 47700.0 + ) + with subtests.test("Haber Bosch: hydrogen consumption"): + assert ( + pytest.approx(73265794.3391077, rel=1e-6) + == h2i.model.get_val("haber_bosch.hydrogen_consumed", units="kg/h").sum() + ) + with subtests.test("Haber Bosch: nitrogen consumption"): + assert ( + pytest.approx(339363089.0568392, rel=1e-6) + == h2i.model.get_val("haber_bosch.nitrogen_consumed", units="kg/h").sum() + ) + with subtests.test("Haber Bosch: electricity consumption"): + assert ( + pytest.approx(182579.15194510925, rel=1e-6) + == h2i.model.get_val("haber_bosch.electricity_consumed", units="MW").sum() + ) + + with subtests.test("NH3 storage ammonia set point (storage)"): + nh3_set_point_max = h2i.model.get_val( + "system_level_controller.nh3_storage_ammonia_set_point", units="kg/h" + ).max() + nh3_set_point_min = h2i.model.get_val( + "system_level_controller.nh3_storage_ammonia_set_point", units="kg/h" + ).min() + assert pytest.approx(100700.0, rel=1e-6) == nh3_set_point_max + assert pytest.approx(47700.0, rel=1e-6) == nh3_set_point_min + + with subtests.test("NH3 storage: ammonia out"): + assert np.all(h2i.model.get_val("nh3_storage.ammonia_out", units="kg/h") == 0.0) + + with subtests.test("Ammonia Storage Command"): + assert np.all(h2i.model.get_val("nh3_storage.ammonia_command_value") == 0.0) diff --git a/h2integrate/core/h2integrate_model.py b/h2integrate/core/h2integrate_model.py index 91acb6e9d..27bda5914 100644 --- a/h2integrate/core/h2integrate_model.py +++ b/h2integrate/core/h2integrate_model.py @@ -593,7 +593,6 @@ def _classify_slc_technologies(self): for k, v in self.tech_control_classifiers.items() if k in upstream_controllable_techs } - # Check if storage models have a controller storage_tech_to_control = {} for tech, classifier in upstream_tech_control_classifiers.items():